Overﬂows and Pyroclastic Density Currents in March-April 2020 at Stromboli Volcano Detected by Remote Sensing and Seismic Monitoring Data

: Between 28 March and 1 April 2020, Stromboli volcano erupted, with overﬂows from the NE crater rim spreading along the barren Sciara del Fuoco slope and reaching the sea along the NW coast of the island. Poor weather conditions did not allow a detailed observation of the crater zone through the cameras monitoring network, but a clear view of the lower slope and the ﬂows expanding in the area allowed us to characterize the ﬂow features. This evidence was integrated with satellite, GBInSAR, and seismic data, thus enabling a reconstruction of the whole volcanic event, which involved several small collapses of the summit cone and the generation of pyroclastic density currents (PDCs) spreading along the slope and on the sea surface. Satellite monitoring allowed for the mapping of the lava ﬂow ﬁeld and the quantiﬁcation of the erupted volume, and GBInSAR continuous measurements detected the crater widening and the deﬂation of the summit cone caused by the last overﬂow. The characterization of the seismicity made it possible to identify the signals that are associated with the propagation of PDCs along the volcano ﬂank and, for the ﬁrst time, to recognize the signal that is produced by the impact of the PDCs on the coast.


Introduction
Rapid changes of the surface morphology often occur in open-conduit basaltic volcanoes that frequently erupt, such as Etna and Stromboli (Italy), Piton de la Fournaise (La Réunion Island), or Kilauea (Hawaii). Cinder-cones~60 m high can form in just one or a few weeks [1,2], large lava flows can spread over roads or villages [3][4][5][6], and summit collapses following major lava withdrawal can involve large areas and result in new calderas [6][7][8][9]. However, as they are rapid in building new reliefs, the often loose and unstable mixture of breccia, ash, and spatter, which accumulates quickly known as the "Lighthouse of the Mediterranean", because of its persistent explosive activity from the summit craters, with bursts occurring every few minutes [29][30][31]. The summit crater of the volcano is a depression ~300 m long in a NE-SW direction (Figure 1c), ~50 m wide, and ~50 m deep, located at ~750 m elevation [14,19]. Three crater areas are located within the summit depression: the NE crater zone (NEC), the Central crater zone (CC), and the SW crater zone (SWC), each of them comprising a variable number of active vents (Figure 1c). The capacity of the uppermost feeder conduit increased after the 2002-2003 and 2007 flank eruptions [20,32], but sudden changes in the magma level may result in a greater magmastatic pressure, which could cause the destabilization and collapse of portions of the summit cone [14,33]. The aim of this paper is to present our study of the eruptive activity occurred at Stromboli between 28 March and 1 April 2020, along with our interpretation and quantification of the eruptive processes that are based on the analysis of monitoring data, comprising time-lapse videos recorded by the camera network, satellite images, GBInSAR, and seismic stations.

Methods
The description of the events, the calculation of velocity for lava flows and PDCs spreading along the SdF and on the sea surface, and the count of the explosions occurring over time from the summit craters were obtained through the analysis of the videos recorded by the network of fixed monitoring cameras maintained by INGV-OE (Istituto Nazionale di Geofisica e Vulcanologia-Osservatorio Etneo). The volcano deformation was measured through two Ground-Based Interferometric Synthetic Aperture Radar (GBInSAR) devices that were installed by the University of Florence. Seismicity was analyzed using data from the broadband seismic network that was installed by INGV-OV (Istituto Nazionale di Geofisica e Vulcanologia-Osservatorio Vesuviano). Lava flow field area and volume, as well as an estimation of the effusion rate, were computed by the TechnoLab of INGV-OE while using multispectral infrared and optical satellite data.

The INGV Cameras Monitoring Network
The INGV cameras monitoring network at Stromboli volcano in March-April 2020 comprised three fixed instruments, two thermals (SCT and SPCT) and a visual (SQV). Their details are listed in Table 1 and their location is shown in Figure 1. SQV acquires at a frequency of one image every two seconds, SCT one image every second, and SPCT two images every second. The difference in acquisition frequency, as well as in the distance from the crater area and viewing angle, result in a different count of the explosions, depending on the camera taken into consideration. The number of explosions is thought to represent an expression of the magma level within the feeder conduit [19,39,49,50]. We have manually counted the total number of events occurring within the whole crater area using only the thermal cameras (SCT and SPCT), because they allowed a comparison between day and night views. However, the presence of clouds and/or dust may limit or hide visibility, as happened during most of the day on 30 March and in the early h on 1 April 2020. The calculated velocity of PDCs and lava flows along the SdF slope, obtained from the images of the monitoring cameras, are average velocities, because they were calculated dividing the whole travelled path by the time. Some of the images were affected by shadows due to fog, clouds, or ash, resulting in a non-well-defined outline of the features. The cumulative error on these spatial measurements, which was due to the poor quality of some frames, is about 2%, and it was obtained from the number and size of uncertain pixels. All of the times are expressed here as UT.

GBInSAR
Measuring surface deformation, exploiting the phase difference between two spaceborne SAR images (differential InSAR, DInSAR; [51]) makes it possible to recognize ground displacements along the satellite line of sight (LOS) direction on a centimeter-scale. Processing a long stack of images using multi temporal (MT) InSAR techniques allows for the detection of millimeter-scale displacements over long time frames through the reduction of error sources [52,53]. GBInSAR has the additional advantage of producing frequent SAR images (on the order of seconds to minutes), resulting in very high frequency slope maps and time series. Moreover, the use of GBInSAR in the Ku-band (17.0-17.1 mm radar) can penetrate dust clouds, abundant especially during collapse events, and can work with variable light and atmospheric conditions [14]. The NE portion of the summit crater terrace at Stromboli and the northern portion of the SdF are monitored by two GBInSAR devices, which are located in a stable area N of the SdF (Figure 1). The first GBInSAR (GBInSAR NE400; Model: GB-InSAR LiSALab, Ellegi srl, Rovello Porro, Italy, http://lisalab.com/home/; Revisiting time; 11 min; [54] Radar images were obtained through sampling techniques; for this reason, particular frequency and spatial steps had to be selected in order to avoid ambiguity in range and cross-range [54]. The system is able to measure line-of-sight (LoS) ground displacement in the time interval between two acquisitions and the displacement is calculated from the phase difference between the back-scattered signals received at different times, through the cross-correlation between two SAR images. Range and cross-range resolution are, on average, 2 × 2 m, with a measurement precision being referred to the displacement of less than 1 mm [54]. The displacement rate is the result of the mathematical division between the displacement measured in an interferogram (referred to the difference between two SAR images) and the elapsed time between the two images, allowing for the identification of very low displacement rates (0.010-0.001 mm/h) related to the creep of the northern sector of the SdF or very fast displacement rates (up to 300 mm/h) associated with effusive vent opening [55]. The capability of InSAR to detect ground displacement depends on the persistence of phase coherence over appropriate time intervals, therefore a SAR coherence mask (threshold = 0.5) was set to mask the noisy areas of the interferogram [54]. The phase values can be affected by ambiguity (unwrapped phase), but, due to the short-elapsed time between two subsequent measurements on Stromboli volcano, the interferometric displacements were usually smaller than half wavelength, so that no unwrapping procedures were needed.
Averaged power (amplitude squared) images produced by the GBInSAR NE400 were used to quantitatively evaluate the changes in the NEC's morphology, as they were the closest devices to the NEC. Each analyzed image was derived-from-48-images averaged (≈ 1 every 4 h) in order to increase the signal to noise ratio and, in doing so, facilitate its interpretation. Because the NEC rim produced a shadow zone corresponding to the crater depression, it was possible to calculate the area of the crater itself as it changed over time [14]. The standard deviation was calculated as equal to 80 m 2 (see [14]). A threshold was set at~50 dB to map the area affected by the NEC widening and narrowing collapse.

Satellite Remote Sensing Monitoring
Multispectral satellite images processing is increasingly demonstrating its potential in providing both timely event detection for volcanic effusive events and, in the case of eruption, extraction of parameters that can help tracking the lava flow [59]. Even if the volcanic features of interest are usually much smaller than the nominal pixel size of the satellite image, moderate spatial resolution sensors (~1 km) can detect emitted radiance in the mid-infrared (MIR) wavelengths, a spectral region in which high temperature events, such as active lava flows, vents, and domes, emit copious amounts of energy. Satellite infrared data represented a useful means to describe the evolution of the eruptive event occurred during 28 March-1 April at Stromboli volcano. In particular, we used the HOTSAT system [60][61][62] to detect the presence of thermal anomalies through the analysis of multispectral infrared images acquired by a variety of satellite sensors with a revisit time of about 12 h per satellite, such as VIIRS (Visible Infrared Imaging Radiometer Suite), providing at-nadir pixel footprint of 375 m for I-bands, SLSTR (Sea and Land Surface Temperature Radiometer) on board of SENTINEL-3, and MODIS (Moderate Resolution Imaging Spectroradiometer), both providing 1 km pixels at-nadir. The combined use of sensors that differed for spatial characteristics (from 375 m to 1 km) and different acquisition times has proved to be a robust and reliable instrument for the thermal monitoring of active volcanoes [63][64][65]. The HOTSAT system locates the thermal anomalies (hotspot), computes the associated radiant heat flux summing up the contribute of each hotspot pixel, and, in the case of effusive eruption, provides the Time Averaged Discharge Rate (TADR) as proportional to the radiant heat flux [66]. The conversion from radiant heat flux to TADR was performed according to Harris et al. [67] using: TADR = Q/(ρ (c p ∆T + c L ∆Φ)), where Q is the total thermal flux obtained summing up the radiative power computed for each hotspot pixel, ρ is the lava density (2600 kg m -3 ), c p is the specific heat capacity (1150 J kg -1 K -1 ), ∆T is the eruption temperature minus temperature at which flow stops (100-200 K), c L is the latent heat of crystallization (3.5 × 10 5 J kg -1 , and ∆Φ is the volume percent of crystals that form while cooling through ∆T (30-54%).
The HOTSAT system was extended with a new module to process data acquired by Landsat 8 OLI and TIRS, Sentinel 2 MSI, and ASTER images in order to exploit higher spatial resolution multispectral images. Besides providing further information on the radiant heat flux, these data can be used to locate eruptive vents and describe the evolution of the lava flow field [68][69][70], based on the spatial and spectral resolution of the available bands and the phase and size of the eruption they catch. For example, in the case of an ongoing eruption, Sentinel-2 MSI, thanks to its bands in the SWIR (bands 11 and 12), can provide the position of an active vent and flow at the spatial resolution of 20 m, whereas Landsat-8 OLI, with its SWIR bands 6 and 7, can provide the same information at 30 m of spatial resolution. If a crusted lava flow is cooling, SWIR bands might not be able to detect it; on the other hand, it could still be visible in the thermal infrared TIR bands, i.e., Landsat 8 TIRS (bands 10 and 11) or ASTER (bands 10-14) at 100 and 90 m of spatial resolution, respectively. Due to the limited temporal resolution of these higher spatial resolution multispectral images, post-eruptive images occur more often than intra-eruptive ones and, in many cases, the flows cool too fast to be visible, even in the TIR bands.
Recently, the high spatial resolution and freely available information coming from the Multispectral Imager (MSI) on-board Sentinel-2 satellite has been used to facilitate the two-dimensional (2D) mapping of lava flows [71] through a new Machine Learning (ML) classifier, which discriminates the recent lava flows from pre-and post-eruptive multispectral images acquired by MSI, combined with pre-eruptive digital topography. Bands 2, 3, 4, and 8 at the spatial resolution of 10 m are used as input to the classifier. This ML approach relies on two steps: (i) a k-medoids unsupervised classifier separating input data in clusters whose pixels have the same properties; and, (ii) a Bayesian neural network mapping recent lava flows. In particular, the first step reveals pixels undergoing similar changes in time between preand post-eruptive images, adopting the correlation distance as a measure of similarity. Subsequently, in the second step, a small representative subset of each cluster is exploited to train the BNN, so that it provides us with the pixels belonging to the recent lava flow.
The advancement of satellite remote sensing techniques also has great potential for what concerns the three-dimensional mapping of volcanic products. Indeed, high spatial resolution data acquired in stereo, tri-stereo, or multi-view configuration (e.g., Pléiades, PlanetScope, ASTER) can be used to frequently update the topography and to estimate volcanic deposits by differencing successive topographies. Such estimates can improve the 2D mapping of lava flows while providing an independent maximum bound to lava flow volume that can be derived from the satellite infrared data [59,64,72,73]. We exploited the Pléiades constellation, which is composed of two optical satellites, Pléiades 1A and 1B, respectively, launched on December 2011 and 2012, in order to retrieve areas, volumes, and thickness distribution of the recent volcanic deposits in Stromboli. These satellites provide images at 50 cm spatial resolution in stereo and tri-stereo mode [74]. The 3D processing of the Pléiades imagery was performed using the free and open source MicMac photogrammetric library (available at http://micmac.ensg.eu), in this way 1-m Digital Elevation Models (DEMs) were obtained.
We derived two 1-m DEMs to further constrain the volume of volcanic deposits: a pre-eruptive one from the tri-stereo optical imagery acquired on 8 October 2019, and a post-eruptive one from the Pléaides-1 images of 7 April 2020. We differentiated them, so to obtain the thickness distributions of volcanic deposits emplaced between October 2019 and April 2020. The two DEMs were first co-registered using the Nuth and Kääb [75] method in order to avoid any errors that could derive from a misalignment.  Stromboli seismicity is typically characterized by explosion-generated signals [77][78][79][80][81] and by persistent volcanic tremor [82]. Seismic signals that are associated with landslides [83,84], rolling blocks, PDCs, and lava flows were also recorded in the period 28 March-1 April 2020 ( The signals that are caused by the typical landslides in loose pyroclastic deposits on the SdF ( Figure 2c) have a fusiform envelope and are characterized by a relatively high frequency (4-15 Hz). An increase in the occurrence of these signals was a short-term precursor of the flank effusive eruption on 2007; therefore, they have been studied in detail in various works, especially from the point of view of their automatic and early detection [83,84]. Such signals are partly due to the morphogenetic processes of the SdF slope. Moreover, they can be related to the explosive activity that throws incoherent materials on the slope, which can be easily re-mobilized, thus generating landslides. These are also favored by effusive activity, as the lava front is often a source of incoherent material, which can move on the steep slope of the volcano's flank [84]. Additionally, the signals that are caused by the rolling of large blocks (Figure 2d), in some cases, can be associated with detachments occurring at the lava flow front. They are similar to the ones associated with landslides in loose pyroclastic deposits, however some impulsive phases can be recognized in the waveform, and the spectrogram shows a less gradual onset when compared to that of the landslides in loose pyroclastic deposits. In the example of Figure 2d, the association of the seismic signal with the rolling of large blocks was based on its comparison with the images of the cameras. But an experienced seismic analyst is able to distinguish these signals, even on the basis of the visual analysis of the seismogram alone. This type of signal can also originate from rock-falls along the cliffs of the Labronzo area (North edge of the SdF). Occasionally, they were recorded on the southern side of the island, where there are steep morphologies. The signals caused by the rolling of large blocks (Figure 2d) are similar to the ones that were associated with landslides, however some impulsive phases can be recognized in the waveform, and the spectrogram shows a less gradual onset as compared to that of the landslides. These signals can also be associated with detachments occurring at the lava flow front. The signals that are linked to hot avalanches or PDCs (Figure 2e) are characterized by frequencies with a range 1-5 Hz and are due to incoherent hot materials that form massive flows. On 31 March 2020, the PDC mechanism generated seismic signals with large amplitude at the stations that are closest to the SdF slope.

STRG
Guralp CMG40T 50 IST3 Nanometrics Trillium120PA 100 Stromboli seismicity is typically characterized by explosion-generated signals [77][78][79][80][81] and by persistent volcanic tremor [82]. Seismic signals that are associated with landslides [83,84], rolling blocks, PDCs, and lava flows were also recorded in the period 28 March-1 April 2020 ( Figure 2). Although the explosion signals are transient events, compared to the recordings of tectonic and volcano-tectonic earthquakes, they ( Figure 2a) are characterized by an emerging onset, so that a clear arrival of the P wave cannot be recognized on the seismogram. Moreover, they show a wide frequency band (0.05-10 Hz) containing a Very Long Period (VLP) event. The volcanic tremor ( Figure  2b), typical of open conduit volcanoes, has a frequency content of between 1 and 3 Hz.

Eruptive Activity between 28 March and 1 April 2020
The eruptive activity taking place at Stromboli between 28 March and 1 April 2020 was studied using the images recorded both by satellites and by the INGV-OE monitoring cameras network. On 30 March and for the first half of the day on 1 April, poor weather conditions limited the visibility of the crater area.
The explosive activity at the summit vents was rather intense on 28 March, with 20-25 explosions/h obtained from both cameras SPCT (West flank) and SCT (East flank of the SdF, Figures 1 and 3a). The eruptive activity took place at the NEC, featuring very intense explosions generating a spherical shape of incandescent spatters, which spread more laterally than vertically and all along the crater rim, fell on the NE outer flank and rolled down the SdF slope. This relationship between magma depth and column shape has been established previously [85], thus this shape of the explosions suggested that the level of magma within the vent was very shallow-estimated in a few tens of meters. Meanwhile, the explosions produced by the SWC were essentially of hot gas and ash with collimated jets extending more vertically than laterally, thus indicating that the level of magma within this vent was rather deep-estimated in a few hundred of meters [85]. At 15:38 the first landslide of hot debris coming from the NE rim of the NEC was observed, and it was followed by a powerful explosion. A similar landslide-the second reported on that day -occurred at 15:42, and it was apparently caused by the instability of the hot debris accumulated by the explosions on the steep slope of the NEC outer crater Remote Sens. 2020, 12, 3010 9 of 26 flank. Between 16:09 and 17:00 a series of small landslides from the NEC crater rim began, whose frequency increased in time until it became almost continuous. At 17:02 a lava flow started overflowing from a vent located at the base of the NE crater rim, and both thermal cameras suddenly recorded a significant decrease of the total number of explosions/h, from the previous 20-25 to 5-15 explosions/h (Figure 3a). At 17:44, the crater outline, as observed from SQV, showed a v-shaped cut in the crater rim, located above the effusive vent, which was caused by a collapse. The area of the missing block, estimated from the images recorded by the SQV camera, was~110 m 2 . When considering a thickness of~10 m for the missing block, the volume of the crater rim eroded by the landslides resulted iñ 1.1 × 10 3 m 3 . The erosion caused by the spreading of the lava flows on the North flank of the cone, forming a channel that was also widening by failures and erosion of the lateral margins, was~1000 m 2 , estimated again by the images of the SQV camera. When considering a thickness of~5 m for the collapsed area, a total volume of~5.0 × 10 3 m 3 can be estimated. At 18:00, several incandescent blocks detaching from the lava flow fronts reached the sea, where formed an apron. After 19:00, the lava flow gradually decreased its output rate, and so did the number of landslides along the SdF, that had been triggered by the failure of incandescent blocks detaching from the flow fronts ( Figure 3b). By midnight, the lava flow was no longer fed and, as soon as it stopped, the number of explosions per h gradually increased during the 29 March (Figure 3a). Poor visibility characterized most of 30 March between 07:00 and midnight, but from the images of SCT at 23:31 we could distinguish incandescent blocks that rolled over the NEC crater rim and down the SdF slope, signaling the start of another overflow. The number and size of the incandescent blocks rolling down the slope increased at 23:49, while two lava flows were spreading, one towards East and another towards North. The North flow was the longest and turned out to have the greater flux. The incandescent blocks detaching from the flow fronts accumulated along the coast and formed a hot talus that could be easily viewed during the early h of 31 March. At 01:46 and 01:49, two PDCs descended the SdF reaching the sea, and expanded as a cloud on the sea surface, followed by several other similar events, most of which are listed in Table 3. It is worth noting that the speed of the PDCs spreading on the sea surface was generally increasing from 6.9 m s −1 to 23.3 m s −1 between 01:50 and 02:51, and decreasing afterwards to 5.9 m s −1 until 03:41. The measured distance that was travelled by the PDC on the sea surface from the coast varied between 108 and 165 m, with speeds between 4.4 and 23.3 m s −1 (Table 3).       Table 4 also shows how the path to the coast decreased with time with the extension of the lava flows down the slope. This happens because many PDCs were starting from the lava flow fronts by detachment of hot blocks at breaks in slope. Once the PDCs reached the coast and spread on the sea surface, most of them formed ash clouds that expanded backwards and upslope to the crater area. Not all of the lava flows and PDCs that actually occurred are reported in Table 4, but only those where the visibility was clear enough to allow for an accurate measurement of the path and speed.
The alternation between lava flows and PDCs that descended the SdF slope indicates the gradual erosion of the summit cone that is caused by the emplacement of the lava flow and by the erosion of the summit cone. After 04:44 the lava flow widened at the coastline, forming a hot apron. Lava flows and PDCs continued during the morning at a decreasing rate corresponding to the gradual decrease of the supply to the lava flows, accompanied by a gradual increase of the explosions number from the summit craters (Figure 3a). The lava flow output decreased even more after 19:00 and, by 22:30 of 31 March, the lava flow was apparently no longer fed. Figure 3 shows a comparison between different parameters measured during the period of interest (28 March-1 April 2020). Figure 3a shows how the number of explosions per h before this second lava flow was~25 explosions/h, it decreased to below 5 explosions/h with the start of the lava flow output, and increased again after the end of the lava flows. This observation is consistent with the erupted volume for the 30-31 March lava flow, which was much greater when compared to that of the 28 March lava flow. After the emplacement of the 30-31 March lava flow, the NE flank of the summit cone outline was significantly modified. Being concave upwards at first, it appeared convex at the end, and the eroded surface was estimated at~730 m 2 . When considering a depth of~10 m, the eroded volume of the summit cone can be estimated to be~7.3 × 10 3 m 3 . This brings the total volume of the summit cone, collapsed from the uppermost NE flank between 28 and 31 March 2020, to 13.4 × 10 3 m 3 . Figure 3b displays the RSAM [86] of the STRE station ( Figure 1) seismic signal (East-West component) filtered at a frequency >10 Hz, calculated within 15-s windows. RSAM stands for Real Time Seismic Amplitude Monitoring and it is based on the moving average of the seismic signal absolute value, optionally filtered in specific frequency bands. This parameter is sensitive to landslides, which generate frequencies >10 Hz in the seismic wave field of Stromboli that is generally dominated by frequencies <10 Hz. These signals can also be associated with lava overflows as the collapsing lava flow fronts can generate landslides. Therefore, the RSAM shown in Figure 3b clearly highlights both of the lava overflows that occurred on 28 and the 30-31 March 2020 as well as the PDCs spreading along the SdF.

GBInSAR Data
The eruptive activity reported between 28 March and 1 April 2020 was not accompanied by a long-term inflation of the summit area (Figure 3c,d and Figure 4), as was recorded for the 2012-2013 and 2014 eruptive activities [14,40,44,58]. Displacements that were recorded by the GBInSAR devices were located around the NEC and were mainly associated to the accumulation and instability of the newly emplaced material (Figure 3a). Abrupt change in deformation behavior (that is, movement away from the sensors) occurred between 01:50 on 31 March 2020 and 05:34 on 1 April 2020 (Figure 3c). Deflation was restricted to the very upper part (Figure 4b

Satellite-Derived Lava Flow Field Retrievals
During the effusive phase that occurred between 30 March and 1 April, we converted the radiant heat flux (Figure 3e) into Time-Averaged Discharge Rate (TADR), which is an estimation of the effusion rate averaged over a certain duration ([87]; Figure 5). Integrating the TADR curve, we obtained an upper and lower bound for the erupted Dense Rock Equivalent (DRE) lava volume that can be placed between 37 and 69 × 10 3 m 3 . This compares to the volume of the NEC eroded during the overflows, which was estimated at 13.4 × 10 3 m 3 .
The mapping of the lava overflow that occurred between 28 March and 1 April was performed through the ML classifier [72], using the Sentinel-2 MSI image acquired on 13 March as representative of the pre-eruptive, and the two images of 7 and 12 April 2020 as post-eruptive. As pre-eruptive topography, we used a 1-m Digital Elevation Model (DEM) that was generated by very high-resolution tri-stereo optical imagery acquired by the Pléiades-1 satellite constellation on 8 October 2019 [88]. Following the same steps that are discussed in [72], pixels with similar spectral properties were grouped in 50 clusters by the k-medoids unsupervised classifier. Subsequently, three pixels for each cluster were labelled to train the BNN, i.e., 150 pixels were used overall. To improve the accuracy, the lava flow field was then refined using the Pléiades image acquired on 7 April 2020, which provides a pixel resolution of 0.5 m. Figure 6 shows the lava flow map resulting from the ML classifier. The area measures 94,500 ± 3380 m 2 . The uncertainty was calculated by multiplying the satellite-derived perimeter (6760 m) by the pixel resolution of the Péiades-1 imagery (0.50 m) [69].
During the effusive phase that occurred between 30 March and 1 April, we converted the radiant heat flux (Figure 3e) into Time-Averaged Discharge Rate (TADR), which is an estimation of the effusion rate averaged over a certain duration ([87]; Figure 5). Integrating the TADR curve, we obtained an upper and lower bound for the erupted Dense Rock Equivalent (DRE) lava volume that can be placed between 37 and 69  10 3 m 3 . This compares to the volume of the NEC eroded during the overflows, which was estimated at 13.4  10 3 m 3 . The mapping of the lava overflow that occurred between 28 March and 1 April was performed through the ML classifier [72], using the Sentinel-2 MSI image acquired on 13 March as representative of the pre-eruptive, and the two images of 7 and 12 April 2020 as post-eruptive. As pre-eruptive topography, we used a 1-m Digital Elevation Model (DEM) that was generated by very highresolution tri-stereo optical imagery acquired by the Pléiades-1 satellite constellation on 8 October 2019 [88]. Following the same steps that are discussed in [72], pixels with similar spectral properties were grouped in 50 clusters by the k-medoids unsupervised classifier. Subsequently, three pixels for each cluster were labelled to train the BNN, i.e., 150 pixels were used overall. To improve the accuracy, the lava flow field was then refined using the Pléiades image acquired on 7 April 2020, which provides a pixel resolution of 0.5 m. Figure 6 shows the lava flow map resulting from the ML classifier. The area measures 94,500 ± 3,380 m 2 . The uncertainty was calculated by multiplying the satellite-derived perimeter (6,760 m) by the pixel resolution of the Péiades-1 imagery (0.50 m) [69].
In the area that was identified thanks to the ML classification ( Figure 6), we found a thickness of volcanic deposits that goes from -14 m (due to the coastal erosion) to 14 m (in proximity of the NE crater) (Figure 7). The volume of the deposits accumulated near the NE crater amounts to 34,600 ± 9,700 m 3 . For calculating the volume of the main lava flow spreading on the SdF, we analyzed the histogram of the thicknesses (inset in Figure 7), finding a peak to 0.8 m. Being the most frequent value, we assigned it to the pixels where the DEM difference was negative, thus estimating a volume of 144,400 ± 79,000 m 3 . Consequently, the total bulk volume of deposits thus amounts to 179,000 ± 89,000 m 3 . The uncertainty was computed by multiplying the areas by the residual vertical accuracy outside the margins of the deposits, i.e., the standard deviation (~1.7 m) of the DEM difference in the area that was not covered by deposits.  In the area that was identified thanks to the ML classification ( Figure 6), we found a thickness of volcanic deposits that goes from −14 m (due to the coastal erosion) to 14 m (in proximity of the NE crater) (Figure 7). The volume of the deposits accumulated near the NE crater amounts to 34,600 ± 9700 m 3 . For calculating the volume of the main lava flow spreading on the SdF, we analyzed the histogram of the thicknesses (inset in Figure 7), finding a peak to 0.8 m. Being the most frequent value, we assigned it to the pixels where the DEM difference was negative, thus estimating a volume of 144,400 ± 79,000 m 3 . Consequently, the total bulk volume of deposits thus amounts to 179,000 ± 89,000 m 3 . The uncertainty was computed by multiplying the areas by the residual vertical accuracy outside the margins of the deposits, i.e., the standard deviation (~1.7 m) of the DEM difference in the area that was not covered by deposits. Figure 6. Google Earth view of the satellite-derived lava overflows from the NE crater rim spreading along the Sciara del Fuoco. The lava flow field (red contour) has been superimposed from the Pléiades image acquired on 7 April 2020.

Seismicity
The observation of large amplitude signals that were associated with PDCs, in particular those that occurred on March 31, is the main peculiarity of the seismic data recorded during the March-April 2020 eruptive crisis at Stromboli. To gain insight into the nature of these signals, we compared the seismic recordings to the images of thermal cameras (Figure 8a), timed with the same reference system as the seismic network (UTC based on GPS). The signal due to a PDC consists of a landslidetype initial part and a near monochromatic phase (peak frequency around 3 Hz) with much larger amplitude. The spectrogram (Figure 8b) clearly highlights the transition between the two phases that occurs exactly at the time of the impact of the PDC on the coast. The landslide-type signal is shown in red in Figure 8c, whereas the near monochromatic phase is drawn in blue. This observation highlights that the PDCs flowing on the ground generate signals that are similar to those typical of

Seismicity
The observation of large amplitude signals that were associated with PDCs, in particular those that occurred on March 31, is the main peculiarity of the seismic data recorded during the March-April 2020 eruptive crisis at Stromboli. To gain insight into the nature of these signals, we compared the seismic recordings to the images of thermal cameras (Figure 8a), timed with the same reference system as the seismic network (UTC based on GPS). The signal due to a PDC consists of a landslide-type initial part and a near monochromatic phase (peak frequency around 3 Hz) with much larger amplitude. The spectrogram (Figure 8b) clearly highlights the transition between the two phases that occurs exactly at the time of the impact of the PDC on the coast. The landslide-type signal is shown in red in Figure 8c, whereas the near monochromatic phase is drawn in blue. This observation highlights that the PDCs flowing on the ground generate signals that are similar to those typical of landslides moving on the SdF. The material accelerates on the slope and then impacts on the coastline (or on the sea surface). The impact generates the large amplitude 1-3 Hz phase.
We also focused on the seismic amplitude of STRE station, near the SdF where the PDCs flow, and STR1 station, which is relatively far from the SdF. The comparison between different seismograms produced in one h at STRE and STR1 stations (Figure 1), including the major PDCs that occurred on 31 March, is displayed in Figure 9 and indicates the explosion (E) and PDC signals. We calculated the ratio between the amplitude of the E3 explosion and PDC4 signals (top of Figure 9) that were recorded at the two stations (STRE: black; STR1: red). In order to evaluate the seismic amplitude, we calculated the average of the absolute values (RSAM [86]) of 30-s windows of both the explosion and PDC signals, starting from the onsets of the transients recorded by the two different stations (vertical component). Subsequently, we calculated the amplitude ratios of the explosion signals and the PDC signals. We obtained E3 ratio (STRE RSAM)/(STR1 RSAM) ≈ 2 and PDC4 ratio (STRE RSAM)/(STR1 RSAM) ≈ 4. This observation highlights the rapid decay of the seismic signal amplitude moving away from SdF slope and confirms that the source of the seismic phase with dominant frequency around 3 Hz is on the surface. landslides moving on the SdF. The material accelerates on the slope and then impacts on the coastline (or on the sea surface). The impact generates the large amplitude 1-3 Hz phase. We also focused on the seismic amplitude of STRE station, near the SdF where the PDCs flow, and STR1 station, which is relatively far from the SdF. The comparison between different seismograms produced in one h at STRE and STR1 stations (Figure 1), including the major PDCs that occurred on 31 March, is displayed in Figure 9 and indicates the explosion (E) and PDC signals. We calculated the ratio between the amplitude of the E3 explosion and PDC4 signals (top of Figure 9) that were recorded at the two stations (STRE: black; STR1: red). In order to evaluate the seismic amplitude, we calculated the average of the absolute values (RSAM [86]) of 30-s windows of both the explosion and PDC signals, starting from the onsets of the transients recorded by the two different stations (vertical component). Subsequently, we calculated the amplitude ratios of the explosion signals and the PDC signals. We obtained E3 ratio (STRE RSAM)/(STR1 RSAM) ≈ 2 and PDC4 ratio (STRE RSAM)/(STR1 RSAM) ≈ 4. This observation highlights the rapid decay of the seismic signal amplitude moving away from SdF slope and confirms that the source of the seismic phase with dominant frequency around 3 Hz is on the surface.

Discussion
Between 28 March and 1 April 2020, the summit craters of Stromboli volcano produced two overflows, on 28 and 30-31 March. The two episodes lasted ~7 h and ~23 h, respectively, and were accompanied by the descent of PDCs down the SdF slope and on the sea surface. Integrating several monitoring data, comprising visual and thermal images from a network of fixed cameras, ground

Discussion
Between 28 March and 1 April 2020, the summit craters of Stromboli volcano produced two overflows, on 28 and 30-31 March. The two episodes lasted~7 h and~23 h, respectively, and were accompanied by the descent of PDCs down the SdF slope and on the sea surface. Integrating several monitoring data, comprising visual and thermal images from a network of fixed cameras, ground deformation from GBInSAR, seismicity, and satellite images, we could reconstruct the sequence of events that occurred between 28 March and 1 April 2020 and understand their eruptive processes, gaining useful insights for hazard assessment.
The explosive activity, in terms of the number of explosions versus time, increased before the 28 March and 30 March lava flows (Figure 3a), which suggested that the magma level was becoming shallower within the feeder conduit [19,39,49,50] prior to the lava flow output. Lava flows were heralded and followed by rock-falls and landslides, which triggered the descent of PDCs down the SdF slope ( Figure 3b). The number and duration of these events, as recorded by the seismic network, appear much greater during the 30-31 March lava flow, which lasted longer than the previous event. The first effusive episode was not preceded by significant deformation of the summit zone as detected by the GBInSAR (Figure 3c,d and Figure 4), whereas the second was accompanied by a sudden widening and narrowing of the NE crater zone (Figure 3c,d), and followed by a deflation of the summit crater terrace (max 17 mm recorded away from the sensors, Figure 4b). It should be noted that deflation of the crater terrace (Figure 3c,d) did not occur after the 28 March lava flow, confirming that its erupted volume was probably much smaller compared to the 30-31 March event. The number of explosions at the summit craters decreased after both the lava flow outputs, but this decline was greater and lasted longer on 30-31 March (Figure 3a), yet again suggesting a greater drainage of the shallow conduit, consistent with a larger erupted volume, which generated this lava flow. Hence, during the eruptive phase that took place between 28 March and 1 April at Stromboli, the shift at the summit craters from the persistent Strombolian explosions to lava flow output was heralded by an increase in the rate of explosions, and by a shallower magma level within the feeder conduit [19,39,49,50]. Conversely, the decrease of explosive activity following the lava flow output suggests that the drainage of the uppermost conduit was efficient, requiring a certain amount of time (of the order of h, Figure 3a) in order to allow the magma level to rise again after drainage in order to restore the persistent Strombolian activity at the summit craters. This is consistent with similar events, which were observed during the much longer 2002-2003 and 2007 effusive phases [19,39,[89][90][91]. Satellite data allowed us to obtain the map of the lava flows expanding on the SdF slope between 28 March and 1 April, which was estimated at 94,500 ± 3380 m 2 (Figure 6), the radiant heat flux over time (Figure 3e), and, consequently, the TADR ( Figure 5), providing an estimation of the cumulative erupted lava flow volume at 37-69 × 10 3 m 3 DRE. This compares to the volume of the NEC flank eroded during the overflows that was obtained from the camera images, resulting in 13.4 × 10 3 m 3 .
The behavior of the lava flow spreading from the craters down to the steep SdF slope was consistent with an increase of the lava viscosity and of the yield strength, caused by a decrease in gas-content and bulk temperature of the flow, as well as its crystallization [92]. As a consequence, the lava flow front fragmented and hot blocks detached and descended the slope, the fragmentation generating several PDCs [19,39,44]. Again, PDCs formation was more relevant during the 30-31 March lava flow than during the previous 28 March lava flow, confirming a lower volume and/or shorter extent of the former lava flow. Poor weather conditions impeded a clear view of the summit crater area from the camera monitoring network, which only detected the mid-lower portion of the slope along which the lava flowed and PDCs were spreading. However, we could obtain the average spreading velocity of several lava flows and PDCs descending down the SdF slope (Table 4) and on the sea surface (Table 3). These data show an increase of the speed of the PDCs spreading on the sea surface from 6.9 m s −1 to 23.3 m s −1 between 01:50 and 02:51 on 31 March, which later decreased to 5.9 m s −1 until 03:41. This is consistent with the greater speed of PDCs and lava flows along the SdF slope that was detected during the same lapse of time (Table 4), which led us to assume that the greater PDC velocity on the sea surface was probably caused by an initial greater thermal efficiency of the PDC starting from a higher elevation and flowing along the slope above the active, well-fed and hot lava flow. This might have increased PDC mobility, allowing for it to reach the coast and spread on the sea surface at high speed. By comparison, the velocity of some PDCs emplaced in 1997 at Montserrat by dome collapses whose velocity was between 8 and 21.9 m s −1 [93,94], whereas much greater values of~100 m s −1 were obtained for the PDC that was emplaced during the 1888 phreatic eruption at Bandai Volcano [95]. Several PDCs on 30 and 31 March at Stromboli spread on the sea surface up to a maximum distance of~165 m from the coast, with an estimated speed up to~23 m s −1 (Table 3). This velocity compares pretty well to the measurement of more than 27. propagating for several tens of meters from the coast at speeds of 5.9 and 9.8 m s −1 [33].
The velocity to which the PDC expanded on the sea surface decreased along with the supply to the lava flow that was descending down the slope; at the same time, the PDC was travelling a shorter distance along the SdF slope (Table 4), thus causing a slower expansion of the PDC along the slope and reaching the coast at a lower speed. The measured distance that was travelled by the PDC on the sea surface from the coast varied between 108 and 165 m, with speeds between 4.4 and 23.3 m s −1 (Table 3). Hence, we understand that also small volume PDCs, like those produced during 28 March-1 April 2020 at Stromboli, can spread on the sea surface for hundred meters distance from the coast, causing a potential hazard for bathers, fishermen, or touristic boats sailing along the North coast of the island. Moreover, it must be considered that the ash cloud spreading backwards from the sea to the craters might produce ash fallout that could possibly reach tourists trekking along the SdF margins, even at lower heights.
PDCs spreading on the sea surface are not uncommon on Stromboli, and they surely represent an underestimated threat. They occurred during the subaerial and submarine landslide that caused a tsunami in December 2002 [16,17,96]; during the initial phases of flank effusive eruptions in [2002][2003]2007, and 2014 [32,36,39]; after paroxysmal explosive eruptions triggered by column collapse [39,97,98]; and every time there is a small failure (with volumes of the order of 10 4 or 10 5 m 3 ) of the summit craters due to overloading or instability [33,34]. The two most impressive episodes occurred in 2019 as a result of the 4-5 km high eruptive column collapse after the paroxysmal explosions of 3 July and 28 August. On these occasions, the PDCs expanded on the sea surface for several hundred meters threatening a boat with tourists onboard (https://www.youtube.com/watch?v=RPKgS3sPP1Y). Lava flow velocities are normally much lower when compared to the mobility of PDCs, with speeds normally below 10 m s −1 ([99], and references therein), because of their greater viscosity (e.g., [100]). Only rarely they pose a threat to people, with the exception of few cases of lava accumulated within a crater and suddenly drained by a fissure opening, such as on Nyiragongo [101,102]. However, at Stromboli, lava flows can represent an indirect threat because they can erode the base of the summit cone and trigger summit collapses [14]. Table 4 shows the very different speed along the same SdF slope between lava flows, which had velocities between 2.8 and 5.2 m s −1 , and the PDCs, which displayed velocities between 12.9 and 40.3 m s −1 . Several authors [103][104][105] observed that abrupt slope changes or variations in magma supply affect the velocity of lava flows. When the slope is greater than 24 • , it can result in the detachment of blocks from the flow front to form talus and, when the slope is even greater, as in the case of the SdF on Stromboli that reaches 30-35 • [44,106], the tensional stresses overcome the tensional strength, so that the lava cannot flow any longer. Instead, it breaks into incandescent blocks [103] rolling down the slope and forming a PDC and a distal pile of talus.
The monitoring of PDCs along the slope of steep volcanoes, like Stromboli, is therefore of crucial importance, because they unexpectedly evolve into flank collapses possibly triggering tsunamis, as happened at Stromboli in 2002 [16,17,96], and more recently at Anak Krakatau in 2018. On that occasion, the tsunami caused over 430 fatalities, injured 14,000 people, and displaced 33,000 more along the Sunda Strait [18,107]. Interesting enough, the precursory phase of the Anak Krakatau flank collapse was characterized by an increase of the eruptive activity that lasted for 175 days, and the collapse was preceded by two seismic signals consistent with minor mass movements as well as a momentary quiescence [18,107]. Luckily, the mass movements in the eruptive crisis of 28 March-1 April 2020 at Stromboli were of modest size, not comparable to the events that were described at Anak Krakatau. However, the PDC seismic recordings of Stromboli share some characteristics with the PDCs recorded on other volcanoes, such as Merapi [108] and Unzen [109]. In general, the PDC seismic signals have frequency content between 2 and 15 Hz, similar to that of signals due to landslides in loose clastic material. When PDCs originate from dome collapse, as in the case of Unzen [109], they may contain a lower frequency component caused by the collapse of the dome. In general, the onset of PDC seismic signals emerges with a gradual increase in amplitude, which then remains constant on average for the entire time of the event. These characteristics lead to the classification of PDC seismic signals as continuous and non-transient signals [110]. In addition to these characteristics, in the Stromboli PDC seismograms we were able for the first time to identify the seismic signals caused by the impact of PDCs on the coastline. They could be due to the development of a T phase that can be generated when a landslide enters the underwater environment ( [89] and references therein). However, the very well defined peak frequency around 3 Hz (Figure 2e), within the tremor frequency band (Figure 2b), can also be attributed to the resonance of the shallow conduit located in the volcano edifice at a small depth below the SdF [77], which generates the volcanic tremor of Stromboli [82].
Processing and combining multispectral infrared images that were acquired by a variety of satellite sensors with different spatial characteristics and acquisition times allowed for us to derive the radiant heat flux from 28 March at 01:10 to 1 April at 20:30, finding a peak of thermal activity of~1. Using a ML approach, we were able to estimate a cumulated area over which the two lava overflows emplaced, which amounts to 94,500 m 2 . From DEM difference, we retrieved a total average bulk volume of~179 × 10 3 m 3 emplaced between 8 October 2019 and 7 April 2020, which reduces to a DRE volume of~136 × 10 3 m 3 when considering an average lava vesicularity of 25% [111]. This value includes four overflows, which all occurred in the same area, on 18 January, 3 February, 28 March and 30-31 March 2020. Comparing it to the satellite-derived estimate for 30-31 March, we found a DRE cumulative volume of~83 × 10 3 m 3 emitted during the three previous events. Due to the comparable duration (a few hours) of the overflows occurred on 18 January, 3 February, and 28 March, it is plausible to divide this volume equally, obtaining~27.5 × 10 3 m 3 per each eruptive episode. Summing up this value with the one derived from multispectral infrared satellite images for the 30-31 March, a lava volume of 80.5 × 10 3 m 3 could have been emitted during 28 March-1 April, eventually providing an average thickness of 0.85 m.

Conclusions
Even if eruption-induced mass-flows at Stromboli volcano are common, the triggering mechanisms are yet to be fully understood, also because of the great diversity of the observed phenomena. In this study, the mass-flows that were associated with overflows occurred between 28 March and 1 April 2020 have been analyzed through the use of remote sensing data, both with ground and satellite based sensors, and deriving from seismic sensors. The analysis of the videos recorded by the network of fixed monitoring cameras allowed for the description of the events, as well as the calculation of the velocity to which the lava flowed and the PDCs descended down the SdF and on the sea surface. These videos also made it possible to count the explosions occurring at the summit craters over time. Two GBInSAR devices detected both ground-deformation and morphological changes induced by the slope instability. Remote sensing data also include multispectral satellite data, used to constrain lava flow field area and volume, as well as an estimation of the effusion rate. Moreover, seismic data made it possible to characterize the various stages of the instability phenomena and, at the same time, to integrate the camera data for the description of the relationship between overflows and PDCs, and they were also useful for detecting the signal of the impact of PDCs on the sea surface.
The main results of this study can be summarized, as follows: before the analyzed phase, the explosive activity at the summit vents was reasonably intense (20-25 explosions/h), with a prevalence of explosions that produced coarse material in the NEC (i.e., shallow magma level in the conduit); -the 28 March 2020 overflow was anticipated by some landslides that involved the material accumulated in the areas around the NEC (total eroded volume~5-6 × 10 3 m 3 ), even if these did not generate a substantial widening of the crater itself; -the first overflow was accompanied by a decrease of the total number of explosions/h (from the previous 20-25 to 5-15 explosions/h); -PDCs were also generated by the crumbling of the overflow front, they reached the sea and formed an apron on the coast; -no ground deformation was recorded before nor after the 28 March event, meaning that the lava flow volume was small; -after the first overflow, the number of landslides detected with the seismic network decreased, while the number of explosions increased again, suggesting a new upward movement of the magma level within the conduit; -the onset of the new overflow phase occurred on 30 March together with a new sharp reduction in the number of explosions, a new increase in the number of landslides, which produced a significant variation in the morphology of the crater and which were associated with the accumulation of incandescent material along the coast line; -the PDCs linked to the initial phase originated from the NEC area (total eroded volumẽ 7.3 × 10 3 m 3 ), whereas, as the effusive phase progressed, the subsequent PDCs were generated directly by crumbling of lava flow front along the steep slope of the SdF; -PDCs reached the sea with variable speed (between 12.9 and 40.3 ms −1 ), partly flowing on the water; -the entry into the sea of these mass-flows is associated with a strong variation in seismic signals, with the disappearance of the typical signal associated with the landslides in Stromboli (high frequency; 4-15 Hz) and the appearance of another one characterized by a large amplitude and lower frequency (1-3 Hz); -this change in the seismic signal could be due to the PDC entrance in the underwater environment, as well as to the resonance of the Stromboli conduit, which is located in the volcano edifice, at a small depth below the SdF; -the lava overflows that were emplaced between 28 March and 1 April covered a total area of 94,500 ± 3380 m 2 ; -the volume of the deposits accumulated from October 2019 to April 2020 near the NE crater amounts to 34,600 ± 9700 m 3 , whereas the volume in the overflows area was of 144,400 ± 79,000 m 3 , for a total amount of 179,000 ± 89,000 m 3 . Thermal satellite data also allowed for constraining the DRE lava volume between 37 and 69 × 10 3 m 3 emplaced from 30 March to 1 April 2020; integrating this result with those that were obtained from DEM difference, a lava volume of~80.5 × 10 3 m 3 could have been emitted during 28 March-1 April.