# Very-long-period seismicity associated with the 2009–2015 reawakening of Cotopaxi Volcano, Ecuador

Indira Molina, Hiroyuki Kumagai, Mario Ruiz, Stephen Hernández, Patricia

Mothes, Gabriela Arias, Joan Andújar

► To cite this version:

Indira Molina, Hiroyuki Kumagai, Mario Ruiz, Stephen Hernández, Patricia Mothes, et al. Very-long-period seismicity associated with the 2009-2015 reawakening of Cotopaxi Volcano, Ecuador. Journal of Volcanology and Geothermal Research, 2024, 453, 10.1016/j.jvolgeores.2024.108150. insu-04676492

HAL Id: insu-04676492

https://insu.hal.science/insu-04676492v1

Submitted on 9 Sep 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
1

2

3

4 Indira Molina¹, Hiroyuki Kumagai², Mario Ruiz³, Stephen Hernández³, Patricia Mothes³, Gabriela Arias⁴, Joan Andújar⁵

5

6

7 ¹Servicio Geológico Colombiano, Bogotá, Colombia cmolinap@sgc.gov.co

8 ²Nagoya University, Nagoya, Japan kumagai@eps.nagoya-u.ac.jp

9 ³Instituto Geofísico, Escuela Politécnica Nacional, Quito, Ecuador mruiz@igepn.edu.ec, pmothes@igepn.edu.ec, shernandez@igepn.edu.ec

10 ⁴Institut Géoazur, Université Côte d'Azur arias@geoazur.unice.fr

11 ⁵ISTO, Orléans, France juan.andujar@cnrs-orleans.fr

13

14 To be submitted to: Journal of Volcanology and Geothermal Research.

15

16 ABSTRACT

17 Cotopaxi is a large, ice-capped stratovolcano located in the Ecuadorian Andes. After 72 years of 18 repose, Cotopaxi erupted on August 14, 2015. The precursory activity included long-period (LP) 19 events followed by volcano-tectonic (VT) earthquakes, very-long-period events accompanying LP 20 signals (VLP/LP events), tremor, deformation and SO₂ emissions. VLP/LP events were first 21 observed at Cotopaxi in 2002, and occurred persistently from 2009 to 2014 and during 22 the 2015 eruptions. Previous studies of the VLP/LP seismicity suggested that these events 23 originated by repetitive volume changes in a crack due to degassing of water from magma at a 24 depth of 2–3 km beneath the NE flank. Based on this interpretation, we estimated the magma 25 volumes related to individual VLP/LP events from 2009 to 2015, which were systematically
extracted from continuous seismic records of the Cotopaxi broadband seismic network. Based on the accumulated magma volume and the VLP/LP activity, our study is divided into seven periods (phases A–G), during which the magma supply rate significantly fluctuated. In phase E (June 1–July 27, 2015), before the eruptions, the magma supply rate increased. Degassing at the VLP source generated gas flows in the conduit and pre-eruptive tremor, gradually drying out a shallow hydrothermal system. In phase F (July 28–September 15, 2015), we estimated the highest magma supply rate, leading to magma fragmentation at the VLP source and eruptions accompanied by tremor. In phase G (September 16–December 29, 2015), the magma supply rate decreased, and overall eruptive activity, VLP/LP events, and tremor gradually waned. These results indicate that the VLP/LP events were likely generated by degassing from magma supplied to the VLP source through an intruded dike before and during the eruptions. The VLP/LP activity provides critical useful information about the magma supply rates that controlled eruptive and gas emission activity at Cotopaxi during this period and may help to constrain magma volumes during future reactivations.

# 1. INTRODUCTION

Cotopaxi (elevation 5876 m) is a large, ice-capped volcano, and one of the highest active stratovolcanoes in the world. Cotopaxi is located south of Quito, Ecuador's capital city of Ecuador, and 40 km north of Latacunga (Fig. 1a). Its eruptive activity has been characterized by periodical andesitic and rhyolitic eruptions (Hall and Mothes, 1995). Over the past 500 years, Cotopaxi has experienced multiple cycles of eruptions in 1532–1534, 1742–1744, 1766–1768, 1853–1855, 1877–1880, 1903–1906, and 1942 (Coordinated Universal Time, UTC, is used throughout this paper). During this period, there were 70 eruptions with Volcanic Explosivity Indices (VEIs) of
2-4, characterized by ash falls, as well as pyroclastic, debris and lava flows (Andrade et al., 2005, Mothes et al., 1998; Hidalgo et al., 2018). Barberi et al. (1995) estimated the average recurrence interval of eruptive episodes of Cotopaxi to be 117 years based on historical and stratigraphic records over the last 2000 years. According to these historical records, smaller and less intense eruptions, lasting months to years, preceded large eruptions (Hidalgo et al., 2018). Identifying pre-eruptive signals before any unrest is crucial because of the population centers lying near the volcano, and especially near lahar channels.

The Cotopaxi seismic network was deployed in 1989 and has been maintained and operated by the Instituto Geofísico at the Escuela Politécnica Nacional (IG-EPN). It has recorded various seismic signals, such as long-period (LP) events, tremor, volcano-tectonic (VT) earthquakes, very-long-period (VLP) events accompanying LP signals (i.e., VLP/LP events), explosions, and ice quakes in the summit glaciers (e.g., Ruiz et al., 1998; Metaxian et al., 2003; Molina et al., 2008; Kumagai et al., 2010; Bernard et al., 2016; Hidalgo et al., 2018). According to the classification scheme of Chouet (1996a), LP and VLP denote seismic events with frequency ranges of 0.5-5 Hz (0.2-2 s) and 0.01-0.5 Hz (2-100 s), respectively. LP events display multiple oscillation cycles, while VLP events often show one cycle, except for some unusual events reported at Mt. Erebus and Kilauea (Waite, 2015). LP events are interpreted as acoustic vibrations of fluid-filled resonators in magmatic and hydrothermal systems (e.g., Chouet, 1996b). VLP events are linked to mass transport processes in magmatic and hydrothermal systems (e.g., Chouet and Matoza, 2013, and references cited therein) related to magma intrusion and eruptive activity (e.g., Arciniega-Ceballos et al., 1999; Molina et al., 2008; Maeda and Takeo, 2011; Park et al., 2020). VLP events exhibit high-frequency signals at event onsets when they are associated with explosions (see Fig. 1d of Ripepe
et al., 2021). However, VLP events are not always associated with explosions or surface activity (e.g., Hill et al., 2002; Marchetti and Ripepe, 2005). VLP/LP events at Cotopaxi show both VLP and LP signatures (Molina et al., 2008).

Alternatively, tremor usually displays emergent onsets and dominant frequencies in the range of 0.066–10 Hz (0.1–15 s) (e.g., McNutt, 1992; Yukutake et al., 2017). The source processes of tremor are interpreted as (i) flow-induced oscillations of conduit walls, (ii) successive cracking of rocks due to injection of fluids, and (iii) hydrothermal interaction with magmatic gas and ash emissions (Li et al., 2022 and references cited therein).

According to Ruiz et al. (1998), the seismic activity at Cotopaxi between 1989 and 1997 was caused by the interaction between magmatic heat and the shallow hydrothermal system. During renewed seismic activity in 2001–2003, VLP/LP events were observed at Cotopaxi for the first time between June and September 2002 (Molina et al., 2008). Later, VLP/LP events from 2002 to 2015 were reported by Molina et al. (2008), Kumagai et al. (2010), Arias (2016), and Anzieta et al. (2019). VLP/LP sources were located below the NE flank of the volcano at a depth of 2–3 km (Molina et al., 2008; Kumagai et al., 2010; Arias, 2016; Anzieta et al., 2019). During the seismic crisis from 2001 to 2003, Molina et al. (2008) reported on an increasing number of LP events in January 2001, which were followed by a VT swarm in November 2001 and subsequent deformation beneath the NE flank in or earlier than November 2001, prior to the appearance of fumarolic activity in early January 2002 (Global Volcanism Program, 2002) (Fig. 2a). These changes led to the interpretation of a magmatic intrusion beneath the NE flank. Continuing into late June 2002, impulsive VLP signals accompanied by LP signals (i.e., VLP/LP events), characterized by a transition from non-
harmonic to harmonic oscillations were observed. VLP signals displayed a broad spectral peak of around 0.5 Hz (2 s), and accompanied LP signals characterized by non-harmonic oscillations in 0.5-5 Hz (0.2-2 s) (see Figs. 3 and 6 of Molina et al., 2008). Molina et al. (2008) performed waveform inversion of VLP signals and suggested that these VLP and LP signals were produced by volume changes in the crack due to degassing of water at the top of an intruded magma and the resonance of a crack filled with ash-laden gases released from the magma, respectively. The crack was located beneath the NE flank at a depth of 2-3 km. Molina et al. (2008) also performed particle motion analyses for the initiation onset (onset) portions of the VLPs events and found that these events were located beneath the NE flank. Using waveform inversion analysis Kumagai et al. (2010) estimated an inclined tensile crack mechanism for a VLP/LP event in 2009. This VLP/LP event displayed impulsive VLP signals at its onset in a frequency band of 0.08 and 0.2 Hz (5-12.5 s), followed by LP signals consisting of non-harmonic oscillations with dominant peaks in 0.5-10 Hz (0.1-2 s) (see Figs. 4 and 5 of Kumagai et al., 2010). The source of this VLP event was located at 2 km depth below the NE flank. In 2009, Cotopaxi's activity was accompanied by an increase in LP events and sporadic SO₂ emissions of around 500 t/d (Bourquin et al., 2009). Small increases in LP events also occurred in 2010 and 2011, and an active fumarole became more evident on the western flank.

In early April 2015, a swarm of LP events was registered (Hidalgo et al., 2018; Gaunt et al., 2016), peaking to 8-12 VLP/LP events per day in early May 2015 (Anzieta et al., 2019). On 4 June 2015, LP activity transitioned to tremor characterized by broadband and non-harmonic oscillations, which occurred almost continuously over the following two months (Hidalgo et al., 2018). This tremor preceded a VT swarm, which then coincided with four phreatomagmatic eruptions with
VEI=2 on August 14 and 15, 2015 (Hidalgo et al., 2018; Bernard et al., 2016). These eruptions were associated with numerous VLP/LP events (Anzieta et al., 2019) and tremor, which notably increased on 21 August 2015 (Hidalgo et al., 2018). Ash emissions characterized the eruptive activity that continued until the end of November 2015 (Bernard et al., 2016; Gaunt et al., 2016; Hidalgo et al., 2018).

From 2001–2002, ground deformation was detected under the NE (Molina et al., 2008) and SW flanks of the volcano, indicating magma migration from the SW to NE flanks along NNE–SSW trending faults (Hickey et al., 2015). Morales Rivera et al. (2017) proposed that the observed deformation between April and August 2015 resulted from an inclined sheet intrusion located a few kilometers SW of the summit. The volume was estimated to be 6.8 × 10⁶ m³, extending from a depth of 12.1 km to 5.5 km below the summit. Based on petrological experiments, Martel et al. (2018) proposed that the 2015 eruptions were fed by a continuous flux of magma with an average water content of 3.5 wt%, originating from depths of 16 km to 7 km below the summit.

Aerial surveys conducted by IG-EPN revealed that the crater lake that existed before 14 August 2015 had disappeared, suggesting that the magmatic heat dried out water in the shallower system (Gaunt et al., 2016). Hidalgo et al. (2018) estimated a magma volume of 65.3 × 10⁶ m³ to supply the emitted SO₂ gases during the 2015 eruptions. According to Bernard et al. (2016), the volume of erupted materials was estimated to be 8.6 × 10⁵ m³ based on the analysis of tephra falls. They found that the amplitude variations of seismic tremors were correlated with the variations of ash fallout mass emitted during the 2015 eruptions. Based on satellite images, Arnold et al. (2018) estimated an erupted volume between 5 × 10⁵ and 1.5 × 10⁶ m³.
In this paper, we perform an analysis of VLP/LP events associated with the reawakening of Cotopaxi Volcano during 2009–2015 (Fig. 2) to estimate the temporal evolution of the VLP/LP activity and its relation to eruptive activity. Waveform inversion, as performed by Kumagai et al. (2010), is the direct approach to quantifying VLP/LP sources. However, the VLP/LP event analyzed by Kumagai et al. (2010) was the largest event during the studied period, and our application of waveform inversion to other smaller events was unsuccessful due to low signal-to-noise ratios in these signals. Therefore, we took another approach to estimate the source volumes of the individual VLP/LP events by using their observed amplitudes. We assume that all the VLP/LP events occurred at the same source location and with the same mechanism as those of the largest VLP/LP event estimated by Kumagai et al. (2010). This assumption is justified by the existence of a unique family of VLP/LP events in 2009–2015 as indicated by Anzieta et al. (2019) and Arias (2016). Anzieta et al. (2019) employed unsupervised machine learning techniques in a two-stage process to analyze seismic events using waveform data at station BREF (Fig. 1). Their findings identified a unique family of VLP/LP events registered throughout 2015. This same study also showed that VLP events disappeared during periods of heightened tremor activity and resurged just before eruptions on August 14, 2015. Additionally, Arias (2016) applied a waveform cross-correlation method to study waveform data at BTAM (Fig. 1). The result was an average cross-correlation value of 0.8 between the largest VLP/LP event reported by Kumagai et al. (2010) and VLP/LP events in 2009–2015, indicating that the same source repetitively triggered these VLP/LP events.

We then estimated the magma volume related to the source of each VLP/LP event caused by H₂O degassing following the interpretation of Molina et al. (2008). Accumulating these estimates, we derived the temporal evolution of magma supply rate in 2009–2015. We found that a sharp increase
in magma supply rate occurred at the beginning of the 2015 eruptive activity. We also note that the
amplitudes and the daily number of VLP events were correlated with tremor amplitudes and SO₂
emissions during the eruptive activity. Based on these results and other observations, we discuss
the source processes of VLP/LP events and tremor and their implications for magmatic activity
before and during the 2015 eruptions at Cotopaxi Volcano.

170

171 **2. DATA**

172 A network of four permanent seismic stations equipped with 1 Hz seismometers was deployed in
173 1989 at Cotopaxi. The network was reinforced in 2001 by three additional permanent stations with
174 1 Hz seismometers and one temporal station featuring a seismometer with bandwidths of 0.2-40
175 Hz (Molina et al., 2008). In 2006, a joint project between the Japan International Cooperation
176 Agency (JICA) and IG-EPN (Kumagai et al., 2007; 2010) installed a seismic network (Fig. 1b) in
177 which each seismic station is comprised of a triaxial broadband seismometer (Güralp CMG-4T:
178 0.02-50 Hz). Seismic waveform data were digitized by a 24-bit data logger (Geotech Smart24D)
179 with a sampling frequency of 50 Hz at each station and then transmitted by a digital telemetry
180 system using 2.4 and 5 GHz wireless LAN to the central office of IG-EPN in Quito through
181 repeaters.

182

183 Waveform data acquired from the Cotopaxi broadband seismic network (Fig. 1b) was continuously
184 sampled and triggered by a detection system using filtered waveform amplitudes (Kumagai et al.,
185 2010). In this system, the maximum filtered waveform amplitudes using the latest one-minute files
186 were estimated. If the maximum amplitude at each station during the latest one-minute interval
exceeded a pre-established threshold, and if the number of triggered stations exceeded a pre-set number, the waveform data were saved as an event for the latest three-minute-long files. Two frequency bands for this detection trigger system were used: 0.3–2 Hz (LP band) and 0.02–0.1 Hz (VLP band). These bands were empirically determined by performing tests for efficient triggering of LP and VLP events at Cotopaxi by Kumagai et al. (2010). In this study, we observe that VLP/LP events were triggered in the LP band because their amplitudes were larger in this band (Fig. A1).

The trigger level was set according to the noise level at each station. We used waveform data from four broadband stations (BREF, BVC2, BTAM, and BMOR; Fig. 1b) to identify a long-term trend in the VLP/LP activity. Due to poor data quality, waveform data at BNAS (Fig. 1b) were excluded.

BREF provided the highest quality data due to its proximity to the VLP/LP events, although this station was occasionally affected by seismic signals from the descent of secondary lahars, mainly occurring on the western flank after the eruption onset in August 2015. To examine the relationship between VLP/LP events and tremor, we used the tremor data from Bernard et al. (2016) during the 2015 eruptive activity, which was derived from continuous waveform data at BREF. They manually filtered out tectonic earthquakes to calculate the median and root mean square (RMS) amplitudes of tremor signals. This process was performed using 10-minute sliding windows (without overlap) in the 0.5–5 Hz band. We also generated similar tremor data in 2009–2015 using 1-day sliding windows (without overlap) in the 0.5–5 Hz band. Using the solution for the moment source time function of the largest VLP/LP event recorded at Cotopaxi on January 14, 2009 at 18:54 (Kumagai et al., 2010), we calculated the magma supply rate as explained below (see Section 3.1).

Since 2008, SO₂ has been measured at Cotopaxi by two permanent scanning Differential Optical Absorption Spectroscopy (DOAS) stations located 2.5 km north (REF) and 5 km west (NASA) of the crater (Fig. 1b). Two additional stations were installed in 2015 (Hidalgo et al., 2018). The
processed data provided the mass emission of SO₂ in tons per month in 2009–2013 and the detailed evolution in tons per day in 2015 (Hidalgo et al., 2018). Since 2010, IG-EPN has conducted geodetic measurements using two surface and two borehole tiltmeters designed to detect small ground deformation changes. From 2010 to 2015, the data acquisition procedures remained unchanged, in which digital data were recorded every 5 minutes and simultaneously transmitted to IG-EPN. Tiltmeters on the NE sector responded very well to magma ascent, and VC1 showed notable changes (Fig. 1b). At VC1, a biaxial surface-mount tiltmeter is set on a cement base within an old, well-settled massive lava. The summit glaciers are located 3.5 km upslope and have no direct contact with VC1, which ensures that this tiltmeter has no relation to the changes in the glacier thickness or melting. The distance from the VLP source to the tiltmeter is 6.1 km. Another tiltmeter above the Refuge (close to REF DOAS station in Fig. 1b) also recorded moderate signals.

# 3. METHODS

## 3.1 VLP/LP EVENT SELECTION AND AMPLITUDE DETERMINATION

We used the catalog of LP and VLP events triggered by the automatic system between 2009 and 2015. From this catalog, we extracted the VLP/LP events by visually inspecting waveform plots with 10-minute durations. According to Peterson (1993) and Webb (1998), the microseism falls within the range of approximately 0.04–0.2 Hz (5–25 s). This microseism band overlaps with the regional and teleseismic earthquakes and VLP events in a frequency band of approximately 0.04–0.2 Hz. To distinguish a microseism from VLP/LP events, we followed the procedure proposed by Dawson and Chouet (2014), which uses the power spectral density (PSD) for noise over 1 hour before a VLP/LP event and the PSD for the same 1-hour length of both noise and VLP/LP signals.
Based on the spectral comparison, we determined that an amplitude larger than approximately 10⁻⁶ m/s for VLP/LP events is a suitable threshold for selecting events in our catalog. To avoid incorporating false VLP/LP events into our database, we eliminated events classified as regional or teleseismic earthquakes by the Ecuadorian National Seismic Network of IG-EPN and USGS catalogs. We then performed the following operations to select VLP/LP events with good signal-to-noise ratios. After instrumental correction, vertical (UD) displacement seismograms were band-passed between 0.08 and 1.0 Hz using a two-pole Butterworth filter to cut a time window of 30 s for each VLP/LP event. This window consisted of 10 s before and 20 s after the maximum peak-to-peak amplitude. We then selected VLP/LP events whose maximum amplitudes were three times larger than the average amplitudes over the first 5 s (noise level). The time window was also applied to the EW and NS components for each event at the same station. We manually deleted the events with poor waveform features (i.e., pulses and noise). We determined the amplitude ratios of the individual horizontal components to the UD component in a frequency band of 0.08–0.5 Hz and estimated peak-to-peak amplitudes in the individual components. We chose the frequency band of 0.08–0.5 Hz because it produces larger amplitude signals within the range of our instrument response (Fig. A1) and removes high-frequency signals affected by structural heterogeneities.

249

250

## 3.2 VLP/LP EVENTS WAVEFORM FEATURES

The vertical velocity seismograms and amplitude spectra observed at BREF for VLP/LP and LP events are displayed in Fig. 3. The VLP/LP seismogram shows impulsive VLP signals with spectral peaks in the frequency band of 0.2–1.6 Hz (0.63–5 s) with dominant peaks between 0.6 and 0.9 Hz (1.11–1.66 s), which are accompanied by LP signals consisting of non-harmonic oscillations and broad spectral peaks in the frequency band of 1.6–10 Hz (0.1–0.63 s) with dominant peaks between

233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255
2 and 4 Hz (0.25–0.5 s). VLP/LP events with similar spectral features occurred in 2002 (Molina et al., 2008), 2009 (Kumagai et al., 2010), 2010–2014 (Arias, 2016), and 2015 (Anzieta et al., 2019).

Vertical velocity seismograms band-passed between 0.08 and 0.5 Hz for selected VLP events observed at BREF between 2009 and 2015, show similar VLP waveforms throughout this period (Fig. 4). This indicates that a single family of VLP events likely occurred at Cotopaxi Volcano, as suggested by Anzieta et al. (2019) and Arias (2016). Higher frequency signals in the 1–5 Hz range (as shown in Fig. 3b) appear in the VLP signals during the increased activity in 2009 and 2015.

Particle motion analysis by Arias (2016) and Anzieta et al. (2019) and waveform inversion analysis by Molina et al. (2008) and Kumagai et al. (2010) revealed that the VLP sources are located at a depth of 2–3 km beneath the NNE flank.

## 3.3 SOURCE VOLUME REPRESENTATION

A seismic source may be represented by the moment tensor, and the n-th component of the displacement field is related to the moment tensor through the following equation (e.g., Aki and Richards, 1980):

$$u_n(\mathbf{x}, t) = M_{pq}(t) * G_{np,q}(\mathbf{x}, t), \qquad (1)$$

where t is time, x is space, Mpq is the pq-th component of the moment tensor, Gnp,q is the spatial derivative of Green’s function, and * denotes the convolution operation. In the frequency domain, Eq. (1) is equivalently written as $u_n(x, \omega) = M_{pq}(\omega)G_{np,q}(x, \omega)$, where $\omega$ is the angular frequency. Based on the waveform inversion results of Kumagai et al. (2010), the most appropriate seismic source to account for the largest VLP event at Cotopaxi is a tensile crack, for which the moment tensor (M) can be written as follows (e.g., Chouet, 1996b):
$$279 \quad \mathbf{M} = s(t) \begin{pmatrix} \lambda/\mu + 2\sin^2\theta \cos^2\theta & 2\sin^2\theta \sin\phi \cos\phi & 2\sin\theta \cos\theta \cos\phi \\ 2\sin^2\theta \sin\phi \cos\phi & \lambda/\mu + 2\sin^2\theta \cos^2\phi & 2\sin\theta \cos\theta \sin\phi \\ 2\sin\theta \cos\theta \cos\phi & 2\sin\theta \cos\theta \cos\phi & \lambda/\mu + \cos^2\theta \end{pmatrix}, \quad (2)$$

where

$$281 \quad s(t) = \mu \Delta V(t). \qquad (3)$$

Here, $\theta$ and $\phi$ are the angles defining the normal direction of a crack, $\lambda$ is Lamé's first elastic parameter, $\mu$ is Lamé's second elastic parameter or shear modulus, $s(t)$ is the source time function, and $\Delta V$ is the volume change of the crack (e.g., Chouet, 1996a). For the same crack angles ($\theta$ and $\phi$) and elastic parameters ($\lambda$ and $\mu$), $u_n$ in a frequency band is thus proportional to $\Delta V$ in this band, and the displacement amplitude ratios of EW/UD and NS/UD only depend on values of the crack angles.

We assume that the volume change in the crack is due to degassing of water from intruded magma (see Section 4.2). By applying the ideal gas law, we can calculate the mass of steam ($M_g$) associated with a VLP/LP event as follows:

$$292 \quad M_g = \frac{\Delta V P_l}{\bar{R} T_m}, \qquad (4)$$

where $\bar{R}$ is the specific gas constant for steam, $T_m$ is the magma temperature, and $P_l$ is the pressure of the magma. $M_g$ is related to magma volume ($V_m$) through the following equation (e.g., Taguchi et al., 2018):

$$296 \quad V_m = \frac{M_g (1-c_o)}{\rho_m (c_o - \sqrt{K_H P_l})}, \qquad (5)$$

where $c_o$ is the initial water content of the magma, $\rho_m$ is the density of the liquid magma, and $K_H$ is Henry's law constant. We assume the lithostatic pressure for $P_l$.
With reference to Section 4.1, the displacement amplitude ratios (EW/UD and NS/UD) of our analyzed VLP/LP events at the individual stations were almost constant, indicating that the crack angles did not change. Since the UD component is indisputable, we estimated the maximum amplitude of the vertical displacement seismogram ($A_{ref}^j$) in a frequency band of 0.08–0.5 Hz for the largest VLP/LP event analyzed by Kumagai et al. (2010) at the j-th station. Using the following approximation, we obtained the magma volume ($V_{ij}^m$) for the i-th VLP/LP event at the j-th station:

$$V_{ij}^m = \frac{A_i^j V_m}{A_{ref}^j}, \qquad (6)$$

where $A_i^j$ is the maximum amplitude of the vertical band-passed displacement seismogram of the i-th VLP/LP event at the j-th station, and $V_m$ is the magma volume estimated from the maximum amplitude of the source time function $s(t)$ for the largest VLP/LP event analyzed by Kumagai et al. (2010). We then summed the $V_{ij}^m$ values to derive the accumulated magma volume and magma supply rate.

# 4. RESULTS

## 4.1 VLP/LP EVENT AMPLITUDE EVOLUTION

We observed maximum displacement amplitudes in the UD, EW, and NS components of the individual VLP/LP events within the frequency band of 0.08–0.5 Hz (2–12.5 s) at BREF, BTAM, BVC2, and BMOR. The temporal evolution of the UD amplitudes (Fig. 5) shows increases around January 14 and November 2, 2009, from early April to November 2011, early February to late July 2012, early June 2015, late August, and mid-October 2015. However, the EW/UD and NS/UD ratios remained stable over the studied period (Fig. 5). This suggests that the crack angles did not
change during the VLP activity and supports our assumption that the displacement of either
component (UD, EW, or NS) is proportional to ΔV, i.e., the volume change. We selected the most
apparent UD value to calculate the magma volume estimates.

No data were recorded in 2009 at BTAM. Temporal changes in the UD amplitudes at BMOR in 2009 were smoother than those at BREF. The seismometer at BMOR is buried in unconsolidated materials on a ridge at a remote location, about 6 km from the summit (Fig. 1b). In contrast, BREF, located closest to the summit, is deployed on massive lava. Amplitude variability among these stations may originate from such local site effects.

## 4.2 MAGMA VOLUME ESTIMATES

According to the VLP source estimates of Molina et al. (2008), Kumagai et al. (2010), Arias (2016), and Anzieta et al. (2019), we assume magma emplacement at the VLP source occurred at a depth of $h = 2500$ m (the average depth estimated by these studies) beneath the NNE flank. To determine the magma volume at the VLP source, we employed Eqs. (3)–(6) as follows: $\mu = 10.6$ GPa for the 2500 m depth, which was derived from an S-wave velocity ($\beta$) of 2000 m/s and a rock density ($\rho_s$) of 2650 kg/m$^3$ (Molina et al., 2008). Incorporating the source time function $s(t)$ estimates presented by Kumagai et al. (2010) for the largest VLP/LP event, we obtained a maximum value of $6.32 \times 10^{12}$ Nm. We then estimated $\Delta V$ for this event to be 596.2 m$^3$. According to Martel et al. (2018), recently erupted magma at Cotopaxi (primarily mafic andesite) was stored at depths ranging from approximately 7 km to $\ge 16$ km (with pressures from ~200 to $\ge 400$ MPa) and contained an average water content of 3.5 wt%. Since water is the most abundant volatile in magma (Wallace, 2005), we assume that $\Delta V$ was caused by the degassing of water. Using a lithostatic
pressure ($P_t$) of 65 MPa, $\bar{R} = 462$ J/(K kg), $T_m = 1223$ K, and $\Delta V = 596.2$ m³, from Eq. (4), we calculated $M_g = 6.85 \times 10^4$ kg for the largest VLP/LP event. Using Eq. (5), we estimated $V_m$ for this VLP/LP event to be $9.54 \times 10^3$ m³, assuming a water content ($C_o$) of 3.5 wt% in andesitic magma with a density ($\rho_m$) of 2500 kg/m³ and Henry's constant ($K_H$) of $1.6 \times 10^{-11}$ Pa⁻¹ (Nishimura, 2004). Using the $V_m$ estimate, we computed $V_{ij}^m$ values from Eq. (6). The $A_{ref}^j$ values are listed in relation to the individual stations (Table 1).

The accumulated amplitude and accumulated magma volume at BREF compared to tremor activity are shown in Fig. 6a. Fig. A2 compares the accumulated magma volumes estimated at BVC2, BTAM, and BMOR with those at BREF. Similar relative temporal changes are observed across all stations, though BMOR shows the largest value (Fig. A2). As mentioned in the previous subsection, these differences are likely attributable to site effects. Subsequently, we focused on estimates at BREF. Based on the temporal variation in the accumulated magma volume (Fig. 6a) and the individual VLP/LP amplitudes (Fig. 7a), we divided the VLP/LP activity into seven phases (A through G) (Table 2). Specifically, we observe gradual increases in accumulated magma volume between January 14 and September 18, 2009 (phase A) and between June 1 and July 27, 2015 (phase E), followed by sharp increases between September 19 and November 25, 2009 (phase B) and between July 28 and September 15, 2015 (phase F), respectively (Figs. 6 and 7). The gradual increase in phase E was associated with an increased number of LP events (Fig. 2c), elevated SO₂ emissions (Fig. 2e), and heightened tremor activity (Fig. 6a). In phase F, a significant number of VLP/LP events occurred in August 2015 (about 200 events). Two phreatomagmatic explosions associated with VLP/LP events took place on August 14–15, 2015, after which tremor reached its peak during the studied period (Figs. 6a and 7a). Subsequently, the number of VLP/LP events
367 decreased in September 2015, and tremor associated with ash venting also decreased (Fig. 7a).

368 From September 7 until the end of November 2015 in phase G, VLP/LP and LP events, tremor, and

369 SO₂ emissions gradually decreased (Figs. 2, 6a, and 7).

370
371 From mid-February to mid-March 2010 and in April 2011, there were minor increases in the
372 accumulated magma volume. Therefore, we defined phase C as lasting from November 26, 2009,

373 to April 2011 (Fig. 6). From mid-September 2014, the accumulated magma volume again slightly
374 increased. Therefore, we defined phase D from April 2011 to May 2015 (Fig. 6a). The magma
375 supply rates varied significantly among the phases, with phase G's production being approximately
376 one order of magnitude smaller than those in other active periods (phases B, E, and F) (Table 2).

377
378
4.3 VLP/LP EVENTS, TREMOR, ASH FALLOUT MASS, SO₂ EMISSION,
379
AND DEFORMATION

380 Tremor amplitudes were correlated with the accumulated amplitudes of VLP/LP events during the
381 pre-eruptive, eruptive, and post-eruptive periods in 2015 (see phases E, F, and G in Fig. 6a, and
382 Table 2), where the VLP/LP accumulated amplitudes (i.e., magma supply rates) increased with
383 expanding tremor amplitudes. Before these periods, from 2009 to 2015, no evident relationship
384 between VLP/LP and tremor amplitudes was observed. However, it is noteworthy that tremor
385 amplitudes showed correlations with SO₂ emissions during this period (Figs. 2e and 6a). Bernard
386 et al. (2016) divided the eruptive activity in 2015 into 4 phases based on ash fallout mass, visual
387 observations, and tremor signal amplitudes. Phase 1 (August 14–15, 2015) was marked by the
388 opening of the conduit, which was associated with a major explosion on August 14, 2015 at 09:02,
389 followed by smaller explosions on the same day, generating plumes reaching 9.3 km above the
crater. Phases 2 (August 15–October 2, 2015), 3 (October 2–November 4, 2015), and 4 (November 4–30, 2015) were characterized by intense ash emissions. The magnitudes of these emissions increased between phases 1 to 2 and later decreased over phases 2 to 4. These phases of Bernard et al (2016) coincide with our phases F and G (Figs. 6b and 7). Pre-eruptive LP activity occurred from early April to June 4, 2015 (Figs. 2a and 2c), followed by a registry of pre-eruptive tremor activity from June 4 to August 13, 2015, during which phase E occurred (Fig. 6b).

Based on tremor activity, we defined peaks 1–4 in phases E and F (Fig. 7a). We also identified peaks E and F based on the VLP/LP activity in phases E and F, respectively (Fig. 7a). We compared the tremor amplitudes with the individual VLP/LP amplitudes and ash fallout mass (Fig. 7a), as well as the daily number of VLP/LP events and SO₂ emissions (Fig. 7b). The number of daily events and the amplitudes of VLP/LP events are highly correlated. Additionally, 45 days after the onset of pre-eruptive tremor, we observe a close correlation between VLP/LP events, tremor, ash fallout mass, and SO₂ emissions (Fig. 6b). This correlation substantiates the inception of phase E, which reaches its peak on June 2, 2015 (Fig. 7a). After peak E, the VLP/LP amplitudes became smaller during phase E until the end of July 2015. A fortnight later, phase F started with large VLP/LP events, reaching peak F, which coincided with the first phreatomagmatic explosion on August 14, 2015 at 09:02.

In summary, the eruptive activity in 2015 (Fig. 7a), particularly during peak F, exhibited the highest daily number of VLP/LP events (Fig. 7b) and the most energetic VLP/LP activity when the conduit system was opened. Subsequently, the VLP/LP amplitudes decreased during phase F, except for an event that occurred on September 7, 2015, at 16:56 (Fig. 7a). In phase G, the VLP/LP activity remained low, although a large event occurred on December 17 at 02:01 (peak G). It is notable that
this event took place after the eruptive activity had already ceased.

Prior to the eruptions, pre-eruptive tremor on June 2, 2015 (peak 1) and VLP/LP activity on June 6, 2015 (peak E) occurred before an SO₂ emission of ~4500 t/d on June 14, 2015 (Fig. 7b). On July 14, 2015, there was another SO₂ emission peak of ~5000 t/d (Fig. 7b). This peak was not linked to either VLP/LP or tremor activity, and SO₂ emission then decreased gradually (Fig. 7b). The high VLP/LP activity on August 13, 2015, occurred before the opening of the conduit (phase 1) with the first explosion on August 14, 2015 at 9:02 (peak F), increased tremor activity, and an SO₂ emission peak of ~7200 t/d on August 15, 2015 (Fig. 7b). Following the first eruption on August 14, 2015, VLP/LP events, tremor activity, and SO₂ emission sharply decreased, although tremor activity peaked again on August 28, 2015 (peak 2) (Fig. 7a). SO₂ emission gradually increased to reach the maximum peak of ~16000 t/d on October 12, 2015, which almost coincided with peak 3 of tremor activity (Fig. 7b). In phases 3 and 4, when the conduit system was open, tremor activity and SO₂ emission showed similar tendencies (Fig. 7b). Therefore, our results suggest that there are correlations among VLP/LP events, tremor activity, and SO₂ emission. Before the eruptions, large VLP/LP events tended to occur before tremor activity peaks, which were followed by SO₂ emission peaks. During the eruptions, particularly in phases 3 and 4, tremor activity and SO₂ emissions occurred simultaneously with a small number of VLP/LP events.

We performed a time-frequency analysis using continuous seismic data at BREF. We divided one day's data into 1080 windows (80 seconds each), over which the PSD was calculated in a frequency range of 0.2–12 Hz. Fig. 6b shows the normalized PSD as a function of time (spectrogram) between January and December 2015. The relatively large spectrum amplitudes in a frequency band of 2.9–7.8 Hz in early April and early June 2015 correspond to pre-eruptive tremor, which contains several
frequency peaks around 2.9, 3.5, 4.6, 5.4, 6.5, and 7.5 Hz. During phase E, we observe low-frequency peaks around 0.5, 0.6, 0.8, 0.9, and 1 Hz of VLP/LP events, as well as the higher frequency peaks of tremor signals. A peak frequency of 2.9 Hz gradually decreased to 2.1 Hz in phase E, but other peaks remained stable (Fig. 6b).

Tiltmeter records (tangential and radial components to the summit crater) at VC1 (Fig. 1b) are compared with the daily number of VLP/LP events in Fig. 8. Cyclic variations with amplitudes of approximately 100 µrad were observed each year from late 2010 to early 2015, attributed to seasonal changes in groundwater and temperature. Significant deformation exceeding the seasonal changes, particularly in the radial component, began to occur from April 10–25, 2015, coinciding with the onset of pre-eruptive LP activity and continued until the peak of VLP/LP activity on June 19, 2015. Then, the deformation gradually decreased in phase F, during which the intensive eruptions occurred. The largest peak of the VLP/LP activity was reached on August 25, 2015 (peak F, Fig. 7a). No clear, convincing deformation signals associated with the VLP/LP activity were detected by seven GPS stations on the volcano in 2015 (Mothes et al., 2017). The depth of the deformation source has not been estimated. Molina et al. (2008) also observed deformation of the NE flank prior to the appearance of the VLP/LP events from 2001 to 2003.

# 5. DISCUSSION

We interpret the VLP/LP and tremor activity at Cotopaxi from 2009 to 2015 as follows. The VLP/LP events were generated by degassing of water from magma supplied at the crack source at a depth of 2–3 km below the NNE flank of the volcano (Molina et al., 2008; Kumagai et al., 2010;
Arias, 2016; Anzieta et al., 2019). Our results indicate that the magma supply rates fluctuated significantly from 2009 to 2015 (Fig. 6b and Table 2). The magma rate during the eruption activity in 2015 (phase F) was higher than during the activity periods in 2009 (phases A and B), when no eruptions occurred (Table 2). The continuous supply of magma through an intruded dike and degassing at a depth of 2–3 km caused volume changes in a crack, leading to its resonances, which were observed as VLP/LP events (Fig. 9a). When the magma supply rate was less than $10^{-2}$ m³/s, degassing and resultant bubble growth were not intensive, and no eruptions occurred (Fig. 9a). Subsequently, the crack was filled with a magma-gas mixture, and gas was gradually released through the partly closed conduit. The degassed dense magma descended and generated magma convection in the dike (e.g., Kazahaya et al., 1994) (Fig. 9a). At this stage, magma intrusion feeding the crack may have caused pressure fluctuations, evidenced by slight increases in LP and VLP/LP events (Figs. 2c and 2d), and resonances in a shallow conduit observed as tremor from 2009 to June 4, 2015, and as pre-eruptive tremor from June 4 to August 13, 2015 (Fig. 9a). Although the shallower magma system may have been partly closed (Fig. 9a), some gases could escape through the conduit and generate periodic variations, seen as tremor fluctuations accompanied by SO₂ emissions from 2009 to 2014 (Figs. 2e and 6a). There was no evident relationship between VLP/LP and tremor activities between 2009 and May 31, 2015 (Fig. 6a), during which no large magma supply rates were estimated (Table 2). This suggests that the crack was not well connected to the conduit during this period. Before, the eruptions degassing became relatively intense. At this stage, large VLP/LP events preceded the increasing tremor activity (see the beginning of phase E in Fig. 7a), and subsequently, the amplitudes of tremor and VLP/LP events were well correlated (see phase E in Figs. 6a and 7a). The NE flank deformation over seasonal variations was detected during this period (Fig. 8). As the magma supply rate increased (phase F), bubble growth became more intensive and resulted in fragmentation (Fig. 9b). At this stage, significant VLP/LP events again
preceded the increase in tremor activity, with tremor and VLP/LP amplitudes showing a strong correlation (see phases 1 and F in Figs. 6a and 7a). The crack was filled with ash and gas, which ascended in the conduit and produced ash emissions. Such flows in the conduit generated tremor associated with eruptions. The fragmentation condition was reached in phase F (Fig. 7), during which the maximum magma supply rate ($1.12 \times 10^{-2} \text{ m}^3\text{/s}$) was estimated (Table 2). The intruded dike was located in the source region of VT events (Ruiz et al., 1998) and may have been connected through NNE-SSW fault systems to a magma reservoir at depths of 4–5 km beneath the SW flank (Hickey et al., 2015). This reservoir was probably fed by a deeper reservoir located at depths from 7 to ≥ 16 km (Martel et al., 2018). Distal GPS data suggest that the deepest reservoir is located at a depth of ~24 km below the SE flank (Mothes et al., 2017).

In 2015, although VLP/LP events were smaller individually, they occurred more frequently, indicating a higher magma supply rate. A VT swarm preceded the 2015 eruption, unlike in 2009. Deformation and seismicity in 2001–2003 were due to a low influx of magma accumulation that stalled beneath the surface. The 2009 activity also indicated a low influx of magma. Magma intrusions can occur without significant surface deformation due to fracture systems accommodating the magma (Puglisi et al., 2008). At Cotopaxi, NNE-SSW trending faults may accommodate magma from deeper parts without causing a VT swarm. Low magma supply rates can result in fewer VT earthquakes (McNutt, 2005; Poland et al., 2008). VT earthquakes, caused by brittle responses to stress from magmatic fluid movement, are linked to dike intrusions (e.g., Toda et al., 2022). LP seismicity represents magma pathways from the dike-VT earthquake depth to the surface. The crack-generating VLP/LP events at Cotopaxi are proposed to be located above this dike intrusion.
Increased pre-eruptive LP activity occurred between April 1 and the end of May 2015 (Figs. 2c and 6b). During this period, significant deformation occurred on the NE flank (Fig. 8). Due to limited data, estimating the geometry and exact location of the deformation source is not possible; however, we interpreted the deformation as being caused by an inflationary source beneath the NE flank, which occurred in or before April 2015. This coincided with a significant increase in the number of LP events (Figs. 2a and 2c) and a rise in VLP/LP events at the end of May 2015 (Figs. 2a and 2d). This behavior of deformation preceding VLP/LP events was also observed by Molina et al. (2008) during the seismic crisis at Cotopaxi in 2001–2003. In May 2015, IG-EPN reported intense fumarolic activity at the summit crater. Large VLP/LP events on June 1–2, 2015 (peak E) were followed by pre-eruptive tremor on June 4, 2015 (Figs. 6b and 7). As the magma supply rate gradually increased in phase E (Fig. 6a), pre-eruptive tremor reached its maximum amplitude on June 6, 2015 (peak 1 in Fig. 7), and the deformation peaked on June 8, 2015, 7 days after the start of phase E (Fig. 8). As the tremor amplitude decreased, its dominant peak frequency changed from 2.9 to 2.1 Hz from June 6 to August 13, 2015 (Fig. 6b), which may be caused by variations in gas flows in the conduit (e.g., Yukutake et al., 2017). Hot gases passing through the partially closed conduit gradually heated the shallow hydrothermal system (Fig. 9a); eventually, the water reached its boiling point, causing the hydrothermal system to dry out gradually. Then, SO₂ gas emission became active, reaching 4500 t/d on June 14, 2015 after peaks E and 1 (Fig. 7b). The peak of VLP/LP events occurred 8 days before the peak of SO₂ emissions. The appearance of sulfide deposits and a green-colored lake in the summit crater during the last week of July 2015 (Gaunt et al., 2016) resulted from the interaction between magmatic gases and hydrothermal water provided by the melting of the summit glacier. Gaunt et al. (2016) found that hydrothermal minerals such as opal and alunite were emitted during the first phreatomagmatic explosions. These explosions were
caused by bubbly magma ascending through the conduit and coming into contact with
hydrothermal water.

The large VLP/LP event on August 14, 2015 at 9:02 (peak F in Fig. 7a) and another event at 9:07 on the same day partially cleared the conduit (Bernard et al., 2016). Plumes of explosions associated with these events reached heights of ~7.8 km above the crater, and these eruptions were followed by a peak in tremor amplitude (~900 a.u.) 6 minutes later (Fig. 7a). SO₂ emissions reached ~7000 t/d on August 15. The peak of VLP/LP events occurred one day before the peak of SO₂ emissions and only six minutes before the onset of tremor related to the resonance of the conduit.
Two other explosions (August 15, 2015 at 15:45 and 21:02) finally cleared the entire conduit, producing plumes that reached 9.3 km above the crater. These explosions on August 14 and 15, 2015 were classified as phreatic due to the scarcity of juvenile material (Gaunt et al., 2016) and were accompanied by a small VT swarm. This opening of the system, evidenced by the large explosions, likely improved connectivity between the crack and the conduit (Fig. 9b). On these days, the peaks of VLP/LP and tremor amplitudes and SO₂ emissions occurred almost simultaneously. Hidalgo et al. (2018) propose that this was due to the hydrothermal system being superheated and decompressed. At this stage, we believe the conduit had a direct link to the crater (Fig. 9b).

After the explosions on August 14–15, 2015, the relative amount of juvenile materials in ash samples increased (Gaunt et al., 2016), suggesting faster magma ascent. This matches our estimate in phase F, which shows the most abundant magma supply rate (Table 2). During this period, tremor activity reached its maximum intensity on August 28, 2015 (peak 2 in Fig. 7a). According to Gaunt
et al. (2016), after the explosions on August 14–15, 2015, a thermal aureole developed around the magmatic column, indicating drying processes in the hydrothermal system. At this stage, we observed the stabilization of deformation, which may reflect the release of pressure within the magma intrusion, possibly due to the opening of a conduit that allowed the escape of gases and magma (see phase F in Fig. 8). As eruptions continued, juvenile grains transitioned to a greater abundance of micro-crystalline materials due to slower magma ascent (Gaunt et al., 2016). This also agrees with our estimate of a low magma supply rate in phase G (Table 2). The continued stabilization of deformation (see phase G in Fig. 8) may be due to a decrease in the influx of magma from greater depths. The variable crystal textures depending on magma supply rates are consistent with those shown by Wright et al. (2012) for ash particles produced by eruptions at Tungurahua volcano, Ecuador.

We estimated the total magma volume during the eruptive activity of phases F and G to be $5.4 \times 10^4 \text{ m}^3$ (Table 2). It is one order of magnitude smaller than magma volumes of $\sim 8.6 \times 10^5 \text{ m}^3$ and $3 \times 10^5 \text{ m}^3$ determined by field measurements of ash fallouts (Bernard et al., 2016) and deformation data analysis (Arnold et al., 2018), respectively, during the eruptive activity. As shown in Fig. 7, VLP/LP events were not always accompanied by tremor, and thus there existed eruptions without VLP/LP events, which resulted in our smaller magma volume estimate (see phase 4 and the end of phase G in Fig. 7a and Table 2). According to Morales Rivera et al. (2017) and Mothes et al. (2016), magma volumes of $6.8 \times 10^6 \text{ m}^3$ and $(42 \pm 26) \times 10^6 \text{ m}^3$ at depths of 12–13 km and 24 km, respectively, caused the deformation of the volcanic edifice before the 2015 eruptions. Hidalgo et al. (2018) concluded that only 2–6 wt% of the magma supply to the reservoir was emitted during the eruptive activity, based on analyses of SO₂ emissions and comparisons with the results from
deformation, geological, and petrological studies (e.g., Mothes et al., 2016; Bernard et al., 2016; Gaunt et al., 2016; Morales Rivera et al., 2017; Arnold et al., 2018). These results indicate that a large amount of magma remains at depths of ~7–24 km.

Although the magma volumes estimated from the VLP/LP activity reflect a partial amount of magma supplied to the conduit, our results indicate that the temporal variation in the accumulated magma volume is useful to infer magma supply processes. During the eruptive activity, tremor is more directly linked to magma volume, as indicated by Bernard et al. (2016). In addition, we found that the increase in tremor activity was triggered by large VLP/LP events during pre-eruptive periods (Fig. 7a). When the magma supply rate was large, the seismic sequence occurred in the following order at Cotopaxi: LP events, VT swarms, VLP/LP events, and tremor. The time difference among VLP/LP events, tremor, and SO₂ emission became shorter when the magma supply rate was larger and the system was more open. When the magma supply rate was lower and the magma system was more closed, tremor was related to the passive degassing of SO₂ emissions (Figs. 2e and 6a).

# 6. CONCLUSIONS

We systematically analyzed broadband seismic signals observed at Cotopaxi volcano and estimated magma volumes related to the individual VLP/LP events from their observed amplitudes in a frequency band of 0.08–0.5 Hz (2–12.5 s) in 2009–2015. The VLP/LP events were hypothesized to be generated by the degassing of water from magma supplied to the conduit. Based on the accumulated magma volume and the VLP/LP activity, our studied period was divided into 7 phases
599 (A through G), during which the magma supply rate largely fluctuated. We found correlations 600 among VLP/LP events, tremor, deformation, SO₂ emissions, and other observations and interpreted 601 the renewed activity as follows. In phases A–D, between January 14, 2009 and May 31, 2015, 602 during which the magma supply rate was low, degassing and resultant bubble growth at the VLP 603 source were not intensive, and no eruptions occurred. In phase E (June 1–July 27, 2015) before the 604 eruptions, the magma supply rate increased, degassing was relatively intensive, and the resultant 605 gas flows in the conduit generated pre-eruptive tremor. Hot gases surging through the conduit 606 gradually heated and dried out the shallow hydrothermal system. In phase F (July 28–September 607 15, 2015), the magma supply rate sharply increased, resulting in magma fragmentation at the VLP 608 source. This produced the large VLP/LP event on 14 August 2015 and associated explosions with 609 tremor, which opened the conduit. The increased VLP/LP and tremor activity and SO₂ emissions 610 occurred almost simultaneously. In phase G (September 16–December 29, 2015), the magma 611 supply rate decreased and eruptive activity gradually waned. We estimated from the VLP/LP events 612 the total magma volume during the eruptive activity of phases F and G to be 5.4 × 10⁴ m³, which 613 is one order of magnitude smaller than magma volumes estimated from field surveys of ash fallout 614 mass and deformation analysis. 615

616 The VLP/LP events that occurred before and during the eruptions provide a comprehensive view 617 of the renewed magmatic activity. VLP/LP events, tremor, and SO₂ emissions at Cotopaxi were 618 related through non-destructive and continuous processes. The time differences among VLP/LP 619 events, tremor, and SO₂ emissions may depend on the connectivity within the dike, crack, and 620 conduit system. The opening of Cotopaxi's conduit was gradually achieved before eruptions, 621 resulting in shorter time intervals. 622
A small magma input reduces the supply of SO₂ to the summit. This condition corresponds to the stage of long-term passive degassing, characterized by periodic fluctuations of SO₂ and tremor between January 2009 and June 4, 2015. Our results suggest that future volcanic unrest at Cotopaxi may begin when the magma supply rate exceeds 10⁻³ to 10⁻² m³/s. VLP/LP activity is a useful indicator to estimate the magma supply rate and monitor future magmatic activity initiated by a large amount of magma remaining in Cotopaxi's reservoir.

We thank the members of the Instituto Geofísico de la Escuela Politécnica Nacional (IG-EPN) for their continuous efforts to maintain the Cotopaxi monitoring network. We thank Benjamin Bernard for providing tremor and ash fallout mass data. We would like to thank Monica Segovia, Silvana Hidalgo, and Caroline Martel for their helpful comments. Comments from two anonymous reviewers helped to improve the manuscript. This work is a contribution under the cooperation between IG-EPN and SGC (reference #018) from 2014.

Aki, K., Richards, P.G., 1980. Quantitative Seismology. Freeman, New York.

Arias, G., 2016. Estudio de las señales sísmicas de muy largo periodo del Volcán Cotopaxi. Bachelor Thesis, National Polytechnic School, Quito, Ecuador. URL: https://bibdigital.epn.edu.ec/handle/15000/13739.

Anzieta, J., Ortiz, H., Arias, G., Ruiz, M., 2019. Finding possible precursors for the 2015 Cotopaxi Volcano eruption using unsupervised machine learning techniques. Int. J. Geophys. 2019, 6526898, 8 pages. https://doi.org/10.1155/2019/6526898.

Arnold, D.W.D., Biggs, J., Wadge, G., Mothes, P., 2018. Using satellite radar amplitude imaging for monitoring syneruptive changes in surface morphology at an ice-capped stratovolcano. Remote Sens. Environ. 209, 480-488. https://doi.org/10.1016/j.rse.2018.02.040.

Bernard, B., Battaglia, J., Proaño, A., Hidalgo, S., Vasconez, F., Hernandez, S., Ruiz, M., 2016.
Relationship between volcanic ash fallouts and seismic tremor: quantitative assessment of the 650 651 2015 eruptive period at Cotopaxi volcano, Ecuador. Bull. Volcanol. 78, 80. 652 https://doi.org/10.1007/s00445-016-1077-5.

Andrade, D., Hall, M., Mothes, P., Troncoso, L., Eissen, J.-P., Samaniego, P., Egred, J., Ramón, P., 653 Rivero, D., Yepes, H., 2005. Los Peligros Volcánicos Asociados con el Cotopaxi. In: Los 654 peligros volcánicos en Ecuador, No. 3. Corporación Editora Nacional. 655

Arciniega-Ceballos, A., Chouet, B., Dawson, P., 1999. Very long-period signals associated with 656 vulcanian explosions at Popocatepetl Volcano. Geophys. Res. Lett. 26, 3013-3016. 657 658 https://doi.org/10.1029/1999GL010756.

Bourquin, J., Hidalgo, S., Arellano, S., Troncoso, L., Galle, B., Arrais, S., Vásconez, F., 2009. First 659 observations of intermittent, non-eruptive gas emissions of Cotopaxi volcano (Ecuador) 660 during a period of heightened seismicity. Eos Trans. AGU 90 (52), Fall Meet. Suppl., Abstract 661 662 V23D-2140.

Chouet, B., 1992. A seismic model for the source of long-period events and harmonic tremor. In: 663 Gasparini, P., Scarpa, R., Aki, K. (Eds.), Volcanic Seismology. Springer, New York, pp. 133- 664 156. https://doi.org/10.1007/978-3-642-77008-1_5.

Chouet, B.A., 1996a. New methods and future trends in seismological volcano monitoring. In: 666 Scarpa, R., Tilling, R.I. (Eds.), Monitoring and Mitigation of Volcano Hazards. Springer, 667 Berlin, pp. 23-97. https://doi.org/10.1007/978-3-642-61313-5_2.

Chouet, B., 1996b. Long-period volcano seismicity: its source and use in eruption forecasting. 669 Nature 380, 309-316. https://doi.org/10.1038/380309a0.

Chouet, B.A., Matoza, R.S., 2013. A multi-decadal view of seismic methods for detecting 671 precursors of magma movement and eruption. J. Volcanol. Geotherm. Res. 252, 108-175. 672 673 https://doi.org/10.1016/j.jvolgeores.2012.11.013.

Dawson, P., Chouet, B., 2014. Characterization of very-long-period seismicity accompanying 674 summit activity at Kilauea Volcano, Hawai'i: 2007-2013. J. Volcanol. Geotherm. Res. 278- 675 279, 59-85. https://doi.org/10.1016/j.jvolgeores.2014.04.010.

Gaunt, H.E., Bernard, B., Hidalgo, S., Proaño, A., Wright, H., Mothes, P., Criollo, E., Kueppers, 677 U., 2016. Juvenile magma recognition and eruptive dynamics inferred from the analysis of 678 ash time series: The 2015 reawakening of Cotopaxi volcano. J. Volcanol. Geotherm. Res. 328, 679 134-146. https://doi.org/10.1016/j.jvolgeores.2016.10.013.
681 Global Volcanism Program, 2002. Report on Cotopaxi (Ecuador) (Wunderman, R., ed.). Bull.
682 Global Volcanism Network 27:12. Smithsonian Institution.
683 https://doi.org/10.5479/si.GVP.BGVN200212-352050.

684 Hall, M., Mothes, P., 1995. Bi-modal nature of the eruptive history of Cotopaxi Volcano, Ecuador.
685 IUGG XXI General Assembly, Boulder, Colorado. Abstract p. A452.

686 Hickey, J., Gottsmann, J.H., Mothes, P., 2015. Estimating volcanic deformation source parameters
687 with a finite element inversion: The 2001-2002 unrest at Cotopaxi volcano, Ecuador. J.
688 Geophys. Res. Solid Earth 120(3), 1473-1486. https://doi.org/10.1002/2014JB011731.

689 Hidalgo, S., Battaglia, J., Arellano, S., Sierra, D., Bernard, B., Parra, R., Kelly, P., Dinger, F.,
690 Barrington, C., Samaniego, P., 2018. Evolution of the 2015 Cotopaxi eruption revealed by
691 combined geochemical and seismic observations. Geochem. Geophys. Geosyst. 19, 2087-
692 2108. https://doi.org/10.1029/2018GC007514.

693 Hill, D.P., Dawson, P.B., Johnston, M.J.S., Pitt, A.M., Biasi, G.P., Smith, K., 2002. Very-long-
694 period volcanic earthquakes beneath Mammoth Mountain, California. Geophys. Res. Lett.
695 29(10), 1370. https://doi.org/10.1029/2002GL014833.

696 Kazahaya, K., Shinohara, H., Saito, G., 1994. Excessive degassing of Izu-Oshima volcano: magma
697 convection in a conduit. Bull. Volcanol. 56, 207-216. https://doi.org/10.1007/BF00302080.

698 Kumagai, H., Nakano, M., Maeda, T., Yepes, H., Palacios, P., Ruiz, M., Arrais, S., Vaca, M.,
699 Molina, I., Yamashina, T., 2007. Enhancing volcano-monitoring capabilities in Ecuador. Eos
700 Trans. AGU 88, 245-246. https://doi.org/10.1029/2007EO230001.

701 Kumagai, H., Vaca, M., Palacios, P., Nakano, M., Maeda, T., Yepes, H., Molina, I., Arrais, S.,
702 Ruiz, M., Yamashina, T., 2010. Broadband seismic monitoring of active volcanoes using
703 deterministic and stochastic approaches. J. Geophys. Res. 115, B08303.
704 https://doi.org/10.1029/2009JB006889.

705 Li, K.L., Bean, C.J., Bell, A.F., Ruiz, M.C., Hernández, S., Grannell, J., 2022. Seismic tremor
706 reveals slow fracture propagation prior to the 2018 eruption at Sierra Negra volcano, Galápagos.
707 Earth Planet. Sci. Lett. https://doi.org/10.1016/j.epsl.2022.117533.

708 Maeda, Y., Takeo, M., 2011. Very-long-period pulses at Asama volcano, central Japan, inferred
709 from dense seismic observations. Geophys. J. Int. 185, 265-282. https://doi.org/10.1111/j.1365-
710 246X.2011.04938.x.

711 Marchetti, E., Ripepe, M., 2005. Stability of the seismic source during effusive and explosive
712 activity at Stromboli Volcano. Geophys. Res. Lett. 32, L03307.
713 https://doi.org/10.1029/2004GL021406.

714 Martel, C., Andújar, J., Mothes, P., Scaillet, B., Pichavant, M., Molina, I., 2018. Storage conditions
715 of the mafic and silicic magmas at Cotopaxi, Ecuador. J. Volcanol. Geotherm. Res. 354, 74-
716 86. https://doi.org/10.1016/j.jvolgeores.2018.02.006.

717 McNutt, S.R., 1992. Volcanic tremor. In: Encyclopedia of Earth System Science. Academic Press,
718 San Diego, California, pp. 417-425.

719 McNutt, S.R., 2005. Volcanic Seismology. Annu. Rev. Earth Planet. Sci. 33, 461-491.
720 https://doi.org/10.1146/annurev.earth.33.092203.122459.

721 Metaxian, J.-P., Araujo, S., Mora, M., Lesage, P., 2003. Seismicity related to the glacier of
722 Cotopaxi Volcano, Ecuador. Geophys. Res. Lett. 30, 1483.
723 https://doi.org/10.1029/2002GL016773.

724 Molina, I., Kumagai, H., García-Aristizábal, A., Nakano, M., Mothes, P., 2008. Source process of
725 very-long-period events accompanying long-period signals at Cotopaxi Volcano, Ecuador. J.
726 Volcanol. Geotherm. Res. 176, 119-133. https://doi.org/10.1016/j.jvolgeores.2007.07.019.

727 Morales Rivera, A.M., Amelung, F., Mothes, P., Hong, S.-H., Nocquet, J.-M., Jarrin, P., 2017.
728 Ground deformation before the 2015 eruptions of Cotopaxi volcano detected by InSAR.
729 Geophys. Res. Lett. 44, 6607-6615. https://doi.org/10.1002/2017GL073720.

730 Mothes, P.A., Hall, M.L., Janda, R.J., 1998. The enormous Chillos Valley lahar: an ash-flow
731 generated debris flow from Cotopaxi Volcano, Ecuador. Bull. Volcanol. 59, 233-244.

732 Mothes, P., Nocquet, J.-M., Morales, A., Jarrin, P., Gaunt, H.E., Yepez, M., Viracucha, G., 2016.
733 Geodetic signature of 2015-16 unrest at Cotopaxi-Ecuador: Modeling of GPS data, a deep
734 magma source, synchronous seismic swarms and petrologic constraints. Cities on Volcanoes-9,
735 At Puerto Varas, Chile, Volume: S3.

736 Mothes, P.A., Ruiz, M., Viracucha, G., Ramon, P., Hernandez, S., Hidalgo, S., Bernard, B., Gaunt,
737 E., Jarrin, P., Yepez, M., Espin, P., 2017. Geophysical footprints of Cotopaxi's unrest and minor
738 eruptions in 2015: An opportunity to test scientific and community preparedness. In: Gottsmann,
739 J., Neuberg, J., Scheu, B. (Eds.), Volcanic Unrest. Advances in Volcanology. Springer, Cham.
740 https://doi.org/10.1007/11157_2017_10.

741 Nishimura, T., 2004. Pressure recovery in magma due to bubble growth. Geophys. Res. Lett. 31,
742 L12613. https://doi.org/10.1029/2004GL019810.
743 Park, I., Jolly, A., Lokmer, I., Kennedy, B., 2020. Classification of long-term very long period (VLP) volcanic earthquakes at Whakaari/White Island volcano, New Zealand. Earth Planets 744 Space 72, 92. https://doi.org/10.1186/s40623-020-01224-z.

746 Peterson, J., 1993. Observation and modeling of seismic background noise. USGS Technical Report 93-322.

748 Poland, M.P., Miklius, A., Sutton, A.J., 2008. A seismic shift in the monitoring of volcanic unrest. Nat. Geosci. 1(7), 539-540. https://doi.org/10.1038/ngeo279.

750 Puglisi, G., Branca, S., Mattia, M., Aloisi, M., Bonaccorso, A., Boschi, E., 2008. Ground deformation patterns at Mt. Etna before, during, and after the July-August 2001 eruption. J. Geophys. Res. Solid Earth 113(B7). https://doi.org/10.1029/2006JB004791.

753 Ripepe, M., Delle Donne, D., Legrand, D., Valadade, S., Lacanna, G., 2021. Magma pressure discharge induces very long period seismicity. Sci. Rep. 11, 20065. https://doi.org/10.1038/s41598-021-99513-4.

756 Ruiz, M., Guillier, B., Chatelain, J.L., Yepes, H., Hall, M., Ramon, P., 1998. Possible causes for the seismic activity observed in Cotopaxi Volcano, Ecuador. Geophys. Res. Lett. 25, 2305-2308. https://doi.org/10.1029/98GL01514.

759 Taguchi, K., Kumagai, H., Maeda, Y., Torres, R., 2018. Source properties and triggering processes of long-period events beneath volcanoes inferred from an analytical formula for crack resonance frequencies. J. Geophys. Res. Solid Earth 123, 7550-7565. https://doi.org/10.1029/2018JB015866.

763 Toda, S., Stein, R., Sagiya, T., 2002. Evidence from the AD 2000 Izu islands earthquake swarm that stressing rate governs seismicity. Nature 419, 58-61. https://doi.org/10.1038/nature00997.

765 Waite, G.P., 2015. Very-long-period seismicity at active volcanoes: source mechanisms. In: Beer, M., Kougioumtzoglou, I.A., Patelli, E., Au, S.K. (Eds.), Encyclopedia of Earthquake Engineering. Springer, Berlin, Heidelberg. https://doi.org/10.1007/978-3-642-35344-4_46.

768 Wallace, P.J., 2005. Volatiles in subduction zone magmas: concentrations and fluxes based on melt inclusion and volcanic gas data. J. Volcanol. Geotherm. Res. 140(1-3), 217-240. https://doi.org/10.1016/j.jvolgeores.2004.07.023.

771 Webb, S.C., 1998. Broadband seismology and noise under the ocean. Rev. Geophys. 36(1), 105-142. https://doi.org/10.1029/97RG02287.

773 Woods, J., Donaldson, C., White, R.S., Caudron, C., Brandsdóttir, B., Hudson, T.S., Ágústsdóttir,
T., 2018. Long-period seismicity reveals magma pathways above a laterally propagating dyke during the 2014-15 Bárðarbunga rifting event, Iceland. Earth Planet. Sci. Lett. 490, 216-229. https://doi.org/10.1016/j.epsl.2018.03.020.

Wright, H.M.N., Cashman, K.V., Mothes, P.A., Hall, M.L., Ruiz, A.G., Le Pennec, J.-L., 2012. Estimating rates of decompression from textures of erupted ash particles produced by 1999-2006 eruptions of Tungurahua volcano, Ecuador. Geology 40(7), 619-622. https://doi.org/10.1130/G32948.1.

Yukutake, Y., Honda, R., Harada, M., Doke, R., Sakai, S., Ueno, T., Morita, Y., 2017. Analyzing the continuous volcanic tremors detected during the 2015 phreatic eruption of the Hakone volcano. Earth Planets Space 69, 164. https://doi.org/10.1186/s40623-017-0751-y.
784

TABLES

785

<caption>**Table 1.** Maximum displacement amplitudes in a frequency band of 0.08–0.5 Hz of the largest VLP/LP event that occurred on January 14, 2009, at the individual stations ($A_{ref}^j$) and total accumulated amplitude and magma volume estimated from VLP/LP events from 2009 to 2015.</caption>

786

787

788

<table><thead><tr><th>Station</th><th>A<sub>ref</sub><sup>j</sup><br>(m)</th><th>Total accumulated<br>amplitude (m)</th><th>Total accumulated magma<br>volume (m<sup>3</sup>)</th></tr></thead><tbody><tr><td>BREF</td><td>5.96 × 10<sup>-5</sup></td><td>1.37 × 10<sup>-3</sup></td><td>2.19 × 10<sup>5</sup></td></tr><tr><td>BVC2</td><td>2.27 × 10<sup>-5</sup></td><td>3.38 × 10<sup>-4</sup></td><td>1.42 × 10<sup>5</sup></td></tr><tr><td>BMOR</td><td>1.33 × 10<sup>-5</sup></td><td>5.49 × 10<sup>-4</sup></td><td>3.94 × 10<sup>5</sup></td></tr><tr><td>BTAM</td><td>2.49 × 10<sup>-5</sup></td><td>3.80 × 10<sup>-4</sup></td><td>1.46 × 10<sup>5</sup></td></tr><tr><td>Average</td><td></td><td>6.59 × 10<sup>-4</sup></td><td>2.25 × 10<sup>5</sup></td></tr></tbody></table>

789

790

791

<caption>**Table 2.** Magma volumes and magma supply rates in the individual phases (A–G) estimated from the amplitudes in a frequency band of 0.08–0.5 Hz at BREF.</caption>

792

793

<table><thead><tr><th>Phase</th><th>Start date</th><th>End date</th><th>Magma volume (m<sup>3</sup>)</th><th>Magma supply rate (m<sup>3</sup>/s)</th></tr></thead><tbody><tr><td>A</td><td>January 14, 2009</td><td>September 18, 2009</td><td>3.25 × 10<sup>4</sup></td><td>1.52 × 10<sup>-3</sup></td></tr><tr><td>B</td><td>September 19, 2009</td><td>November 25, 2009</td><td>4.73 × 10<sup>4</sup></td><td>8.17 × 10<sup>-3</sup></td></tr><tr><td>C</td><td>November 26, 2009</td><td>April 26, 2011</td><td>1.64 × 10<sup>4</sup></td><td>3.67 × 10<sup>-4</sup></td></tr><tr><td>D</td><td>April 27, 2015</td><td>May 31, 2015</td><td>3.63 × 10<sup>4</sup></td><td>2.81 × 10<sup>-4</sup></td></tr><tr><td>E</td><td>June 1, 2015</td><td>July 27, 2015</td><td>2.38 × 10<sup>4</sup></td><td>4.92 × 10<sup>-3</sup></td></tr><tr><td>F</td><td>July 28, 2015</td><td>September 15, 2015</td><td>4.72 × 10<sup>4</sup></td><td>1.12 × 10<sup>-2</sup></td></tr><tr><td>G</td><td>September 16, 2015</td><td>December 29, 2015</td><td>6.32 × 10<sup>3</sup></td><td>7.32 × 10<sup>-4</sup></td></tr></tbody></table>

794

795
796

797

-0.62

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

Colombia

Quito

Cotopaxi

Latacunga

Ecuador

Peru

(a)

-80

-78

-76

798

-78.5

-78.4

799

800

**Fig. 1.** (a) Location of Cotopaxi Volcano in Ecuador. (b) Topographic contours of Cotopaxi Volcano, every 200 m, showing the locations of the seismic, DOAS, and tiltmeter stations.
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Temporal variations in VLP/LP, VT, and LP activity and SO₂ gas emission. (a) Monthly numbers of LP, VT, and VLP/LP events from 2002 to 2015. (b, c, d, e) Monthly numbers of VT, LP, and VLP events and SO₂ emission (t/month), respectively, from 2009 to 2015. Vertical arrows in (e) display periods when SO₂ increased.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Vertical waveforms and their spectra of (a, b) VLP/LP and (c, d) LP events at BREF. An LP band of 1.6–10 Hz (0.1–0.63 s) is indicated by gray in each spectrum, and a VLP band of 0.2–1.6 Hz (0.63–5 s) is indicated by white in (b).</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Fig. 4. Vertical velocity seismograms of VLP signals in VLP/LP events band-passed in a frequency band of 0.08–0.5 Hz at station BREF. Amplitude (m/s) is indicated by a vertical bar at the left of each seismogram. The first waveform is the largest VLP/LP event analyzed by Kumagai et al. (2010), which is our reference event to determine magma volume.</figcaption></figure>
Frequency band 0.08–0.5 Hz

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

**Fig. 5.** Temporal variations in UD amplitudes and amplitude ratios of EW/UD and NS/UD of VLP signals in VLP/LP events in a frequency band of 0.08–0.5 Hz from 2009 to 2015. (a, b, c) BREF. (d, e, f) BTAM. (g, h, i) BVC2. (j, k, l) BMOR. The averages and standard deviations of the EW/UD ratios were 1.1 ± 0.7, 1.5 ± 0.4, 1.7 ± 0.5, and 1.6 ± 0.6 at BREF, BTAM, BVC2, and BMOR, respectively. Those of the NS/UD ratios were 0.9 ± 0.4, 1.3 ± 0.7, 1.7 ± 0.9, and 2.2 ± 0.8 at BREF, BTAM, BVC2, and BMOR, respectively.

818

819

820

821

822

823

824

825

826

827

828

829

830
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** (a) Accumulated amplitude and magma volume estimated of VLP signals in VLP/LP events in a frequency band of 0.08–0.5 Hz at BREF from 2009 to 2015 (see the text for details), overlapped with tremor amplitudes in arbitrary unit (a.u.). (b) Spectrogram of continuous vertical velocity waveform data at BREF between January and December 2015. Phases A–G were defined through the accumulated magma volume and VLP/LP activity. Phases 1–4 in (b) were defined by Bernard et al. (2016) based on tremor activity.</figcaption></figure>

831

832

833

834

835

836

837
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Comparison of VLP/LP activity (red triangles) with tremor (black vertical bars), SO₂ emission (blue lines), and ash fallout mass (green squares). Tremor amplitudes in arbitrary units (a.u.) and ash fallout mass are from Bernard et al. (2016), while SO₂ emission data are from Hidalgo et al. (2018) (online dataset at https://ecl.earthchem.org/view.php?id=1165). (a) Temporal variations in VLP amplitudes in VLP/LP events in the 0.08–0.5 Hz band, tremor amplitudes, and ash fallout mass. (b) Temporal variations in the daily number of VLP/LP events, tremor amplitudes, and SO₂ emission.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Temporal variations in ground deformation observed by the tiltmeter at VC1 (radial and tangential components to the crater shown by green and blue lines, respectively) and the daily number of VLP/LP events (red vertical bars). The period of pre-eruptive tremor is shown by gray area, and the first eruption on August 14 is indicated by an arrow. Peaks of VLP/LP activity on June 19 and August 25 are also shown by arrows. Vertical dashed lines emphasize the pre-eruptive LP events and phases E through G (see text for details).</figcaption></figure>

847

848

849

850

851

852

853

854
<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Conceptual model for source processes of the VLP/LP events observed at Cotopaxi Volcano. See Section 5 for details.</figcaption></figure>

855

856

857

858
859

SUPPORTING INFORMATION

860

**Appendix A.** Amplitude spectra of a VLP/LP event and accumulated amplitudes and magma volumes.

861

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

862

863

864

865

866

**Fig. A1.** Amplitude spectra of a VLP/LP event at Cotopaxi, band-passed with a Butterworth filter at 8 distinct pass bands, with the low-frequency corner fixed in all instances (see legend). This event occurred on January 14, 2009 at 18:54 and is the largest VLP/LP observed at Cotopaxi (see text for details).
Frequency band 0.08–0.5 Hz

<figure><img src="image_13.png" /><figcaption>**Fig. A2.** Accumulated amplitudes of VLP signals in VLP/LP events in a frequency band of 0.08–0.5 Hz at the individual stations from 2009 to 2015 and accumulated magma volumes estimated from the accumulated amplitudes at the individual stations.</figcaption></figure>

867
868
869
870
871