# Rapid pre-explosion increase in dome extrusion rate at La Soufrière, St. Vincent quantified from synthetic aperture radar backscatter

E. W. Dualeh, S. K. Ebmeier, T. J. Wright, M. P. Poland, R. Grandin, A. J. Stinton, M. Camejo-Harry, B. Esse, M. Burton

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

E. W. Dualeh, S. K. Ebmeier, T. J. Wright, M. P. Poland, R. Grandin, et al. Rapid pre-explosion increase in dome extrusion rate at La Soufrière, St. Vincent quantified from synthetic aperture radar backscatter. Earth and Planetary Science Letters, 2023, 603, 10.1016/j.epsl.2022.117980. insu-04155771

HAL Id: insu-04155771

https://insu.hal.science/insu-04155771

Submitted on 7 Jul 2023

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Rapid pre-explosion increase in dome extrusion rate at La Soufrière, St. Vincent quantified from synthetic aperture radar backscatter

<figure><img src="image_4.png" /></figure>

E.W. Dualeh<sup>a,*</sup>, S.K. Ebmeier<sup>a</sup>, T.J. Wright<sup>a</sup>, M.P. Poland<sup>b</sup>, R. Grandin<sup>c</sup>, A.J. Stinton<sup>d,e</sup>, M. Camejo-Harry<sup>f,g</sup>, B. Esse<sup>g</sup>, M. Burton<sup>g</sup>

## ARTICLE INFO

**Article history:**

Received 22 September 2022

Received in revised form 20 December 2022

Accepted 27 December 2022

Available online 13 January 2023

Editor: C.M. Petrone

**Keywords:**

volcano remote sensing

Synthetic Aperture Radar (SAR) backscatter

dome growth

La Soufrière

St. Vincent

## ABSTRACT

The extrusion rate of a lava dome is a critical parameter for monitoring silicic eruptions and forecasting their development. Satellite radar backscatter can provide unique information about dome growth during a volcanic eruption when other datasets (e.g., optical, thermal, ground-based measurements, etc.) may be limited. Here, we present an approach for estimating volcanic topography from individual backscatter images. Using data from multiple SAR sensors we apply the method to the dome growth during the 2021 eruption at La Soufrière, St. Vincent. We measure an average extrusion rate of $1.8 \text{ m}^3\text{s}^{-1}$ between December 2020 and March 2021 before an acceleration in extrusion rate to $17.5 \text{ m}^3\text{s}^{-1}$ in the 2 days prior to the explosive eruption on 9 April 2021. We estimate a final dome volume of 19.4 million $\text{m}^3$, extrapolated from the SAR sensors, with approximately 15% of the total extruded volume emplaced in the last 2 days. A possible explanation for the acceleration in extrusion rate could be the combined emptying of a conduit and reservoir of older material before the ascent of gas-rich magma in April 2021.

© 2022 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

## 1. Introduction

After 41 years of quiet, La Soufrière Volcano, St. Vincent, began a new phase of extrusive dome growth on 27 December 2020. By the beginning of April 2021 the new dome had grown steadily around the southwest section of the 1979 dome (Fig. 1). Starting on 9 April 2021 the eruption suddenly transitioned into an explosive phase that lasted until 22 April 2021, completely reshaping the summit crater (Joseph et al., 2022). Understanding and forecasting such transitions from extrusive to explosive phases of eruptions is critical for hazards assessment and mitigation (Fink and Anderson, 2000; Joseph et al., 2022) and relies on robust monitoring data (e.g., Griffiths and Fink, 1997; Pallister et al., 2013; Scharff et al., 2014).

Transition from an effusive to explosive eruption occurs when the overpressure in a conduit exceeds confining pressure and has

been attributed to both shallow changes in a volcano's edifice ('top-down,' depressurisation) and events deeper in the magmatic system ('bottom-up,' increased overpressure). Some top-down triggers of explosive eruption, such as changes in the shallow stress regime caused by dome collapse (e.g., Pallister et al., 2013), high rainfall (e.g., Carn et al., 2004), or changes in extrusion direction (e.g., Watts et al., 2002) are detectable from ground- or satellite-based observations (e.g., Ryan et al., 2010; Pallister et al., 2013). Deeper triggers (Metcalfe et al., 2022), such as the crystallisation or degassing of an ascending magma, or an increase in conduit permeability, may not result in diagnostic changes at the surface. However, the resulting increases in magma ascent rate (e.g., Castro et al., 2012; Wang et al., 2015) can be inferred from (1) changes in seismicity (e.g., Nakada et al., 1999; Roman et al., 2016) or (2) increases in extrusion rate, if measurements are sufficiently frequent.

High-temporal-resolution measurements of dome growth are difficult to obtain from the ground, being limited by accessibility, resources, cloud coverage and rapidly changing hazards. Although optical aerial and satellite imagery can have sub-metre pixel resolution (e.g., Pléiades, WorldView), and one can make Digital El-

<footnote><sup>a</sup> COMET, School of Earth and Environment, University of Leeds, Leeds, LS2 9JT, UK</footnote>

<footnote><sup>b</sup> USGS, Yellowstone Volcano Observatory, Vancouver, WA, United States</footnote>

<footnote><sup>c</sup> Université de Paris, Institut de Physique du Globe de Paris, CNRS, Paris, France</footnote>

<footnote><sup>d</sup> Montserrat Volcano Observatory, Flemmings, Montserrat</footnote>

<footnote><sup>e</sup> Department of Earth Sciences, University of Oxford, Oxford, OX1 3AN, UK</footnote>

<footnote><sup>f</sup> The University of the West Indies, Seismic Research Centre, St. Augustine, Trinidad and Tobago</footnote>

<footnote><sup>g</sup> Department of Earth and Environmental Sciences, University of Manchester, Oxford Road, Manchester, M13 9PL, UK</footnote>

<footnote>* Corresponding author.
E-mail address: edna.dualeh@bristol.ac.uk (E.W. Dualeh).</footnote>
<figure><img src="image_5.png" /><figcaption>**Fig. 1.** Location and overview of La Soufrière. (A) map of the Lesser Antilles arc and the Atlantic and Caribbean plate subduction zone, (B) Sentinel-2 (S2) image (29-12-2019) of the island of St. Vincent with extent of SAR tracks shown in green and the summit of La Soufrière in red, (C) optical and (D-H) SAR images in radar geometry showing the dome growth in La Soufrière crater between December 2020 to April 2021. SAR satellite abbreviation: CSK, COSMO-SkyMed; TSX, TerraSAR-X; S1, Sentinel-1. (For interpretation of the colours in the figure(s), the reader is referred to the web version of this article.)</figcaption></figure>

elevation Models (DEM) directly with stereo optical imagery, these datasets depend on daylight and low cloud coverage. Similarly, the bistatic approach using TanDEM-X data to calculate measurements of elevation change over time (e.g., Poland, 2014; Albino et al., 2020) require a pair of images acquired at the same time. Historical time series of dome extrusion (e.g., Zharinov and Demyanchuk, 2008; Ryan et al., 2010; Harnett et al., 2019) therefore tend to have large gaps or rely on measurements averaged over weeks to months.

Synthetic Aperture Radar (SAR) sensors have the benefit of being independent of the time of day and meteorological conditions. SAR backscatter has previously been used to observe dome stability and growth (e.g., Pallister et al., 2013; Plank et al., 2019) and to track changing morphology (e.g., Wadge et al., 2011; Walter et al., 2015). Structural information such as dome heights and volumes have been extracted from SAR backscatter through the use of radar shadows (El Reventador, Ecuador, Arnold et al., 2017), feature tracking (Mount Cleveland, Alaska, Wang et al., 2015) and simulating SAR backscatter signals to construct synthetic domes (Shinmoe-dake, Japan, Ozawa and Kozono, 2013). The recent study by Angarita et al. (2022) demonstrates a method similar to that presented here to quantify changes at Shishaldin Volcano.

Here, we present an analysis of the 2021 dome growth at La Soufrière volcano on the island of St. Vincent using radar backscatter from multiple sensors. We describe an approach to extract topographic profiles from SAR backscatter images. Using data from multiple SAR satellite systems (i.e., TerraSAR-X, COSMO-SkyMed and Sentinel-1), we construct a time series of 32 dome volume estimates between December 2020 and April 2021. We show that the extrusion rate remained relatively constant until two days before the explosive eruption on 9 April 2021.

## 1.1. Soufrière St. Vincent

La Soufrière (1,234 m) is the current active volcanic centre on the main island of St. Vincent and is situated in the northern section of the island (Fig. 1). La Soufrière is considered one of the most active volcanoes in the Caribbean, with at least five explosive eruptions since the first recorded eruption in 1718 (Robertson, 1995).

La Soufrière is monitored by the University of the West Indies Seismic Research Centre (UWI-SRC), assisted by the Soufrière Monitoring Unit (SMU) of National Emergency Management Organisation (NEMO) of St. Vincent and the Grenadines. Historical eruptions were characterised by both explosive activity and effusive dome growth.

In December 2020, La Soufrière entered a new period of activity with a new dome growing on the edge of the 1979 dome first observed on 27 December 2020 (Joseph et al., 2022). The Government of St. Vincent raised the volcanic alert level in December to 'Orange' in response to the extrusive activity (NEMO, 2020a). By mid-February 2021 the dome had grown to approximately 90 m in height and had transitioned into an elongated coulee, curling around the southwest edge of the 1979 dome (Fig. 1). During the extrusive phase, 12 measurements of the 2021 lava dome were made by UWI-SRC and SMU through photogrammetry using ground-based (i.e., from the crater rim) and aerial (i.e., UAVs, helicopters and fixed-wing aircrafts) imagery. A volcano-tectonic (VT) earthquake swarm was observed on 23–25 March 2021, marking a change in eruptive style and suggesting movement of magma towards the surface (NEMO, 2020c; Joseph et al., 2022). The rate of earthquakes increased on 5–6 April 2021 with the first banded tremor appearing on the 8 April (Joseph et al., 2022); the volcanic alert level was raised to 'Red' (at 18:00 UTC on 8 April), triggering an evacuation order for sections of the island (NEMO, 2020b; Joseph et al., 2022).

At 12:41 UTC on 9 April 2021 La Soufrière started to erupt explosively, with multiple explosions occurring over the day. Over the following week the volcano continued to erupt, with pyroclastic flows, continuous venting, and ashfall over St. Vincent and on neighbouring Barbados (11 April 2021). The newly emplaced dome that had been growing since the beginning of the year, the 1979 dome, and much of the previously existing crater floor were completely destroyed. These changes were confirmed by the first observations of the volcano summit from SAR images acquired on 10 April 2021 from ICEYE at 02:03 UTC (Fig. S12) and Capella at 14:02 UTC (Joseph et al., 2022). Explosive activity ceased on 22 April 2021 at La Soufrière and the seismicity has remained low since then. The alert level was lowered to 'Orange' on 6 May 2021 (NEMO, 2020b; Joseph et al., 2022) and to 'Green' in March 2022.
<caption>**Table 1**</caption>

SAR acquisitions used to monitor St. Vincent Dome growth. Images are acquired by three sensors (COSMO-SkyMed, CSK; TerraSAR-X, TSX; Sentinel-1, S1) for five tracks (indicated by subscript).

<table><thead><tr><th>Acquisition Date</th><th>Sensor</th><th>Acquisition Date</th><th>Sensor</th><th>Acquisition Date</th><th>Sensor</th></tr></thead><tbody><tr><td>31 Dec. 2020 22:19</td><td><sup>1</sup>S1</td><td>24 Jan. 2021 22:19</td><td>TSX<sub>a</sub></td><td>09 Mar. 2021 22:19</td><td>TSX<sub>a</sub></td></tr><tr><td>02 Jan. 2021 22:19</td><td><sup>2</sup>TSX<sub>a</sub></td><td>30 Jan. 2021 22:18</td><td>S1</td><td>13 Mar. 2021 22:18</td><td>S1</td></tr><tr><td>08 Jan. 2021 22:10</td><td><sup>3</sup>TSX<sub>b</sub></td><td>04 Feb. 2021 22:19</td><td>TSX<sub>a</sub></td><td>19 Mar. 2021 22:18</td><td>S1</td></tr><tr><td>12 Jan. 2021 22:19</td><td>S1</td><td>05 Feb. 2021 22:19</td><td>S1</td><td>20 Mar. 2021 22:19</td><td>TSX<sub>a</sub></td></tr><tr><td>13 Jan. 2021 22:19</td><td>TSX<sub>a</sub></td><td>10 Feb. 2021 22:10</td><td>TSX<sub>b</sub></td><td>25 Mar. 2021 22:19</td><td>S1</td></tr><tr><td>17 Jan. 2021 21:38</td><td><sup>4</sup>CSK<sub>a</sub></td><td>11 Feb. 2021 22:18</td><td>S1</td><td>31 Mar. 2021 22:18</td><td>S1</td></tr><tr><td>18 Jan. 2021 21:38</td><td>CSK<sub>a</sub></td><td>15 Feb. 2021 22:19</td><td>TSX<sub>a</sub></td><td>31 Mar. 2021 22:19</td><td>TSX<sub>a</sub></td></tr><tr><td>18 Jan. 2021 22:18</td><td>S1</td><td>17 Feb. 2021 22:18</td><td>S1</td><td>06 Apr. 2021 22:19</td><td>S1</td></tr><tr><td>19 Jan. 2021 22:10</td><td>TSX<sub>b</sub></td><td>23 Feb. 2021 22:18</td><td>S1</td><td>07 Apr. 2021 21:38</td><td>CSK<sub>a</sub></td></tr><tr><td>22 Jan. 2021 21:32</td><td><sup>5</sup>CSK<sub>b</sub></td><td>01 Mar. 2021 22:18</td><td>S1</td><td>08 Apr. 2021 21:38</td><td>CSK<sub>a</sub></td></tr><tr><td>24 Jan. 2021 22:19</td><td>S1</td><td>07 Mar. 2021 22:18</td><td>S1</td><td></td><td></td></tr></tbody></table>

## 2. Data and methods

We develop a method that allows us to make measurements of new volcanic topography from single SAR backscatter images on the assumption that relative brightness derived from SAR backscatter range is proportional to the local topographic gradient. We make empirical estimations of the scaling factors and offsets required to map backscatter onto local slope and integrate along range lines to retrieve a relative height for each range pixel. We then apply a triangular mesh interpolation in radar geometry, including a smoothing factor and edge constraints, to estimate the full lava dome shape, before returning to geographic coordinates to calculate the dome volume and extrusion rates.

We exploit 32 SAR images from three sensors (TerraSAR-X, COSMO-SkyMed and Sentinel-1) acquired from January 2021 to April 2021 (Table 1) to observe the dome growth at St. Vincent (Fig. 1). The acquisition of high-resolution X-band data (i.e., CSK and TSX) was coordinated and supported by the CEOS Volcano Demonstrator, and output products were provided to UWI-SRC during the eruption to support monitoring efforts. For each sensor geometry we co-register the image to the first in the dataset and multi-look to give the pixel range dimensions shown in Supplement Table 1. Multi-looking reduces the contribution of speckle, which would otherwise introduce noise to our estimations of topography. We use the GAMMA-remote sensing software (Werner et al., 2000) to produce the backscatter images. We use a Digital Elevation Model (DEM) constructed from three 2014 Pleiades images (2 m pixels) and 2018 Copernicus DEM (30 m pixels) (Grandin and Delorme, 2021) to calibrate our estimates of relative brightness and act as reference points for the height profiles extracted from SAR cross sections.

SAR backscatter ($\sigma_0$), describes the radar pulse scattered back towards the sensor by the ground surface and is affected by local gradient, surface roughness and dielectric properties.

Surface scatterers are the main variables that determine the backscatter signal of the ground surface. Slopes facing towards the satellite will appear brighter than ones facing away (due to the effect of differences in local incidence angle), although this contribution can be corrected through radiometric terrain corrections if accurate DEMs are available (Small, 2011; Meyer et al., 2015). The roughness of a surface is dependent on instrument wavelength. Smooth surfaces will demonstrate specular reflection, but as the surface roughness increases, there will be more non-coherent or diffuse scattering, resulting in lower magnitude backscatter. Lastly, the dielectric properties of the ground surface also affect backscatter. However, in many applications moisture content is considered the main influencing factor, with higher moisture content reducing the proportion of the wave that is absorbed and increasing

the measured backscatter (Woodhouse, 2017; Flores-Anderson et al., 2019). The resulting backscatter is therefore controlled by the constructive and destructive interaction among these factors (i.e., surface roughness, local incidence angle and dielectric properties), which can lead to a complex signal that can be difficult to interpret. Further, these surface scatterers are dependent on the satellite instrument properties (e.g., polarisation, wavelength and incidence angle), making direct comparisons among sensors and tracks more challenging. Different polarisation, (e.g., VV or VH) will produce different signals over the same type of surface depending on orientation of the scattering target. For example, cross-polarised sensors are more sensitive to vegetation changes as a result of the depolarising effect from volume scattering (Woodhouse, 2017; Flores-Anderson et al., 2019). Cross-polarisation may therefore better distinguish between different degrees of vegetation coverage than co-polarisation, which can help better distinguish volcanic deposits (Solikhin et al., 2015).

La Soufrière's crater prior to the explosive eruption was dominated by grass and small shrubs, with more rock exposure found on the crater wall and the rougher 1979 dome. The area between the 1979 dome and the crater wall was relatively flat with a small crater lake in the eastern section.

We use calibrated SAR backscatter intensity images ($\sigma_0$), rather than the raw intensity images ($I_0$) in our analysis. These two quantities are related by $I_0 = k\sigma_0$, where k is a calibration constant of the SAR system. For our dataset, we assume that local incidence angle is the dominant contributor to the backscatter signal. Slopes facing towards the satellite are foreshortened, meaning there are more returns per ground unit, whereas slopes facing away are lengthened and have fewer ground returns. When the local incidence angle of slopes (i.e., angle between the normal to the slope and the satellite line-of-sight) is greater than the satellite incidence angle (i.e., angle between the satellite line-of-sight and the normal to the surface of the ellipsoid, $\theta$ in Fig. S1a), we observe radar shadows for slopes facing away and layovers for slopes facing towards the sensor. In a SAR image a radar shadow appears as a gap in data where no returns are recorded. Regions of layover do return radar signals, but as the sensor will image the top of the slope before the base, each radar unit receives returns from multiple locations at the same radar range, producing a layered image that results in ambiguity in the location of the reflectors.

### 2.1. Topographic change from SAR backscatter

#### 2.1.1. Estimation of relative brightness

SAR backscatter is dependent on sensor parameters that vary between sensors. We therefore estimate scaling factors that describe how relative brightness relates to topography. We estimate

<footnote n="1"><sup>1</sup>Sentinel-1 Descending; <sup>2</sup>TerraSAR-X Descending 085; <sup>3</sup>TerraSAR-X Descending 039</footnote>

<footnote n="4"><sup>4</sup>COSMO-SkyMed S2_17 Ascending; <sup>5</sup>COSMO-SkyMed S2_25 Ascending</footnote>
<table><caption>**Table 2**
SAR sensor parameters used to examine the St. Vincent dome growth.</caption><thead><tr><th>Sensor</th><th>Band</th><th>Incidence Angle</th><th>Polarisation</th><th>Multilook Resolution Factor</th><th>Range × Azimuth [m]</th></tr></thead><tbody><tr><td>COSMO-SkyMed</td><td>X (3.1 cm)</td><td>44.2-44.4</td><td>HH</td><td>4×4</td><td>1.5×2.8</td></tr><tr><td>TerraSAR-X</td><td>X (3.1 cm)</td><td>50.9</td><td>HH</td><td>4×4</td><td>3.6×5.1</td></tr><tr><td>Sentinel-1 (Asc)</td><td>C (5.5 cm)</td><td>43.0</td><td>VV & VH</td><td>4×1</td><td>9.3×14.0</td></tr><tr><td>Sentinel-1 (Desc)</td><td>C (5.5 cm)</td><td>31.0</td><td>VV & VH</td><td>4×1</td><td>9.3×14.0</td></tr><tr><td>ALOS-2</td><td>L (23.6 cm)</td><td>36.3</td><td>HH & HV</td><td>2×2</td><td>4.3×3.8</td></tr></tbody></table>

these factors for each image by comparing an area where topography is unchanged to the simulated backscatter ($\sigma_s$) that we calculated from the pre-eruptive DEM using a Lambertian backscatter model. We consider a 'stable' area to be one where there has been no change in topography and minimal change in scattering properties, ideally where the scaling factor is dominantly influenced by the local incidence angle. In the SAR images of St Vincent, we exclude the site surrounding the new dome and the surrounding area of vegetation die-back and alteration (approximately 27% of the St. Vincent crater in the TSX image) from estimation of scaling factors. Using an unstable area will distort the simulated SAR, which will then not correlate well with the real SAR signal and produce inaccurate scaling factors. For most datasets we used the 1979 dome as the 'stable' area used to calculate the scaling factors (Fig. 1). However, for the lower-resolution Sentinel-1 we incorporated areas outside the 1979 dome to increase the number of data points used to estimate the scaling factors. We then apply an iteratively reweighted least squares method to determine the best linear scaling factors (m, c) between the simulated backscatter, $\sigma_s$, and the real backscatter values, $\sigma_0$. These scaling factors vary among sensors due to different (1) levels of noise present, (2) calibrations of sensors, and (3) satellite geometry and parameters (i.e., wavelength). Finally, we apply the scaling factors to the whole image (including where topography has changed) to estimate a relative backscatter ($\sigma_r = m\sigma_0 + c$).

Although no radiation is scattered back to the sensor from shadows, they still have non-zero $\sigma_0$ values introduced by thermal noise and residual radiation from sidelobes. We therefore mask out the data gaps caused by radar shadows using a threshold estimated from the maximum backscatter within a known shadow in each SAR dataset (e.g., from La Soufrière crater wall, Table 2).

## 2.1.2. Retrieval of topography along range lines

Based on the assumption that the local incidence angle is the dominant scattering property of the radar backscatter signal, we can express this in terms of the gradient of topography, $h$, in radar (slant range) geometry,

$$ \sigma_r \propto \frac{dh}{d\delta} \qquad (1) $$

where $\delta$ is the range pixel spacing. A relative brightness of 1 implies flat topography, whereas values > 1 are caused by slopes that face the satellite and < 1 are caused by slopes facing away (Fig. S1).

To calculate $h$ from $\sigma_r$ we need to integrate along each range line. Following Taud and Parrot (1995) and Paquerault and Maitre (1998), we calculate relative height ($\epsilon_i$) with respect to the first point, set to 0, for the $i^{th}$ point along a line in range by;

$$ \epsilon_i = \epsilon_{i-1} + \delta(\sigma_r - 1) \qquad (2) $$

This equation is valid provided there is no layover or radar shadows (Supplementary Text). To mitigate the gaps in data caused by radar shadows, we use the width of the shadow in the SAR range line and calculate the corresponding shift in height, $\Delta_H$, for each radar shadow present in the cross-section:

$$ \Delta_H = W_{\text{range}} \cos \theta \qquad (3) $$

where $W_{\text{range}}$ is the width of the shadow in radar geometry and $\theta$ the incidence angle, (Fig. S1d). We then apply this correction $\Delta_H$ to all of the heights following that radar shadow.

Because we apply the summation independently on each range line, inconsistencies between consecutive range lines can give rise to apparent 'jumps,' typical of topographic maps derived from radarclinometry methods (Wildey, 1986; Guindon, 1990). Differences between the shape of our relative height profile and the true topography can be introduced by layover, inadequate estimate of scaling factors, speckle, or by changes in backscatter caused by differences in roughness or dielectric constant. We therefore estimate a best-fit vertical adjustment by masking out the 2021 dome and comparing the shape of our relative height profile with the pre-eruption DEM along the same line in range. Until this point, we have not made any adjustment for the influence of layovers in our estimated profiles. To accurately calculate the vertical shift, we apply a layover mask constructed from the pre-eruption DEM. We then use the DBSCAN clustering algorithm (Ester et al., 1996) to determine the different clusters between the two profiles and then calculate the linear gradient, such that a value of 1 indicates a perfect match between our estimations and DEM topography. We select the cluster with a linear gradient closest to 1 and use the corresponding y-intercept, $y$, to shift our relative height profiles, such that best estimates of the new topography, $d_{\text{topo}}$, takes the form $d_{\text{topo}(i)} = \epsilon_{i-1} + \delta(\sigma_r - 1) + y$, neglecting corrections for any shadows (Fig. 2).

## 2.1.3. Inversion for dome shape

To construct a complete dome surface from the multiple topographic profiles along range lines we divide the outline of the SAR dome into a regular triangular mesh and solve for the height of each node, assuming that the topography changes linearly across each triangle (Fig. 3). The resolution of the constructed mesh is lower than the original heights estimated from the SAR images such that multiple height estimations fall into a single triangle (Table S2). This smooths the small-scale variations originating in backscatter noise and bridges areas where there are gaps in height measurements (e.g., radar shadows). Using the topographic profiles in range, $d_{\text{topo}}$, and the DEM height values for the edge of the dome, $d_{\text{edges}}$, the system of equations for this model can be written as,

$$ [h] = (G^T W^{-1} G)^{-1} G^T W^{-1} d \qquad (4) $$

$$ G = \begin{bmatrix} G_{\text{sar}} \\ K^2 \nabla^2 \\ G_{\text{edges}} \end{bmatrix} \qquad W_{\text{diag}} = [r_1, r_2, \dots, r_n, 0] \qquad d = [d_{\text{topo}} \ 0 \ d_{\text{edges}}] \qquad (5) $$

where $h$ are the unknown height observations for the triangular vertices; $G$ is a design matrix constructed from $G_{\text{sar}}$ (number of SAR points × total number of vertices), the matrix of interpolated kernels representing the shape function (Wang and Wright, 2012) for each height observation (Fig S4); $\nabla^2$ is the Laplacian smoothing operator (number of interior vertices × total number of vertices) with an umbrella smoothing factor, $K^2$; and $G_{\text{edges}}$ is the matrix
<figure><img src="image_6.png" /><figcaption>**Fig. 2.** Topography extracted from radar backscatter range line. TerraSAR-X images (A) 13 January 2021 and (B) 31 March 2021 with cross-section lines over the 2021 dome. The growth of the 2021 dome can be seen in (C) the radar backscatter range line cross-sections with the corresponding simulated Synthetic Aperture Radar (SAR) signal produced from the Digital Elevation Model (black line) and (D) the retrieved topography for all cross-sections shown in (C). Zero values (C) and gaps in data (D) correspond to radar shadow locations.</figcaption></figure>

that contains a weighted edge constraint (number of edge vertices × total number of vertices). For the smoothing factor, $K^2$, used in the Laplacian operator, we apply a small value ($K^2 = 0.5$) to minimally smooth the structure but remove any large artifacts otherwise introduced over data gaps. Further, to avoid smoothing of the dome beyond where we know its edges to be, we include edge constraints ($G_{edges}$) that indicate where we are confident about location the dome edge ($d_{edges}$), and elevations can be obtained from the DEM in radar geometry. We reduce the contribution from range lines that were particularly distorted by radar layovers, which can cause errors to propagate along the range line (e.g., the negative and positive elevation difference over stable portions of topography, Fig. 3d-f), by applying a weighting matrix (W) to our system to reduce the influence from these particular lines.

The weighting values used for our triangular mesh interpolation are the RMSE between each range line and the corresponding pre-existing DEM line (with the new dome masked out). The diagonals of the weighting matrix ($W_{diag}$), are the respective RMSE values associated with each SAR data point ($r_1, r_2, ..., r_n$) and we assign additional zeros to the end of the diagonal that corresponds to the total number of vertices to account for the Laplacian ($K^2\nabla^2$) and edge constraint ($G_{edges}$) components in the design matrix (G). Profiles that deviate from the known topography will produce larger values, which are then accordingly down weighted in Eq. (4).

After solving for a smooth continuous mesh of dome heights in radar geometry, we translate our new elevation maps back into

geographic coordinates (Fig. 3g-i). To do this, we return our triangular mesh back into range lines and convert the coordinates of each line into Cartesian geometry (Text S1, Eq.1). We use a look-up table to geographically locate the start and end of each line.

Because we used multi-looked images for the SAR height extractions, when we translate an individual range line back into its geographic coordinates, the translation results in gaps in heights between lines. We apply a second triangular mesh in geographic coordinates to bridge these data gaps and fully construct the dome geometry without any additional weighting, smoothing, or constraints. By comparison to the original pre-eruption DEM, we can then calculate the dome heights and volumes for each SAR image (Fig. 3g-i, S6, 4). For a more complete description of the method see supplementary text and figures S1–S3.

## 2.2. Estimation of extrusion rates and uncertainties

Using our dome volumes estimates, we calculate extrusion rates between pairs of images that have the same geometry. We do this separately for each sensor to avoid introducing uncertainties associated with combining different sensor resolutions and geometries. We consider the largest sources of uncertainty in our estimates to be the systematic effects of radar layover, which introduce distortions to the topographic profiles retrieved, and which vary between satellite geometries. We calculate the volume uncertainties by estimating the percentage of the dome that is affected by layover.
<figure><img src="image_7.png" /><figcaption>**Fig. 3.** (A) Sentinel-1 (18 January 2021), (B) TerraSAR-X (15 February 2021) and (C) COSMO-SkyMed (07 April 2021) examples of the range line retrieval of topography. (D-F) show the difference between the extracted height from SAR and the pre-eruption DEM. (G-I) show the final dome shapes shown as shaded relief maps. (J-L) show the original SAR backscatter images used to estimate the dome topography. Note (A-F) are in radar geometry and (G-I) are in geographic coordinates. The green outline in E indicates the area used to calculate average difference between SAR topography and the pre-eruption DEM.</figcaption></figure>

This estimate depends on the accuracy and resolution of the pre-
eruption DEM, which we need to identify where these distortions
occur. We then use these percentage errors to calculate the up-
per and lower limits for our dome volumes. It is important to
note that while our estimations of uncertainties relate solely to the
systematic, higher-magnitude errors associated with the satellite
geometry, there are other sources of error associated with noise in
backscatter (e.g., speckle, changes to surface scatterers). We mit-
gigate the contribution of these sources of noise by multi-looking
our data, applying speckle filters, and using the triangular mesh
interpolation. However, even after these corrections, some uncer-
tainties remain, especially in areas where the properties of surfaces
scatterers vary over time (e.g., crater floor and wall). For example,
the differences between TSX-retrieved topography and the DEM for
areas where topography is unchanged vary between 10 and 24 m
(area used is outlined in green, Fig. 3e).

2.3. *SAR pixel tracking*

We also make pixel-offset tracking estimates of surface dis-
placements in both range (i.e., perpendicular to the satellite track
direction) and azimuth (i.e., parallel to the satellite track direction)
using pixel-based cross-correlation. Although Interferometric SAR
(INSAR) methods are able to measure surface displacements be-
tween acquisitions with centimetre accuracy (e.g., Delgado et al.,
2017; Moore et al., 2021), they require stable ground surface scat-
terers, otherwise the data become completely decorrelated. In con-
trast, pixel-offset tracking can measure larger surface movements
(i.e., metre-scale) without the use of SAR phase (e.g., Wang et al.,
2015; Schaefer et al., 2017) using pixel-based cross-correlation be-
tween the backscatter signals of the two images. The two SAR
images used are required to be precisely co-registered, as the ac-
curacy of the pixel-offset measurements is strongly dependent on
the image co-registration. This method only works if changes are
relatively small.
We apply pixel-offset tracking to two 1-day X-band COSMO-SkyMed (CSK) pairs to detect short-term dome growth during the early phase and immediately before the eruption transitioned from extrusive to explosive. Based on the large perpendicular baseline between our 1-day CSK image pairs (564 m and 2411 m), we apply a co-registration that makes use of a look-up table formed between the SAR images and the DEM in radar geometry based on the primary SAR image. For higher-resolution systems with large baselines, this method has the advantage that it incorporates topography-related offsets in the lookup table and applies them during co-registration. We use a small cross-correlation threshold of 0.05 to remove some of the obvious outliers from the offset results. We found that a large cross-correlation window with steps of about a third of the window size (i.e., 100×100 pixels window and 32×32 step size) maximised the spatial coverage of the dome offsets for CSK and returned appropriate offsets across the dome. Nevertheless, measurements were not possible over the centre of the dome, where surface changes were too large. Pixel-offset maps produced from 7 pairs of TerraSAR-X (TSX) images were sparser, with a less-distinct offset pattern (Fig. S11) and higher uncertainties. The difference between the TSX and CSK offset tracking quality results mainly from the longer time intervals for the TSX data (11 days), which capture a longer period of dome growth.

### 3. Results: dome growth

The gradual growth of the 2021 dome, including changes in morphology and surface texture, is observable directly from individual SAR images (Fig. S7-9, 13-16). The 2021 dome gradually changed from uniform central growth to laterally spreading around the 1979 dome. In February 2021 we observe the development of fractures at either end of the 2021 dome (Fig. S7-8, S14) and in April 2021 visible bulging in the center (Fig. S7, S14).

Our analysis has enabled us to estimate dome volumes and extrusion rates from the 32 CSK, TSX and Sentinel-1 SAR images spanning the dome-growth phase (Fig. 4). This allows us to examine patterns in extrusion rates for the dome over 3.5 months in early 2021 before the eruption transitioned from extrusive to explosive on 9 April 2021.

#### 3.1. Lava dome growth

To characterise the main dome-growth phase at La Soufrière from December 2020 to March 2021, we use a total of 29 SAR backscatter images (Fig. 4). The dome at La Soufrière showed steady growth with extrusion rates of ~0.5–5 m³s⁻¹ from December 2020 to March 2021 (Fig. 4b). Our data indicate that extrusion rate varied over short time scales with shorter timespan image pairs (e.g., CSK 1-day pair dates) recording the highest rates. However, the extrusion rate recorded over intervals > 10 days is remarkably consistent (Fig. 4). Our data are well explained with a constant rate of ~1.8 m³s⁻¹ from the first appearance of the dome until 7 April 2021, two days prior to the eruption. Our independent measurements of extrusion rate are consistent with those presented by Stinton (submitted for publication), who found extrusion rates that varied between 0.95 and 2.65 m³s⁻¹, with a long term average of 1.85 m³s⁻¹ from 12 sets of photogrammetric measurements made from locations along the rim of the summit crater, UAVs, and helicopter and fixed-wing aircraft overflights.

#### 3.2. Pre-explosive phase increase in extrusion rate

At the beginning of April 2021, the dome continued to grow laterally around the older 1979 dome without any significant change in rate. However, from the SAR data, we observed a major change

in extrusion rate in the days prior to the transition from extrusive to explosive. There was a sudden increase in extrusion rate from ~1.8 to 17.5 m³s⁻¹ some time between 21:38 UTC on 7 April and 21:38 UTC on 8 April (Fig. 4a). The rapid increase in extrusion rate is confirmed by images of rapid inflation acquired by a remote camera on the crater rim (Joseph et al., 2022). If we assume this new, higher rate was constant during the few days prior to the onset of explosive activity, then the transition to higher-rate extrusion would have occurred on 7 April 2021, just after the onset of rapid dome inflation as observed by a remote camera (Joseph et al., 2022). Our final dome volume calculated from a CSK image acquired at 21:38 UTC on 8 April, approximately 15h before the first explosions were recorded (12:41 UTC on 9 April 2021), is 18.5 million m³, slightly higher than the estimate of final dome volume of 16 million m³ made by Stinton (submitted for publication) by extrapolating a constant extrusion rate of 1.85 m³s⁻¹. Using the increased extrusion rate from the 8 April CSK image, we extrapolate the dome volume from 8 April to the time of the first explosive eruption (12:41 UTC on 9 April 2021) to estimate a final dome volume of 19.4 million m³. Because, there is no information about the extrusion rate in the hour before the eruption. Thus, we assume that the extrusion rate remained constant, which may not have been the case. The first post-explosion SAR image, acquired on 10 April at 02:03 UTC by ICEYE, showed the complete destruction of the new dome (Fig. S12).

The 1-day CSK pixel-offset tracking over the 2021 dome showed displacements within the dome itself and no significant changes elsewhere in the summit crater (Fig. 5). During the initial dome growth phase from 17-18 January 2021, we observed localised movements (maximum offset of ~9m) centred around 13°19'56"N 61°11'08"W, which we estimated to be the location of the vent. The image from 8 April 2021 showed displacements away from the same centre point extending over whole dome (maximum measurable offsets of ~15m, Fig. 5f). The order-of-magnitude increase in extrusion rate could not be accommodated by continued growth of the dome along the central axis of the coulee (e.g., Fig. S8) but can be seen bulging around the vent location on 7-8 April (Movie S14). Displacement in this area exceeds the threshold measurable with pixel-offset tracking. From our measurements, the dome surface rose by 10 m in its central section during this 24 hour period.

### 4. Discussion

#### 4.1. Transition to explosive eruption at St Vincent

The sudden increase in dome volume and extrusion rate seen at St. Vincent in the days prior to the start of the explosive phase indicates that the transition was driven by changes in the magmatic system (a 'bottom-up' process). We see no evidence in the pre-eruptive SAR imagery of dome collapse up to 15h before the eruption. However, the sudden increase in dome volume between 7 and 8 April 2021 does correlate with other evidence of rapid dome inflation in backscatter images: (1) increased fracturing, (2) surface displacements > 15 m, and (3) a bulge surrounding the vent (Fig. 5).

Independent evidence of a magmatic change that preceded the transition to explosive activity on 9 April includes (1) a change in volcano-tectonic earthquakes beneath the volcano on 5-6 April (Thompson et al., 2022), (2) the first detection of SO₂ on 8 April from TROPOMI, a coastal traverse made using UV Spectrometer and MultiGAS instruments (e.g., Joseph et al., 2022; Esse et al., in press), (3) banded tremor and incandescence, both first observed on 8 April (Joseph et al., 2022), and (4) ash venting from fractures in the rapidly inflating region on the 2021 dome.

The high volume of extruded material (16.4 million m³, red star in Fig. 4) up to 7 April and constant extrusion rate indicates the
<figure><img src="image_8.png" /><figcaption>**Fig. 4.** Time series of (A) cumulative volumes and (B) time-averaged extrusion rates of the dome at La Soufrière from December 2020 and April 2021 extracted from SAR backscatter data. The SAR backscatter data extend the photogrammetry observations showing the sudden increase in dome volume and extrusion rate in the days prior to the transition from extrusive to explosive eruption. (C) Enlarged time series of the two weeks prior to the eruption, with important changes seen in other datasets (i.e., seismicity, SO₂ emissions, visual observations superposed). SO₂ measurements are time-averaged fluxes over periods of 5–20 hours for the corresponding days (Esse et al., in press). Date format is shown as 'Day Month' starting in December 2020 until April 2021. VT = Volcano-Tectonic.</figcaption></figure>

steady removal of older material from a combined conduit and established magma reservoir, where the loss of material does not result in a discernible pressure decay. The subsequent increase in extrusion rate to $17.5 \text{ m}^3\text{s}^{-1}$ is consistent either with the addition of new magma to the system in April 2021, or, as proposed by Joseph et al. (2022)'s conceptual model, to a gas-rich magma injected earlier in the eruption (between December 2021 and April 2021) that had been slowly pushing up older material intruded in 1979. The increase in extrusion rate may have started as a gradual response to a decreasing overburden of older, degassed magma and a corresponding increase in gas exsolution rate. An alternative may be that a sudden increase in gas pressure was caused by internal processes within the reservoir, allowing gas-rich magma to break through a remaining degassed cap.

By 7 April 2021 we observed changes in the SAR backscatter that indicated a sudden increase in dome volume as the last of the old, degassed material (we estimate 3 million m³ between the 7-9 April, 15% of the total extruded volume) was erupted. Shortly after, the detection of increased SO₂ on 8 April (Esse et al., in press) indicated that SO₂ reached a critical threshold to force its way to the surface, ultimately triggering the 9 April eruption. The timing of the increase in extrusion rate estimated from SAR backscatter measurements (7 April) is consistent with the detectable increase in SO₂ emissions early on 8 April (Esse et al., in press). SO₂ emissions were reported beginning 1 February 2021 through MultiGAS measurements. However, they remained below the detection threshold

for TROPOMI and DOAS until the increase on the 8 April (Esse et al., in press).

## 4.2. Potential for application to future eruptions

Our method for extracting topographic profiles and thus extrusion rates from SAR backscatter has potential to be a powerful monitoring tool. Extrusion rate is typically one of most challenging geophysical parameters to obtain safely during an eruption, requiring either in-situ measurements near vents or clear conditions for overflights. Volume estimates made from radar satellites have potential to increase the temporal resolution of an extrusion-rate time series, with photogrammetry augmenting when conditions allow (as for the 2021 St Vincent eruption). Photogrammetry measurements shown in Fig. 4a were reported during the progression of the eruption and were based on simplifying assumptions to obtain rapid estimations of dome height, area and volume. We assume that the discrepancies between the photogrammetry and SAR dome volumes are related to the assumptions of both methods.

The usefulness of SAR-retrieved topography depends on (1) the accuracy of our underlying assumption that $\sigma_0$ is primarily controlled by topography, (2) satellite radar image resolution and geometry relative to the new volcanic topography, and (3) the quality and resolution of the pre-eruption DEM used to identify areas of shadow and layover.
<figure><img src="image_9.png" /><figcaption>**Fig. 5.** COSMO-SkyMed (CSK) pixel tracking maps. Radar backscatter images from (a) 17 January, (b) 18 January, (d) 7 April and (e) 8 April with dome outlines in blue (1st date) and red (2nd date). (c) and (f) Pixel tracking maps for (c) January and (f) April overlain on cross-correlation standard deviations. Arrows show the direction of surface displacement in range and azimuth direction, and grid represents the step window used (CSK, 32×32). Location of images is shown in Fig. 1H.</figcaption></figure>

Our assumption that local topographic gradient is the primary control on backscatter is justified only where surface roughness and dielectric proprieties are relatively steady. We can correct for these effects where variations are uniform across the structure. For the lava dome at La Soufrière, based on ground and satellite imagery, we could confidently assume that the surface roughness of the new dome was uniform across the dome and through time; however, this was not the case for other parts of La Soufrière's crater, or other domes elsewhere. For example, we observed that $σ_0$ on sections of the 1979 dome adjacent to the 2021 dome decreased over the course of the eruption, a change that we attribute to the removal of vegetation caused by the emission of gases from the growing dome (Fig. S10).

To explore the accuracy of our method, we use synthetic topographic cross-sections to test the impact of varying incidence angle and the local incidence angle (Fig S2). We found that radar layovers had a greater impact on the accuracy of the retrieved topography than radar shadows, which can be easily identified from backscatter images. The edges of regions of layover are harder to identify in backscatter, and since topographic change is unknown for real data, cannot be reliably predicted. In the synthetic models, a single section of radar layover would skew the final profile shape, but our correlations meant that the overall height estimations were similar to the original profile. When a layover interacted with a radar shadow, the estimated topography contained distortions and differ significantly from the original synthetic cross-section (Fig. S2). The distortion introduced by the layovers in our estimations of topography is under-constrained and cannot be corrected. In this dataset, layovers are most significant during the early period of dome growth, when the edges are relatively steep and not constrained by the surrounding topography. Lower incidence angles produced larger areas affected by radar layovers, and we found that an ALOS-2 ascending track and the Sentinel-1 descending track (that ultimately did not contribute to our results) were less successful in retrieving the dome heights because of smaller incidence angles, 37° and 31° respectively (Table 2). We consider

layover to be the largest source of uncertainty in our results, especially for the descending tracks as the 2021 dome grew into the layover of the crater wall, concealing the backscatter signal from the dome (Fig. S3).

### 4.3. Optimising SAR acquisitions

The different SAR datasets covering the 2021 dome growth at St Vincent allowed us to test how different sensor polarisations and resolutions affected our estimates of topography. We found that unmitigated speckle greatly increased the noise in our retrievals of topography. Higher-resolution sensors allowed us to multi-look and apply speckle corrections to the images while retaining high spatial resolution (CSK 1.5×2.8 m pixels at 4×4 looks, range by azimuth and TSX 3.6×5.1 m at 4×4 looks, range by azimuth). Although the multi-looked Sentinel-1 images (9.3×14.0 m at 4×1 looks, range by azimuth) were much lower resolution, their frequent acquisitions and multiple tracks covering the dome-growth period provided a denser temporal spacing for extrusion rate estimates and allowed a much denser time series than possible from X-band alone. We also found that the Sentinel-1 cross-polarisation images (i.e., VH) produced less accurate topographic estimates, presumably due to their greater sensitivity to volume scattering and vegetation (Patel et al., 2006).

Based on this, we found the ideal SAR imagery for estimating topography from backscatter has the following characteristics: (1) co-polarization, (2) a geometry that minimises shadow and layover, and (3) is acquired as frequently as possible. In practice, these characteristics are most likely to be achieved for any individual eruption by analysis of a constellation of different SAR instruments. The most appropriate geometry for tasking is especially hard to predict in advance of an eruption, as the shape of new topography determines to the locations of shadow and layover. Volcanic craters, domes with steep sides, and flows that fill steep valleys are all likely to be bounded by zones of shadow and layover. Al-
though this source of uncertainty can be minimised by choice of geometry, it is unlikely to be eliminated.

## 5. Conclusions

We successfully estimated topography from 32 SAR backscatter for the 2020–2021 eruption of La Soufrière, St. Vincent volcano using five different geometries from three sensors (TerraSAR-X, TSX; COSMO-SkyMed, CSK; and Sentinel-1). The topographic estimates allowed us to estimate a time series of dome volume and thus extrusion rate. The processing timeframe from raw SAR data to final retrieved topography was about done within a day, and depending on the resolution of the data and size of the target feature. We found that high-resolution, co-polarisation imagery with small areas affected by layover (i.e., larger incidence angles) resulted in the most reliable dome-volume estimates.

Extrusion rate at La Soufrière was relatively stable between December 2020 and March 2021, with an average rate of $1.8 \text{ m}^3 \text{ s}^{-1}$. By 7 April 2021, extrusion rate increased by an order of magnitude to $17.5 \text{ m}^3 \text{ s}^{-1}$, before the explosive phase on 9 April 2021, 15 hours after our last SAR measurement. The final dome volume measured from SAR data was 18.5 million $\text{m}^3$ with a maximum height of ~140 m. This increase in extrusion rate is consistent with an incursion of gas-rich magma into the conduit, which then led to the explosive phase of the eruption.

### CRediT authorship contribution statement

**E.W. Dualeh:** Conceptualization, Investigation, Methodology, Writing – original draft. **S.K. Ebmeier:** Conceptualization, Funding acquisition, Methodology, Supervision, Writing – review & editing. **T.J. Wright:** Conceptualization, Methodology, Supervision, Writing – review & editing. **M.P. Poland:** Conceptualization, Writing – review & editing. **R. Grandin:** Resources, Writing – review & editing. **A.J. Stinton:** Resources, Writing – review & editing. **M. Camejo-Harry:** Conceptualization, Resources, Writing – review & editing. **B. Esse:** Resources, Writing – review & editing. **M. Burton:** Resources.

### Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

### Data availability

Data will be made available on request.

### Acknowledgements

This work was supported by the European Space Agency (ESA Living Planet Fellowship to SKE), the NERC-BGS Center for Observation and Modeling of Earthquakes, Volcanoes and Tectonics (COMET) and the University of Leeds. SKE is funded by a NERC Independent Research Fellowship (NE/R015546/1). COSMO-SkyMed (CSK) and TerraSAR-X (TSX) imagery were provided by the CEOS Working Group on Disasters Volcano Demonstrator, and ICEYE images through SKE's ESA TPM account. CSK (provided by ASI), and TSX (provided by DLR) can be requested through the CEOS Volcano Demonstrator Program, by contacting either the demonstrator leads or the WGDisasters secretariat (https://ceos.org/ourwork/workinggroups/disasters/). Sentinel-1 data are freely available to download from https://scihub.copernicus.eu/dhus/%23/home. We would like to thank the staff at the University of the West Indies Seismic Research Centre (UWI-SRC), the Montserrat Volcano Observatory (MVO), and Soufrière Monitoring Unit (SMU) for the

work they did during and after the eruption of La Soufrière. We also thank C. Harnett for providing useful background information on dome growth literature, R. Contreras-Arratia for the valuable discussion on the seismic monitoring at St. Vincent during the 2021 eruption, and R. Robertson for confirming the rapid dome growth in April 2021. Thanks to Larry Mastin, Rick Wessels, and two anonymous reviewers for their constructive comments that helped improve the manuscript. And to Jon Major and Chris Magirl for taking time to read through the final manuscript. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government.

## Appendix A. Supplementary material

Supplementary material related to this article can be found online at https://doi.org/10.1016/j.epsl.2022.117980.

### References

Albino, F., Biggs, J., Escobar-Wolf, R., Naismith, A., Watson, I.M., Phillips, J.C., Marroquin, G.A.C., 2020. Using TanDEM-X to measure pyroclastic flow source location, thickness and volume: application to the 3rd June 2018 eruption of Fuego volcano, Guatemala. J. Volcanol. Geotherm. Res., 107063.

Angarita, M., Grapenthin, R., Plank, S., Meyer, F., Dietterich, H., 2022. Quantifying large-scale surface change using sar amplitude images: crater morphology changes during the 2019-2020 shishaldin volcano eruption. J. Geophys. Res., Solid Earth, e2022JB024344.

Arnold, D.W., Biggs, J., Anderson, K., Vallejo Vargas, S., Wadge, G., Ebmeier, S.K., Naranjo, M.F., Mothes, P., 2017. Decaying lava extrusion rate at El Reventador volcano, Ecuador, measured using high-resolution satellite radar. J. Geophys. Res., Solid Earth 122, 9966-9988. https://doi.org/10.1002/2017JB014580.

Carn, S., Watts, R., Thompson, G., Norton, G., 2004. Anatomy of a lava dome collapse: the 20 March 2000 event at Soufrière Hills Volcano, Montserrat. J. Volcanol. Geotherm. Res. 131, 241-264.

Castro, J.M., Cordonnier, B., Tuffen, H., Tobin, M.J., Puskar, L., Martin, M.C., Bechtel, H.A., 2012. The role of melt-fracture degassing in defusing explosive rhyolite eruptions at volcán Chaitén. Earth Planet. Sci. Lett. 333, 63-69.

Delgado, F., Pritchard, M.E., Ebmeier, S., González, P., Lara, L., 2017. Recent unrest (2002–2015) imaged by space geodesy at the highest risk Chilean volcanoes, Villarrica, Llaima, and Calbuco (Southern Andes). J. Volcanol. Geotherm. Res. 344, 270-288.

Esse, B., Burton, M., Mayer, C., Contreras-Arratia, R., Christopher, T., Joseph, E., Varam, M., Johnson, C., in press. SO2 emissions during the 2021 eruption of La Soufrière St. Vincent, revealed with back-trajectory analysis of TROPOMI imagery. Geological Society of London Special Publication.

Ester, M., Kriegel, H.P., Sander, J., Xu, X., et al., 1996. A density-based algorithm for discovering clusters in large spatial databases with noise. In: kdd, pp. 226-231.

Fink, J.H., Anderson, S.W., 2000. Lava domes and coulees. In: Encyclopedia of Volcanoes, pp. 307-319.

Flores-Anderson, A.I., Herndon, K.E., Thapa, R.B., Cherrington, E., 2019. The SAR handbook: comprehensive methodologies for forest monitoring and biomass estimation. SERVIR Global. https://doi.org/10.25966/nr2c-s697.

Grandin, R., Delorme, A., 2021. La Soufrière volcano (Saint Vincent) – Fusion of Pleiades (2014, 2 m) and Copernicus (2018, 30 m) digital elevation models. Dataset on Zenodo. https://doi.org/10.5281/zenodo.4668734.

Griffiths, R.W., Fink, J.H., 1997. Solidifying Bingham extrusions: a model for the growth of silicic lava domes. J. Fluid Mech. 347, 13-36.

Guindon, B., 1990. Development of a shape-from-shading technique for the extraction of topographic models from individual spaceborne sar images. IEEE Trans. Geosci. Remote Sens. 28, 654-661.

Harnett, C.E., Thomas, M.E., Calder, E.S., Ebmeier, S.K., Telford, A., Murphy, W., Neuberg, J., 2019. Presentation and analysis of a worldwide database for lava dome collapse events: the Global Archive of Dome Instabilities (GLADIS). Bull. Volcanol. 81, 1-17.

Joseph, E., Camejo-Harry, M., Christopher, T., Contreras-Arratia, R., Edwards, S., Graham, O., Johnson, M., Juman, A., Latchman, J., Lynch, L., et al., 2022. Responding to eruptive transitions during the 2020-2021 eruption of La Soufrière volcano, St. Vincent. Nat. Commun. 13, 1-15.

Metcalfe, A., Moune, S., Komorowski, J.C., Moretti, R., 2022. Bottom-up vs top-down drivers of eruption style: Petro-geochemical constraints from the holocene explosive activity at La Soufrière de Guadeloupe. J. Volcanol. Geotherm. Res. 424, 107488.

Meyer, F.J., McAlpin, D.B., Gong, W., Ajadi, O., Arko, S., Webley, P.W., Dehn, J., 2015. Integrating SAR and derived products into operational volcano monitoring and decision support systems. ISPRS J. Photogramm. Remote Sens. 100, 106-117. https://doi.org/10.1016/j.isprsjprs.2014.05.009.
Moore, C., Wright, T., Hooper, A., 2021. Rift focussing and magmatism during late-stage rifting in afar. *J. Geophys. Res., Solid Earth*, e2020JB021542.

Nakada, S., Shimizu, H., Ohta, K., 1999. Overview of the 1990–1995 eruption at Un-zen Volcano. *J. Volcanol. Geotherm. Res.* 89, 1–22.

NEMO, 2020a. Effusive Eruption at the La Soufrière Volcano alert level raised to orange. http://nemo.gov.vc/nemo/index.php/news-events/news-release/568-volcano-alert-level-raised-to-orange.

NEMO, 2020b. La Soufrière Bulletin #50 April 08, 2021 5:30 PM. http://nemo.gov.vc/nemo/index.php/news-events/news-release/638-la-soufriere-bulletin-50-april-08-2021-5-30-pm1.

NEMO, 2020c. Scientists have noted a change in seismic activity associated with the ongoing eruption of the La Soufrière Volcano. http://nemo.gov.vc/nemo/index.php/news-events/news-release/628-scientists-have-noted-a-change-in-seismic-activity-associated-with-the-ongoing-eruption-of-the-la-soufriere-volcano.

Ozawa, T., Kozono, T., 2013. Temporal variation of the Shinmoe-dake crater in the 2011 eruption revealed by spaceborne SAR observations. *Earth Planets Space* 65, 527–537.

Pallister, J.S., Schneider, D.J., Griswold, J.P., Keeler, R.H., Burton, W.C., Noyles, C., Newhall, C.G., Ratdomopurbo, A., 2013. Merapi 2010 eruption-chronology and extrusion rates monitored with satellite radar and used in eruption forecasting. *J. Volcanol. Geotherm. Res.* 261, 144–152. https://doi.org/10.1016/j.jvolgeores.2012.07.012.

Paquerault, S., Maitre, H., 1998. New method for backscatter model estimation and elevation map computation using radarclinometry. In: *SAR Image Analysis, Modeling, and Techniques*. SPIE, pp. 230–241.

Patel, P., Srivastava, H.S., Panigrahy, S., Parihar, J.S., 2006. Comparative evaluation of the sensitivity of multi-polarized multi-frequency SAR backscatter to plant density. *Int. J. Remote Sens.* 27, 293–305.

Plank, S., Walter, T.R., Martinis, S., Cesca, S., 2019. Growth and collapse of a littoral lava dome during the 2018/19 eruption of Kadovar Volcano, Papua New Guinea, analyzed by multi-sensor satellite imagery. *J. Volcanol. Geotherm. Res.* 388, 106704.

Poland, M.P., 2014. Time-averaged discharge rate of subaerial lava at Kīlauea Volcano, Hawaiʻi, measured from TanDEM-X interferometry: Implications for magma supply and storage during 2011–2013. *J. Geophys. Res., Solid Earth* 119, 5464–5481.

Robertson, R.E., 1995. An assessment of the risk from future eruptions of the Soufrière volcano of St. Vincent, West Indies. *Nat. Hazards* 11, 163–191.

Roman, D.C., Rodgers, M., Geirsson, H., LaFemina, P.C., Tenorio, V., 2016. Assessing the likelihood and magnitude of volcanic explosions based on seismic quiescence. *Earth Planet. Sci. Lett.* 450, 20–28.

Ryan, G., Loughlin, S., James, M., Jones, L., Calder, E., Christopher, T., Strutt, M., Wadge, G., 2010. Growth of the lava dome and extrusion rates at Soufrière Hills Volcano. *Montserrat, West Indies: 2005–2008*. Geophys. Res. Lett. 37.

Schaefer, L., Wang, T., Escobar-Wolf, R., Oommen, T., Lu, Z., Kim, J., Lundgren, P., Waite, G., 2017. Three-dimensional displacements of a large volcano flank move-

ment during the May 2010 eruptions at Pacaya Volcano, Guatemala. *Geophys. Res. Lett.* 44, 135–142.

Scharff, L., Hort, M., Gerst, A., 2014. The dynamics of the dome at Santiaguito volcano, Guatemala. *Geophys. J. Int.* 197, 926–942.

Small, D., 2011. Flattening gamma: radiometric terrain correction for SAR imagery. *IEEE Trans. Geosci. Remote Sens.* 49, 3081–3093. https://doi.org/10.1109/TGRS.2011.2120616.

Solikhin, A., Pinel, V., Vandemeulebrouck, J., Thouret, J.C., Hendrasto, M., 2015. Mapping the 2010 Merapi pyroclastic deposits using dual-polarization Synthetic Aperture Radar (SAR) data. *Remote Sens. Environ.* 158, 180–192. https://doi.org/10.1016/j.rse.2014.11.002.

Stinton, A., submitted for publication. Growth and evolution of the lava dome and coulée during the 2020-21 eruption of La Soufrière, St Vincent. *Montserrat Volcano Observatory*.

Taud, H., Parrot, J.F., 1995. Qualitative digital elevation model from sar data. In: *Synthetic Aperture Radar and Passive Microwave Sensing*. SPIE, pp. 226–233.

Thompson, J.O., Contreras-Arratia, R., Befus, K.S., Ramsey, M.S., 2022. Thermal and seismic precursors to the explosive eruption at La Soufrière Volcano, St. Vincent in April 2021. *Earth Planet. Sci. Lett.* 592, 117621.

Wadge, G., Cole, P., Stinton, A., Komorowski, J.C., Stewart, R., Toombs, A.C., Legendre, Y., 2011. Rapid topographic change measured by high-resolution satellite radar at Soufriere Hills Volcano, Montserrat, 2008–2010. *J. Volcanol. Geotherm. Res.* 199, 142–152. https://doi.org/10.1016/j.jvolgeores.2010.10.011.

Walter, T., Subandriyo, J., Kirbani, S., Bathke, H., Suryanto, W., Aisyah, N., Darmawan, H., Jousset, P., Luehr, B.G., Dahm, T., 2015. Volcano-tectonic control of Merapi's lava dome splitting: the November 2013 fracture observed from high resolution TerraSAR-X data. *Tectonophysics* 639, 23–33.

Wang, H., Wright, T., 2012. Satellite geodetic imaging reveals internal deformation of western Tibet. *Geophys. Res. Lett.* 39.

Wang, T., Poland, M.P., Lu, Z., 2015. Dome growth at Mount Cleveland, Aleutian Arc, quantified by time series TerraSAR-X imagery. *Geophys. Res. Lett.* 42, 10,614–10,621. https://doi.org/10.1002/2015GL066784.

Watts, R., Herd, R., Sparks, R., Young, S., Druitt, T., 2002. Growth patterns and emplacement of the andesitic lava dome at Soufriere Hills Volcano, Montserrat. *Mem. Geol. Soc. Lond.* 21, 115–152.

Werner, C., Wegmüller, U., Strozzi, T., Wiesmann, A., 2000. Gamma SAR and interferometric processing software. In: *Proceedings of the Ers-Envisat Symposium. Gothenburg, Sweden, p. 1620*. Citeseer.

Wildey, R., 1986. Radarclinometry. *Earth Moon Planets* 36, 217–247. https://doi.org/10.1007/BF00055161.

Woodhouse, I.H., 2017. Introduction to Microwave Remote Sensing. CRC Press.

Zharinov, N., Demyanchuk, Y.V., 2008. The growth of an extrusive dome on Shiveluch Volcano, Kamchatka in 1980–2007: geodetic observations and video surveys. *J. Volcanol. Seismol.* 2, 217–227.