Since we are interested in long-term trends, the SO₂ flux measurements presented here are the monthly average values. Campion (2014) showed that a tight linear correlation exists between the monthly averaged SO₂ mass around the volcano and the monthly average of the SO₂ fluxes calculated on a daily basis over the same month using the traverse method. The R ² correlation coefficient of these correlations is usually superior to 0.8, and their slope coefficient is specific for each volcano, depending on the plume transport and atmospheric conditions. The SO₂ fluxes presented in this study were calculated from the monthly averaged SO₂ mass after calibrating the mass to flux relation for each volcano over a period of several months of fluctuating SO₂ flux (see calibration details for each volcano in Supplementary Material). When the SO₂ emissions are too small to cause a monthly SO₂ anomaly that stands consistently out of the sensor noise, we average the images over a year and calibrate with the flux data of Carn et al. (2017). We used the literature values in the cases where no OMI data were available (Nyiragongo during 2002 and 2003)

In this work, we used the MODIS data processed by the MIROVA algorithm (Coppola et al., 2016b; https://www.mirovaweb.it/). This system processes the infrared images at 1 km of resolution acquired daily by the two MODIS sensors on board the polar satellites Terra and Aqua (EOS-NASA missions). The algorithm is built to detect and locate sub-pixel thermal anomalies sourced by the various types of high-temperature volcanic activity, including lava lakes (Coppola et al., 2016a). Once a thermal anomaly is detected the associated VRP is calculated by using the “MIR method” (Wooster et al., 2003), an approach based on the analysis of the mid-infrared (MIR) radiance in excess to the background (see Coppola et al., 2016a for details). The VRP is calculated in Watts (W) and represents a combined measurement of the area of the volcanic emitter (i.e., the lava lake area) and its effective radiating temperature. The MODIS MIR band (Band 21 or 22) is centered at 3,959 μm, with a bandwidth ranging from 3,929–3,989 μm. The CO₂ absorption band in the MIR is centered at 4.2 μm and does not produce a significant absorption in the MODIS band 21. H₂O does not absorb in that wavelength range either. SO₂ does have an absorption band (ν1 +ν2, stretching from 3.95 μm to 4.07 μm) that overlaps marginally on Band 21, but this is a weak absorption band (HITRAN database, https://hitran.org/lbl/, Rothman et al., 2013). Our calculations show that a 200 m thick plume with 200 ppm should produce a band-integrated absorbance of less than 5%. Band 21 can therefore be considered quite transparent to volcanic gases and the atmosphere. Aerosols are more of a concern. The plumes emanating from lava lakes have a quite variable optical thickness that depends strongly on the atmospheric humidity and mixing pattern (see Matsushima and Shinohara, 2006). Several volcanoes in this study are located in the tropical regions where meteorological clouds or water droplets condensing in the plume often mask, partially or totally, the thermal radiance emitted by the lake. To cope with this variable underestimation factor, the VRP values that are used in this study are the maximum values recorded each month. This approach is based on the assumption that most of the data variability over a month is due to the partial masking of the IR radiation by the clouds or the condensed plume so that the maximum VRP of each month is recorded when the volcano is observed on the optimal viewing conditions. We checked this hypothesis by comparing the variability of the monthly means, weekly maxima and monthly maxima calculated for a time series of VRP from Nyiragongo volcano during year 2010, where field observations showed a very stable activity and constant lake dimension. However, we acknowledge that each lava lake’s activity has its own variability and that the abovementioned approach may lead to a slight overestimation of the “true” mean VRP radiated by the lava lake.

Erta Ale (Ethiopia, 13.6°N, 40.67°E) is a broad tholeiitic shield volcano located on the spreading axis of the Danakil Rift. Its 0.7 × 1.6 km, elliptical summit caldera hosts two pit craters. The northern pit crater contained a lava lake until the mid-seventies and has been the seat of permanent fumarolic activity and occasional lava spattering or small intra-crater flows since then (Oppenheimer et al., 2004). The smaller southern pit crater has contained a lava lake discovered in 1969 (Le Guern et al., 1979). The southern lava lake of Erta Ale has shown large variations of activity, size, and level, as shown in Figures 2A, B and sometimes overflowed the rim of its pit, feeding lava flows that have covered the caldera floor (e.g., Field et al., 2012). In January 2017, after a year-long gradual increase in the lake level, with several overflowing events, a dike-fed fissure eruption on the southern flank produced copious lava flows for more than 3 years. The eruption covered ∼30 × 10⁶ m² with >1 × 10⁸ m³ of lava but it did not lead to the disappearance of the lava lake (Global Volcanism Program, 2018) and did not cause any major collapse in the caldera. Instead, Xu et al. (2020) reported a very wide deflation, extending well beyond the limits of the volcanic edifice, which they attribute to the extraction of magma from a mid-crustal storage area. The long duration of the eruption, the absence of a caldera collapse and the persistence of the lava lake indicate that it was not a draining event, although the lava level rise prior to it and drop during it, indicate some degree of hydraulic connection with the lava lake.

The VRP time series, in Figure 2C, shows a relatively weak thermal emission which shows strong variations, ranging between 10⁷ and 1.5 × 10⁸ W. The stronger VRP spikes that exceed the background values correspond to short-lived lake overflows, either inside the pit crater or on the bottom of the caldera. During a small period in 2005, the thermal emission disappeared completely. The SO₂ flux shows small values, between 0.3 and 1.2 kg/s and evolves grossly in a similar way as the VRP. The relatively low values of VRP and SO₂ flux emitted by Erta Ale are in agreement with the previous field studies (Le Guern et al., 1979; Oppenheimer et al., 2004)

Nyiragongo, (1.52°S, 29.25°E, Democratic Republic of Congo) is a foiditic stratovolcano located in the Western Branch of the East African Rift. It has hosted a permanent lava lake that was discovered in 1948 (Tazieff 1974) and disappeared in a catastrophic drainage event in 1977 (Tazieff 1977). Hundreds of people were killed by the fast-moving lava flows rushing downslope at initial velocities estimated close to 30 m/s. The crater suffered a massive collapse as a result of the drainage of the upper magmatic system. Intra-crater eruptions in 1982 and 1994–1995 gradually refilled the crater but failed to form a permanent lake (Tazieff, 1994). In January 2002, a new lateral eruption triggered by a regional dyking and rifting event drained the upper magmatic system (Komorowski et al., 2002). The fast-flowing lava flows invaded the city of Goma, again killing hundreds of people. In the summer of 2002, a new lava lake gradually formed, and has since then enlarged considerably to a maximum size of ∼240 m diameter (Burgi et al., 2014; Valade et al., 2018). Figures 3A, B illustrate the growth and rise of the lava lake between 2003 and 2017. The level of the lake followed a step-wise increasing trend, to reach an altitude comparable to that of 1972–1977, considered as critical for the risk of draining (Burgi et al., 2020). Once again in June 2021, , fissures opened suddenly on the south and southeast flank of the volcano, releasing lava flows that rapidly reached the outskirt of Goma (Global Volcanism Program, 2021; Smittarello et al., 2022). The lava lake drained completely and has not yet reformed at the time of writing. The OMI data, which started in late 2004, do not cover the SO₂ emissions in the first 2 years of the lake existence. In order to fill this gap, the SO₂ flux dataset was completed with the ground-based DOAS measurements of Sawyer et al. (2008), flux measurements from three ASTER (Advanced Spaceborne Thermal Emission and Reflection) images with the band ratio method of Campion et al. (2010), and the data obtained by Carn (2004) with the TOMS (Total Ozone Mapping Spectrometer) spaceborne sensor.

The SO₂ data (Figure 3C) show a 2 year period of extremely high SO₂ emissions (100–250 kg/s) at the beginning of the lava lake’s existence, followed by an initially rapid and then slower decrease to reach, by 2009, a relatively stable value of ∼10 kg/s that has persisted until the disappearance of the lake. This decrease in the SO₂ flux coincides with a continuous increase of the VRP that materializes the growth of the surface area of the lake and the rise of its level (Burgi et al., 2014). Spikes in the thermal power graph are due to the frequent overflowing events that accompanied the lava lake rise from ∼800 m below the crater rim in 2002 to 400 m in 2011 (Burgi et al., 2014). A drastic drop of both SO₂ flux and VRP marks the May 2021 drainage. The activity since then has not completely disappeared but the produced VRP has been weak and discontinuous so far. SO₂ fluxes have stayed much lower than before the lava lake drainage, and two orders of magnitude lower than in 2002 and 2003, when the former lake started to form. We conclude that the activity since 2021 has not yet reformed a permanent lava lake.

Nyamuragira, (Democratic Republic of Congo, 1.408°S, 29.2°E), is a basanitic shield volcano located 15 km northeast from Nyiragongo. It has hosted a lava lake between an uncertain date at the beginning of the 20th century and a massive flank eruption that drained the lake in 1938 (Verhoogen, 1948). From 1948 to 2012, Nyamuragira produced frequent voluminous flank eruptions (e.g., Wadge and Burt, 2011). In 2012, after a major lateral eruption that partially drained the upper magmatic system of the volcano, a new lava lake started to form in a new collapse pit in its summit caldera (Campion 2014; Coppola et al., 2016b). The formation of the lava lake started with a 2.5 years period of intense SO₂ emission during which only very weak radiant heat emissions were infrequently detected (Figure 4C). A fountaining lava lake was observed from the ground for the first time in June 2014 (Figure 4A) and starting in October 2014, the SO₂ flux diminished and the radiant heat emissions increased and became more frequent, as the lava lake had become larger and more stable (Figure 4B). The lava lake has passed through two periods of very reduced activity, from June 2016 to December 2016 and from June 2017 to April 2018, which has shown its surface becoming completely crusted over (Global Volcanism Program, 2017, Global Volcanism Program, 2019), although the emission of SO₂ never stopped. The scarce cloud-free images that ASTER took since the end of these reduced activity periods show weak SO₂ emissions, below or barely above detection limit, which yield an SO₂ flux inferior to 3–5 kg/s, which is dwarfed by the emissions of neighboring Nyiragongo. Burgi et al. (2021) documented a regime of piston-like rising and sinking of the mostly crusted lava lake surface since 2016, possibly differing from the classical convection scheme that is thought to feed lava classical lakes.

Kilauea (19.421°N, 155.287°W, Big Island of Hawaii, United States of America), is a tholeiitic shield volcano. It hosted a permanent lava lake in the Halemaumau pit crater of its summit caldera until 1924, when it disappeared during a collapse event and explosive eruption, possibly associated with a submarine eruption. Since 1983, Kilauea was in a state of constant effusion of lava through the Puu O’o vent on its East Rift zone. After years of summit inflation, a new lava lake formed gradually, starting in March 2008, in a new pit crater located in the Halemaumau summit caldera (Poland et al., 2010). During its first month of existence, the only surface expressions of the lava lake were the extraordinary degassing and fresh fluidal basaltic glass fragments that were lofted by the gas plume. By the end of 2008, the gas emissions started to decrease and the lake level started to rise, becoming more frequently visible. It often showed rapid filling and draining cycles, or sometimes had the aspect of a cascade. The pit containing the lava lake enlarged over the years by frequent collapses of its rim, resulting in a continuous increase of the lava lake size (Patrick et al., 2013). Figures 5A, B illustrate the spectacular growth of the lake between 2008 and 2017. The large changes of the level of the lava lake were caused by changes in the pressure of the magmatic system of the volcano (Patrick et al., 2015), while the smaller short-term changes were caused by the cycles of foam accumulation below the lake skin and outgassing through spattering (Nadeau et al., 2015). Our measurements with OMI (Figure 5C), confirms extremely high SO₂ flux values during the first year of the existence of the lava lake in Halemaumau, reported by Beirle et al. (2014). The monthly averaged flux peaked at ∼300 kg/s, whereas the VRP remained very low. The VRP data show a strong long-term increase that started shortly after the extreme SO₂ emissions of the initial phase had declined, and stabilized around 2016. A short but very prominent spike in VRP appeared in April 2018, without any increase in the SO₂ flux. It corresponds to a series of large lake overflows that were later considered as a precursor to the massive May-August 2018 eruption which drained the lava lake and summit magma reservoir into Kilauea’s Lower East Rift Zone. Although two effusive eruptions have partly refilled the collapse caldera left by the 2018 draining eruption, they have not reformed a permanent lava lake.

Ambrym Volcano, (Republic of Vanuatu, 16.25°S, 168.12°E) is a large island arc basaltic to dacitic volcano crowned by a 12 km wide caldera that contains several cones and craters. Two of these cones, Marum and Benbow (Figures 6A, B), have been hosting lava lakes semi-permanently for at least 3 decades (Allard et al., 2016; Coppola et al., 2016c). Short-lived lateral eruptions in the caldera, like those of 1988 and 2015, or in the rift zones, like those of 1914 (Nemeth and Cronin, 2011) and 2018 have partially or completely drained the lava lakes. The degassing of Ambrym’s lava lakes was extremely vigorous (e.g., Bani et al., 2009; Allard et al., 2016), maintaining its surface permanently bubbling and agitated, and placing Ambrym as the long-term strongest permanent emitter of Sulfur Dioxide. However, following a major dyking episode and submarine eruption in December 2018, the lava lakes of Ambrym have disappeared (Shreve et al., 2019). The time series of total VRP emitted from Ambrym, combining the heat emission from both Marum and Benbow, shows a waxing and waning trend with a clear increase since 2012, leading to a fissure eruption in February 2015, and a marked decrease immediately after it (Figure 7C). The time series of the SO₂ flux shows very high values, with a prominent peak in 2004–2005, followed by ups and downs usually uncorrelated with the thermal emissions. Interestingly, a general decrease in the SO₂ flux is observed during the 3-year-long period of VRP increase preceding the 2015 eruption (Coppola et al., 2016c).

Masaya volcano (Nicaragua, 11.98°N, 86.16°W), is a large complex of calderas in central Nicaragua formed by powerful basaltic plinian eruptions (Perez and Freundt, 2006). A small shield volcano crowned by a complex of juxtaposed cones and pit-craters has built since the last caldera-forming eruption (Rymer et al., 1998). It has hosted lava lakes episodically since before the Spanish conquest, and also produced several large effusive effusions (Viramonte and Incer-Barquero, 2008). Masaya has been known as a strong permanent emitter of SO₂ and other gases even when no lava lake was visible at the surface (Delmelle et al., 2002). Starting in December 2015, a lava lake appeared in the Santiago crater of Masaya volcano and gradually enlarged until May 2016, when it reached its maximal size of 40 × 30 m (Figure 6A). Aiuppa et al. (2018) reported the changes in gas flux and composition that accompany its formation, most notably an increase in the SO₂ flux and the CO₂/SO₂ ratio. The low altitude and small size of the Masaya plume make it difficult to measure with OMI on a monthly base. We, therefore, averaged images over 1-year periods to produce the long-term trend of SO₂ fluxes. Figure 7B shows the SO₂ flux and VRP from Masaya. SO₂ emission is permanent but shows a pronounced waxing and waning pattern, a behavior that has characterized the volcanoes since at least the seventies (Delmelle et al., 1999). A strong and rapid increase of both the SO₂ flux and the VRP marked the appearance of the lava lake in 2016, followed by a gradual decrease over the next 4 years. The detection of low thermal anomalies is occasional before the appearance of the lava lake. These anomalies occurred predominantly during the period of moderately high SO₂ flux that prevailed during 2007–2010 and have low VRP values.

Villarrica volcano (Chile, 39.42°S, 71.93°W), is a medium-sized ice-clad stratovolcano built on two recent calderas associated with the basaltic plinian eruptions and pyroclastic flows (Lohmar et al., 2007). Its narrow summit crater is permanently degassing and is often the seat of mild strombolian and lava lake activity (Palma et al., 2008; Moussallam et al., 2016), but has also produced several stronger eruptions (VEI 2-3, Aiuppa et al., 2019). The daily variability of VRP from Villarrica is the highest among the volcanoes studied in this work. Additionally to the cloud-related variability, it might also result from the frequent bursting of large bubbles at the surface of the lake that eject abundant lava spatters out of the lava lake, and/or rapid cycles of building and collapse of cooled scoria benches that partially covered the lava lake. In the long term, the VRP time series (Figure 8C) highlight the intermittency of Villarrica’s lava lake. The periods of heightened activity lasting 2 to 3 years and radiating up to 3 × 10⁸ W are separated by intervals of nearly complete quiescence also lasting 2 to 3 years. The typical SO₂ flux from Villarrica is too low to be measured on a daily or even on a monthly basis by OMI. So, in this work we used the time series from Carn et al. (2017) and extended it until 2020. It follows a temporal trend that is strikingly parallel to the VRP time series. This cyclic behavior of both the VRP and the SO₂ flux corresponds to periods of alternating periods high and low surface, level, and activity of the lava lake (Figures 8A, B).

Comparing quantitatively lava lakes together is not a straightforward task because their heat and gas emissions span more than two orders of magnitudes (see Table 1). One approach is to compare intensive parameters (i.e. that are not dependent on the lake size). One of them is the equivalent radiating temperature T _AR , which is calculated using the inverse of the Stefan–Boltzmann law: T E R = V R P / A . σ . ε , 4 (1) where A is the lake area, ε its emissivity (0.95 for basaltic lava), and σ is the Stefan–Boltzmann constant (5.668 10⁻⁸ W/m²K⁴). The distinction made by Lev et al. (2019) is particularly apparent in Table 1, which clearly shows two populations of T _ER , which correspond to the chaotic lakes (T _AR between 1,100 and 1,200 K) and the ordered lakes (T _ER between 800 and 900 K). In Figure 9, we plot the apparent radiating temperature against the SO₂ flux per unit area of the lake. This shows an expected positive correlation between these two parameters, demonstrating that the chaotic lakes are not only hotter but also richer in gas. The explanation for this is that the T _ER is dominantly dependent on the percentage of the lake that is covered by a cooling skin and not by the absorption of the gas plume overlying the lake. The gas-poor lakes do have such a crust while the permanent bubbling at the surface of gas-rich lakes hampers its formation.

More information can be gained by comparing extensive parameters. The VRP vs. SO₂ flux diagram in Figure 10 shows that two leagues can be distinguished, those with small SO₂ and heat emissions (including Villarrica, Erta Ale, Masaya and likely Erebus although it was not included in our study), and those whose SO₂ or VRP emission are higher than either 25 kg/s or 5 × 10⁸ W (Nyamulagira, Nyiragongo, Halemaumau, and Ambrym). Different behaviors seem to characterize these two leagues.

The group defined as small lakes contains volcanoes that are known for their permanent degassing over very long periods. Their lava lakes, when visible, have surface areas that vary between a few tens and a few thousands of square meters. They show coupled variations of their SO₂ flux and VRP, which maintain a certain linear correlation over time. In all cases, these fluctuations correspond to large changes of the height of the magmatic column, sometimes to the point of disappearing from sight in the case of Masaya and Villarrica, or overflowing in the case of Erta Ale. These lakes do not seem to suffer transitions from the chaotic type (to which Masaya and Villarrica belong) to the organized type (to which Erta Ale belongs). The most likely explanation for the long-term coupling of their VRP and SO₂ flux is that both are linked to the flow rate of the magmatic convection that feeds the lake. If this is the case, the slope of the regression line in the VRP versus SO₂ plots should depend on the initial sulfur content of the convecting magma and the efficiency of the convective process to dissipate heat.

Instead, large lava lakes, especially at the early stages of their existence show a strong anti-correlation between their heat and gas emissions. Within large lava lakes, two groups can clearly be distinguished. The first group is characterized by very high gas emissions and small to moderate VRP. It regroups the early-stage lava lakes of Kilauea, Nyiragongo, and Nyamuragira and the lava lakes of Ambrym. Visually, these lakes have a small size, brightly incandescent and vigorously bubbling surface, are stirred by very fast currents, or may even look like a small permanent lava fountain (e.g., Smets et al., 2014). They belong to the “Chaotic” category of Lev et al. (2019). The second group is characterized by high heat emissions and comparatively low gas emissions, and contains only mature lava lakes from rifts and hotspot volcanoes [“Organized” category according to Lev et al. (2019)]. The visual observations show that more than 80% of their surface is covered by a cooling skin or crust that is stirred by relatively slow convection currents until it sinks back into the lake.

The strong decoupling of SO₂ flux and radiant flux during the early stage of all three lava lake formations at Nyiragongo, Kilauea, and Nyamuragira is at a first glance surprising result, as most of the previous studies on lava lakes have reported or assumed a correlation between SO₂ emission and lava lake surface (e.g., Oppenheimer et al., 2004; Sweeney et al., 2008). It is, however, easily explained if we consider that the upwelling material feeding the lake is a magmatic foam. We show in the next section that the thermodynamic and phenomenological implications of this hypothesis match all of the observations reported on lava lakes.

A lava lake can be considered, over a large period of time, as a steady state system, where the thermal power brought by the magmatic convection Q is equal to the sum of three following loss terms (Harris et al., 1999). Q = V R P + C L C + C L A + Q g . (2)

The VRP, the radiated thermal power from the lake surface, is the dominant heat loss mechanism at these very high temperatures (Harris et al., 1999). The conductive heat loss through the wall of the conduit (CLC) can be considered constant and negligible for a well-established steady-state system (Henley and Hughes, 2016). The conductive heat loss to the atmosphere (CLA) amounts to 10 (for hot lakes) to 30% (for colder lakes) of the VRP (Harris et al., 1999) and is also proportional to the lake surface area. Q_g is the heat lost by the transfer of hot gases from the lake to the atmosphere, and its magnitude depends on the total gas flux emitted by the lava lake.

The short-term excursions of the system away from this equilibrium can occur when the magma flux F_m feeding the system at depth varies. When F_m decreases, the heat transported by the magmatic convection is insufficient to balance the heat loss and a portion of the lake’s surface solidifies, reducing the exposed radiating surface until the heat loss is again balanced. This leads to the formation of a bank of solidified lava. When F_m increases, the convection gets more vigorous, and the gas flux also increases so that more bubbling and more incandescent surfaces are exposed to the atmosphere. As a result, the heat loss increases until it balances the heat transported by the convection. Lake may enlarge by thermal erosion of its walls, its level tends to rise (because of the increased volume occupied by the rising magma and its gas) possibly until the point of overflowing, as often observed at the lava lakes of Erta Ale, Nyiragongo, Nyamuragira, and Kilauea.

If the material feeding the lava lake is considered as a bubble-rich liquid, or even as foam, composed of a volumetric fraction of magmatic gas X_g and a fraction of silicate melt equal to 1-X_g. The transported thermal energy flux can be calculated as Q = F m . 1 − X g . ρ l . C p l . T u − T d + X S · ∆ H C S + X g · ρ g . C p g · T u − T a , (3) where F_m is the volumetric flux of rising magmatic foam, ρ _l and ρ _g are the densities of the melt and gas, respectively, Cp_l and Cp_g are their specific heat, X_S is the proportion of melt that solidifies before sinking back, ΔHC_S the latent heat of solidification of the melt and T_U, T_D, and T_A are the temperature of the upwelling magma, of the downwelling lava, and of the atmosphere, respectively. Because silicate melts have a much higher density than gases, the thermal energy transported by the melt fraction of the foam is typically one order of magnitude higher than what its gaseous phase can transport. When magmatic gases escape into the atmosphere, they cool rapidly due to their mixing with air. Because of their very low emissivity at the wavelengths used in the MIROVA algorithm (HITRAN database, https://hitran.org/lbl/, Rothman et al., 2013), the gas, once emitted, does not contribute to the radiated power. Additionally, the rise from depth of a bubbly magmatic liquid implies that part of its thermal energy is consumed by the adiabatic expansion of the gas phase during its ascent (Stevenson and Blake 1998). Although the related temperature decrease may be compensated by the latent heat released by the partial crystallization of the foam’s liquid phase (La Spina et al., 2016), the total energy available for radiation into the atmosphere is anyway decreased.

The consequence of all this is that the more gas the foam contains at the bottom of the convection cell, the less thermal energy it can transport to the surface and radiate into the atmosphere. This point is the cornerstone to understand the anticorrelation of heat and gas flux during the early stage of lava lake formation. The gas-rich lava lakes are fed by a material that can be considered true foam, which likely contains more than 80% of gas. The consequences of this are as follows:

1) The density contrast between the undegassed upwelling foam and the degassed downwelling lava is extremely high, which results in a very dynamic convection, with fast superficial currents, a lower exposure time to the atmosphere, and less surface cooling for the downwelling lava.

If our hypothesis that the SO₂ vs. VRP anticorrelation is due to the changes in the gas/melt ratio of the material feeding the magmatic convection, then why do small lakes show a correlation between their SO₂ and VRP? The answer to this question lies probably in the presence or absence of a large and shallow magma storage system, where SO₂-rich gas can exsolve and separate from the magma to accumulate as foam. We propose that large lava lakes have a large shallow magma storage system and small lava lakes do not. The presence of a relatively large shallow magma storage area beneath large lava lakes is also the reason why they are affected by catastrophic draining events. These events release massive lava flows with colossal effusion rates and empty, totally or partially, the shallow magma storage associated with the lava lake, to the point of generating massive collapse and smothering for the next several years of the volcanic activity, as happened in Kilauea in 2018 (Patrick et al., 2020), Nyiragongo in 1977 (Tazieff, 1977) and 2021 (Smittarello et al., 2022), and Ambrym in 2018 (Shreve et al., 2019). Small lakes may have storage areas underneath them (Stephens and Wauthier, 2018 for Masaya; Xu et al., 2020 for Erta Ale) but they must be either small or located deeper than the exsolution depth of the SO₂ in their magma.

In this section, we highlight a four-stage sequence of processes leading to the formation of a lava lake at the rift and hotspot volcanoes. This sequence is schematized in Figure 11.

The first stage may start years to months before the first visible signs of lava lake formation appear. This stage is the previous accumulation of magma in the plumbing system of the volcanoes, when the upper conduit system of the volcano is still closed and does not allow the release of high-temperature magmatic gases to the surface. Several seismic and/or deformation studies have suggested that significant magma accumulation had taken place before the formation of a new lava lake at Nyiragongo, Kilauea, and Nyamuragira. Intrusion of a massive dyke occurred, after a year of seismic unrest, in January 2002 at Nyiragongo volcano (Wauthier et al., 2012), and preceded by 6 months the appearance of the new lava lake. A pulse of increased magma supply and accumulation was documented by deformation, seismicity, and gas measurement at Kilauea during 2003–2007 (Poland et al., 2012; Wauthier et al., 2013). Large-scale inflation was detected before the 2011 and 2012 eruption of Nyamuragira volcano (Albino et al., 2013) which marked the beginning of the formation process of Nyamuragira volcano. Simultaneously to this magma accumulation, a magmatic foam layer rich in gas is likely to accumulate at the top of the magma storage area (Jaupart and Vergniolle, 1989; Wallace 2001; Allard et al., 2005).

The second stage is the development of an open conduit system that allows magmatic degassing from the summit. In all three cases, the opening of the system was achieved through the migration of magma through laterally propagating dikes, the partial drainage of the upper plumbing system of the volcano, and the formation of a pit crater. At Nyiragongo, this second stage was materialized by the catastrophic eruption of January 2002 and the massive collapse of the caldera floor. At Kilauea, it was much less dramatic, being just materialized by a small fissure eruption in June 2007, generalized deflation (Poland et al., 2012), and the formation of an initially very narrow pit crater. This is probably because magma was already injecting continuously through its rift zone to feed the long-lived Pu’u O’o eruption. At Nyamuragira, this second stage was characterized by the large 2011–2012 eruption, one of the most voluminous in its history, and the collapse of a pit crater in its caldera (Campion, 2014; Coppola et al., 2016b).

The third stage is characterized by a strong degassing with little to no magma yet reaching the surface. This stage may initiate simultaneously to the lateral migration of magma and formation of the collapse pit, as was probably the case at Nyamuragira (Coppola et al., 2016b), or be delayed by a few months as in Nyiragongo and Kilauea. During this stage, the SO₂ emission rate typically reaches values of more than 100 kg/s, whereas the lava lake remains very deep in its newly formed collapse pit. If visible, the proto lava lake may look, at this stage, like a small lava fountain (Smets et al., 2014) or like an endless lava cascade (HVO photo archive, 2010), or be subject to huge level variations. All these phenomena result from the dynamics of the gas-rich magmatic foam that is feeding the lake. The duration of this stage of high gas/low heat emission probably corresponds to the time needed for the accumulated foam to be outgassed by the proto lava lake.

The fourth stage, characterized by a high and stable ratio of radiated heat to emitted gas corresponds to a fully mature lava lake fed by magmatic convection. Although the transition from stage three to four is usually progressive, it seems that the drop in SO₂ flux is more rapid (a few months) than the increase in VRP, which spans several years. This materializes the exhaustion of the pre-accumulated foam, the resulting increase of the melt/gas ratio in the material feeding the lake, and possibly thermal erosion of the conduit feeding the lava lake. The convection speed slows accordingly because the density contrasts between the ascending magma and the descending degassed lava decrease. The bubble bursting at the top of the magma column decreases and a cooling skin can start to form. It is likely that the magma feeding these mature lava lakes has already partly degassed during stage 3 and/or is a mixture between the partly degassed magma from stage 3 and the juvenile magma that is continuing to enter the system.

The formation sequence described in the previous section applies to the rift and hotspot volcanoes, which are fed by a volatile-poor basaltic magma. Primitive basaltic magmas from arcs are richer in volatiles and more commonly produce explosive basaltic eruptions. However three of the volcanoes studied here, Villarrica, Ambrym, and Masaya, do host most of the time a lava lake. Interestingly, they all share the following three characteristics:

1) They are associated with large calderas (>5 km) that have produced large plinian basaltic eruptions in recent times (<10 ky)

We can, thus, reasonably suspect that arc lava lakes, even though they all belong to the gas-rich type, are fed by a previously degassed magma or more likely a mix between a large batch of degassed magma and small repeated injection of juvenile magma, as suggested by the strong variations in our time series of SO₂ flux from Ambrym, Villarrica, and Masaya and by recent gas composition time series at these volcanoes (Allard et al., 2016; Aiuppa et al., 2017; Aiuppa et al., 2018, respectively).