Frequency-magnitude distribution of earthquakes at etna volcano unravels critical stress changes along magma pathways

ABSTRACT The high seismic productivity of volcanic areas provides the chance to investigate the local stress conditions with great resolution, by analysing the slope of the

frequency-magnitude distribution of earthquakes, namely the _b-_value. Here we investigated the seismicity of Mt. Etna between 2005 and 2019, focusing on one of the largest known episodes of

unrest in December 2018, when most of the intruding magma aborted, rather oddly, its ascent inside the volcano. We found a possible stress concentration zone along magma pathways, which may

have inhibited the occurrence of a larger eruption. If the origin of such hypothetical loaded region is related to tectonic forces, one must consider the possibility that geodynamic

processes can locally result in such rapid crustal strain as to perturb the release of magma. Strong _b-_value time-variations occurred a few days before the unrest event, suggesting new

possibilities for investigating the volcano state and impending eruptions. SIMILAR CONTENT BEING VIEWED BY OTHERS INCREMENT IN THE VOLCANIC UNREST AND NUMBER OF ERUPTIONS AFTER THE 2012

LARGE EARTHQUAKES SEQUENCE IN CENTRAL AMERICA Article Open access 17 November 2021 MONTHS-LONG SEISMICITY TRANSIENTS PRECEDING THE 2023 _M_W 7.8 KAHRAMANMARAŞ EARTHQUAKE, TÜRKIYE Article

Open access 28 November 2023 STATISTICS OF SEISMICITY TO INVESTIGATE THE CAMPI FLEGREI CALDERA UNREST Article Open access 30 March 2021 INTRODUCTION Measuring changes in crustal stress

associated with the ascent of magmatic bodies may be crucial in understanding a volcano’s state and could make the assessment of incoming volcanic eruptions more effective. On 24 December

2018, Mt. Etna (eastern Sicily, Italy) experienced one of the strongest known unrest episodes, related to the intrusion of about 30 × 106 m3 of magma1 inside the upper crust below the

volcano, as evidenced by different volcanological and geophysical information2,3,4,5,6. Nevertheless, this major unrest event resulted in a relatively modest, short-lived eruption. In

investigating such anomalous behaviour of the magma supply system, we tested using a crustal stress-meter, i.e., the slope _b_ of the _Gutenberg & Richter_ Frequency-Magnitude

Distribution7 (FMD; log_N_ = _a−bM_, where _N_ is the cumulative number of seismic events with magnitude above or equal to _M_ and _a_ represents the productivity), otherwise acknowledged as

the _b_-value. It describes the average size distribution of earthquakes and is sensitive to differential stress8,9,10,11,12,13,14 and material properties15,16,17, which in turn are also

conditioned by fluid pressure and thermal gradients18,19. In particular, we investigated the seismicity and space-time variation of the _b_-value, considering about 13,700 volcano-tectonic

earthquakes, which occurred at Mt. Etna between 2005 and 2019 (Fig. 1). This methodology enables investigating large crustal volumes, wherever they are “illuminated” with a sufficient number

of earthquakes. We speculate that stress concentration zones along magma pathways may temporarily inhibit the occurrence of large eruptions. Our study suggests that for better outlining and

understanding plumbing systems and the volcanoes’ state, spatial and temporal analysis of seismicity and _b_-values could be applied. CONSIDERATIONS ON THE PHYSICAL MECHANISMS CONTROLLING

THE _B_-VALUE The average size distribution of earthquakes, described by the slope of the FMD (i.e. the _b_-value) is inversely dependent on the mean magnitude, and thus to some extent on

the mean length of seismogenic structures20,21. However, it has been confirmed many times in both laboratory8,10,11,13 and field studies9,12,14 that it depends not only on the material

properties and fractal dimension of fractures and fault systems15,16,17 but also closely on the local stress conditions. The _b_-value averages around 1 on a global scale20,22, which means a

nearly tenfold increase in the number of earthquakes for successively smaller magnitudes, but local conditions may lead to remarkable anomalies. Pockets of _b_-values considerably greater

than 1, embedded in average valued crust23, are generally observed in volcanic and geothermal areas24, and an almost single school of thought16,19,21,24,25,26,27,28,29,30,31,32,33,34 trace

them to the neighbourhood of active magma storage and conduits, corresponding to volumes repeatedly intruded by magmatic fluids, and thus intensely fractured by mechanical and thermal

stress18,19. Consequently, a high fractal dimension of seismicity characterizes these crustal regions16,17, which means that a dense net of relatively small-scale fractures almost completely

fills the seismogenic volume, with the preferential occurrence of frequent low-magnitude earthquakes (i.e. high _b_-values). Slip along fault surfaces is further favoured by the high fluid

pressure, which has the effect of locally reducing the normal effective stress9,26,27,31,34. Conversely, values of _b_ considerably lower than 1, which result from a statistically reduced

number of small-energy earthquakes, are related to stress concentration zones8,9, as asperities along fault planes14,35, where a relatively large volume is available for faulting9,20.

Temporal variations of the _b_-value must be interpreted in the framework of possible time-changes in the relevant stress conditions10,11,13,14,36. Experimental and numerical results show

that increasing differential stress is accompanied by the closing of microcracks and their propagation and coalescence into larger ones, resulting in higher-magnitude earthquakes and thus in

decreased _b_-values11. A rupture is likely to connect one sub-volume to another in a highly stressed environment because the system contains more energy8,12. Conversely, an increase in the

_b_-value is believed to be the reflection of decreased differential stress, which can derive from the injection of pressurized fluids or the sudden release of seismic energy soon after

large earthquakes14,36, resulting in a temporary higher proportion of small seismic events. ETNA VOLCANO ACTIVITY FROM 2005 TO 2019 From 2005 to the beginning of 2008, the volcanic activity

at Mt. Etna was confined to the summit area (within the Etna summit caldera, Fig. 1), with occasional powerful lava fountains or strombolian eruptions, and associated modest lava flows37

(Fig. 2). Between 13 May 2008 and 6 July 2009, a major effusive eruption took place from vents close to the summit craters38. During the following year and a half, specifically until the end

of 2010, eruptive activity diminished and ~40 × 106 m3 of magma is thought to have accumulated inside the volcano (according to the steady-state balancing law proposed by ref. 2; Fig. 2a).

Summit eruptions vigorously resumed from January 2011 until the first half of 2016, with intense and frequent paroxysmal events that probably ensured the release of the formerly stored

magma2 (Fig. 2a). Except for a few episodes (e.g. on 19 December 2009), the cumulative number of earthquakes and related seismic strain-release (Fig. 2b; seismic energy is computed following

the relation reported in ref. 39) shows a linear trend from 2005 to the first half of 2016, despite the remarkable volumetric gap of erupted magma between mid-2009 and the end of 2010, also

accompanied by intense inflation of the volcano40. Afterward, a major change in the seismic pattern took place at Mt. Etna from the second half of 2016; in fact, exponential growth in the

cumulative number of earthquakes can be observed (black-dashed arrow, Fig. 2b). Most of the earthquakes were located down to 8 km below sea level (all depths are relative to sea level),

(Fig. 2c). Together, a decrease in erupted volumes2 and an almost continuous long-term inflation3 were recorded. During this time, up to the end of 2019, eruptions took place exclusively

from the summit craters, through moderate to strong strombolian activity, associated with small lava flows; a further ~50 × 106 m3 of magma is believed to have accumulated within the supply

system at the end of 20192 (yellow-dashed arrow, Fig. 2a). Between 24 and 27 December 2018, Mt. Etna underwent one of the strongest known unrest phases, evidenced by a marked release of the

seismic strain (Fig. 2b), a fast increase in the volcanic tremor, strong deformation6,41, and intense degassing2. This unprecedented unrest was produced by the intrusion of an estimated ~30 

× 106 m3 of magma, which aborted its ascent at about sea level1. It was only followed by a short-lived eruption (red-dashed line, Fig. 2) that mainly occurred through a fissure opened within

the Etna summit caldera, during which only ~3 × 106 m3 of magma were emitted to the surface1,2,3. Geochemical data42 indicate a long-lasting prelude stage to the eruption begun in 2017 and

involved magma-fluid accumulation at depths >7 km, followed by pressure build-up at 2–5 km deep, 6–7 months before the eruption; Etna was extraordinarily overpressurized before the onset

of the December 2018 event42. RESULTS SEISMICITY AND _B_-VALUES AT ETNA VOLCANO BETWEEN 2005 AND 2019 To investigate seismicity patterns and stress conditions beneath Mt. Etna, with

particular regard to the changes induced by the December 2018 major unrest, we performed a spatial analysis of earthquake foci distribution and related _b_-values over two time-spans; only

the earthquakes that occurred between 2005 and 24 December 2018 (right before the unrest) were included at first, then we considered the entire 2005–2019 catalogue (seismicity and _b_-value

sections and maps in Figs. 3 and 4). The two models were subtracted to highlight relevant changes after the intrusive episode (percent _b_-value difference sections and maps in Figs. 3 and

4). If on the one hand, a certain degree of approximation subsists since the two models are mediated over two long and largely overlapping periods, on the other hand, this approach ensures a

proper number of earthquakes to calculate reliable _b_-values over wide crustal volumes23. Earthquakes were accurately relocated using the “tomoDDPS” software43 and a 3D velocity model

derived by the integration of passive44 and active45 seismic data, revised and optimized for this work (_Patanè, D_., personal communication; see “Methods” for details). The resulting

earthquake foci distribution within the crustal volume beneath Etna volcano is shown in Fig. 3a, b. The _b_-values were estimated in those nodes of a 3-dimensional grid spaced 250 m, for

which a minimum number of events (_N_min) above the magnitude threshold of complete recording (say “magnitude of completeness”46,47, Mc), calculated locally for each of the nodes, was found

over a maximum search radius of 3 km. To make the results less dependent on an arbitrary a priori choice of the _N_min, different _b_ estimations were computed using four _N_min values (i.e.

50, 100, 150, and 200); the mean of such estimates was then assigned to the _b_-value of the considered node (for further details on computing see “Methods”). Given the high degree of

complexity in computing _b_-values and estimating the associated uncertainties, only marked spatial and temporal _b-_differences were ascribed to effective variations in the physical

condition of the crust22,23. Such differences were statistically evaluated using the _p_-test by _Utsu_48, which estimates the probability that two samples may come from the same population;

differences between _b_-values can be considered significant whenever they are different at the 95% confidence level (i.e. _p_ ≤ 5%) or higher23 (see “Methods” for further details). The

following provides insights into different crustal sectors beneath the volcano. MID TO SHALLOW ETNA CRUSTAL SECTOR An alignment of earthquakes, mostly localized above 10 km depth and

trending nearly WNW-ESE, was detected in the south-western lower flank of the edifice; it appears to extend the San Gregorio - Aci Trezza fault (SG-AT, Fig. 1a) to the west. This arrangement

(along the grey-dashed line, Fig. 1b) bounds a region that occupies the southernmost sector of the study area, where seismicity is rather low (Sector A, Fig. 1b); therefore, it was not

possible to calculate _b_-values therein. Some authors49 inferred the existence of a regional dextral shear zone that controlled the emplacement and distribution of the oldest Etnean

volcanic centres, running along this same belt. Seismic occurrence noticeably increases on moving to the north; the earthquakes confined within the first 10 km appear to be clustered along

preferential NNW to NW directions, up to the Mt. Etna summit area (along the green-dashed lines within Sector B, Fig. 1b). The mapped faults follow the same trends (e.g. the Ragalna,

Trecastagni, and the offshore Acireale Lineament continuing landward with the S. Venerina, S. Tecla, and Fiandaca system; R, TC, AL, SV, ST, F, respectively, Fig. 1a). In this sector, values

of _b_ less than or around 1 have generally been found down to almost 1.5 km (Fig. 4a, b, f, g); this may be due to local asperities that sometimes produce relatively large earthquakes when

triggered by magmatic intrusions (e.g. the 27/10/2002 event, with local magnitude, _M_L = 4.7, which produced extensive surface faulting along the SV fault50, or even the moment magnitude,

_M_w = 4.9 earthquake that occurred on 26/12/2018 between the F and ST faults, Fig. 1a; see also refs. 4,41). At 3–4 km of depth, the _b_-values locally are up to over 2 (Fig. 4c, d, h, i),

suggesting a considerable degree of interaction between magmatic fluids and faults also at a distance from the central sector of Etna volcano, due to intense crustal fracturing and

consequent elevated fluid-phase mobility. Deeper than 4 km, a general reduction of _b_ is observed again (Fig. 4e, l). The above-discussed seismic clusters do not continue further north of

the Etna summit caldera, on the sides of which seismicity is generally between about 2 and 8 km depth and appears to be arranged along preferential E-W to NE-SW directions (above the

purple-dashed line, within Sector C, Fig. 1b). On the western side, this seismogenic belt extends from the southern summit area to over 10 km further west (“Western Seismogenic Zone”, WSZ);

on the other, it is shifted to the north, running for almost 15 km, from immediately east of the summit caldera, along the northern sidewall of the “Valle del Bove” to the town of Riposto

(“Eastern Seismogenic Zone”, ESZ). The analysis of this crustal sector highlights the occurrence of a “very high _b-_value zone” beneath the Etna summit, mainly between 1 and 6 km depth,

which likely corresponds to an “intermediate Etna magma storage” (Figs. 3a, c and 4b–e). Considerable high _b_ anomalies at similar depths were found in other volcanoes worldwide and

likewise interpreted24, because the repeated use of the plumbing systems produces heterogeneous, shatter, and fluid-intruded crustal volumes, characterized by a high fractal dimension of

seismicity16,17,19, or rather by a larger proportion of small earthquakes and vice versa. Moreover, this interval is consistent with the depth at which most of the volatiles evolve from

magma51, with the effect of locally increasing fluid pressure, and where open cracks may effectively exist in the host rock24. Within this sector, a narrow volume of “relatively low

_b_-values” (Figs. 3c and 4c, d) stands right beneath the Etna summit caldera, at depths comprised between 2.5 and 4.5 km. Differences in _b_ between the “intermediate Etna magma storage”

and the enclosed volume with “relatively low _b_-values” can be considered statistically significant, since the _p_-value test by _Utsu_48 gives more than 99% confidence (i.e. _p_ ≪ 1%).

After the December 2018 dyke intrusion, a marked decrease in _b_ (down to −70% difference; _p_ ≪ 1%) is observed within the “intermediate Etna magma storage” (Figs. 3b, d and 4g–i; region 2,

in Figs. 3e and4n–p), due to reduced fluid pressure after the eruption; however, at the resolution of our data, the enclosed volume with “relatively low _b-_values” appears to persist to

some extent, since it generally undergoes a much less marked reduction in _b_ (locally around −15%, but also down to −45%; _p_ ≪ 1%), (Figs. 3d and 4h; region 3, Fig. 3e). A considerable

increase (up to 50% difference; _p_ ≪ 1%) is observed along the “central plumbing system” (Fig. 3c, d; region 4, Figs. 3e and4o–q), breaking around the region with “relatively low

_b_-values” (black-dashed arrow, Fig. 3d). Finally, slightly increased _b_-values (up to around 20% difference; _p_ ≈ 5%) are observed up to above sea level, just north-east of the Etna

caldera (region 1b, Fig. 4m, n). As well as temporal differences, spatial differences of _b_-values between the above regions of interest (Figs. 3 and 4) in both the considered time

intervals (i.e. from 2005 to 24 December 2018 and 2005 to 2019) are statistically significant at more than a 99% confidence level (_p_ ≪ 1%). Beneath the summit area, a NNW-trending seismic

cluster (within the red-dashed oval, Fig. 1b), mostly confined within 1 km depth, with some events localized down to 5 km deep, marks the shallower portion of the Etna plumbing system. Low

_b_-values (Figs. 3c, d and 4a, f) are locally related to dyke propagation, which episodically produced earthquakes with _M_L up to 3.9, within the observed period. North of ESZ (Fig. 1b),

very low _b_-values (around 0.5; Fig. 4a, b, f, g) characterize an isolated, W-E-trending seismic cluster extending down to 2 km depth, related to the activity of the Pernicana fault system

(P, Figs. 1a and 4a), where, despite the shallow depths, earthquakes with _M_L above 4 (e.g. _M_L = 4.1, 08/01/2019 event4) sometimes occur. DEEP ETNA CRUSTAL SECTOR AND PLUMBING SYSTEM At

greater depths, a roughly NE-SW striking seismogenic volume steeply sinks to the northwest, reaching 10–20 km deep in correspondence to the “Parmentelli-Ragalna deep cluster”52, (blue-dashed

oval, Figs. 1b and 3a) and north of the “Valle del Bove”, and over 30 km beneath the town of Bronte (Fig. 1a). It is characterized by several thick clusters of earthquakes with _M_L up to

4.8, but most frequently around 1.5, and _b_-values around 1, both before and after the December 2018 intrusion (Fig. 3a, b). Regardless of whether depth is a universal factor controlling

_b_53,54 or not12, the increase of frictional resistance to failure with confining pressure makes qualitative sense; under these conditions fractures tend to grow larger, producing a

systematic lowering of the _b_-value11,12. At a global scale9,53,54, this was found to hold true generally down to nearly 15 km, while at greater depths an increase of _b_ has been observed.

Beneath the north-western slope of Mt. Etna, no considerable increase in _b_-values occurs with increasing depth (Fig. 3a, b), likely due to the presence of an allochthonous seismogenic

crustal slab with good mechanical properties, which keeps the _b_-value low15,19,35 down to 25 km depth. Considering the eastern Sicily collisional geodynamic setting55, such a “NW sinking

seismogenic volume” (Fig. 3a, b) should correspond to the “Hyblean foreland carbonate slab”, nearly cropping out offshore of the town of Acireale49, (Fig. 1a) and steeply subducting

north-westward beneath the volcano. Relatively low _b_-values below 20 km depth may also be favoured by inverse faulting style along NE-trending planes44 since thrust faulting has shown to

globally hold smaller _b_ than that of normal or strike-slip events12,22,35. From 30 km up to around 10 km, a deep “aseismic cusp” (Fig. 3a, b) occurs between the above-discussed seismogenic

zones, narrowing upwards towards the central part of the volcano. Seismic tomographic studies report a low-velocity, tri-axial ellipsoid region extending beneath the entire volcanic area

between 25 and 15 km and passing up to around 10 km to a narrow high p-wave velocity body56. Such an aseismic cusp might be an expression of the partially molten, “deeper Etna plumbing

system”, where deformation mainly occurs aseismically, confined to the north-west by the subducting “Hyblean foreland carbonate slab”. A poorly earthquake-populated volume overlies the deep

“aseismic cusp” in the central-southern portion of the volcano, almost up to sea level; it matches with a well-known “High p-wave Velocity Body” (HVB, Fig. 4f–l), generally interpreted as

high-density cumulates, fractionated by magma during its ascent, stocked and solidified at depth56. _B_-VALUE MONITORING AT ETNA VOLCANO Using seismology to investigate magma movement within

the crust is one of the most intriguing long-term challenges of geosciences. From a monitoring perspective, the _b_-value time analysis (e.g.57,58,59) must be viewed as a proxy of the

average crustal stress variations, and could play an important role in tracking magma; nevertheless, several authors26,35,36,60 confine its temporal variations to a secondary role. This is

largely due to several issues that may derive from artifacts linked to changes in network configurations, reporting, or analysing procedures (which here can be almost neglected; see

“Methods”), and also because the _b_-value is believed to primarily depend on stationary properties of the crust23,26. We point out that spatial _b_ mapping is more resolute and accurate

than its time variations because by taking into account the whole catalogue and sampling constant, statistically significant _N_min above Mc around each close-spaced node of a

three-dimensional grid, distinct _b_-values can be calculated for relatively small adjacent sub-volumes; in this way, anomalously seismogenic zones, even of limited size, can be

distinguished. Conversely, temporal _b_ variations may not be diagnostic if calculated in too narrow regions; as a matter of fact, several earthquakes occurring in neighbouring sectors,

characterized by diverse seismic patterns, must be included for reliable computing, with the effect of introducing non-local trends. Such issues can be bypassed to some extent if _b_ time

series are calculated within large volumes, which overall show relatively uniform characteristics in the spatial _b_ analysis. For monitoring purposes, it is more reasonable to investigate

the crustal sector of the volcano directly involved in the final stages of magma transfer up to the surface. Therefore, we calculated two different time series (TS) of the _b_-value,

averaging the results from selected points within that sector of the shallower portion of the Etna plumbing system showing increased _b_ after 24 December 2018 (region 1b, Fig. 4m, n; TS-1b

in Fig. 5) and the “intermediate Etna magma storage” (region 2, Figs. 3e and 4n–p; TS-2 in Fig. 5). We used an event-by-event sliding window of 200 earthquakes that occurred within a maximum

radius of 3.3 km, with _N_min above Mc (corrected by a factor of +0.247) equal to 50; each calculated _b_ was associated at the time of the last earthquake. The most interesting results

were achieved by computing the percent differences between each calculated _b_ and the former one, then by plotting the cumulated series of such variations (see “Methods” for details).

Before the end of 2016, the _b_-value shows an overall regular trend, except for an abrupt fall on 13 May 2008 (Fig. 5a, b); it relates to stronger than normal earthquakes during the opening

of the dykes that had drained a major lava effusion until 6 July 2009. The more pronounced drop in TS-1b, compared to TS-2, is linked to the propagation of fractures mainly at shallow

depths. During the following year and a half, the remarkable volumetric gap of erupted magma led to a gradual slight rise of the _b_-value in the time series (Fig. 5b). A similar trend is

also observed during the subsequent resumption of summit paroxysmal explosive activity in January 2011, which discontinuously lasted until mid-2016 and is believed to have drained most of

the stored magma2. Such slight rise in _b_ in the time series can be related to a certain increment of fluid pressure within the plumbing system, which had the effect of diminishing the

differential stress, thus favouring the preferential occurrence of frequent, relatively small earthquakes during both the volumetric gap of erupted magma and the subsequent summit paroxysms.

However, the observed variations are very small (on the order of 0.1 percent _b_ differences), which suggests that only minor changes in the stress conditions perturbed seismicity patterns

during Etna volcanic activity. The main change instead occurred after the middle of 2016; the relevant decrease in erupted magma volumes (Fig. 2a) is accompanied by the exponential growth of

the cumulative number of earthquakes (black-dashed arrow, Fig. 2b) and strongly increased _b_-values (Fig. 5a, c). The onset of such remarkable growth in _b_ is observed at the beginning of

2017, corresponding very well with the long-lasting magma-fluid pressure build-up that preluded the December 2018 eruption42. The magmatic volatiles injected within the plumbing system and

the net of fractures cutting through its surroundings were not efficiently released this time, building up a huge fluid pressure inside the volcano, which in turn increases the _b_-values. A

very important event is observed on 01 December 2018 (Fig. 5d), when the _b-_values in TS-2 (i.e. within the “intermediate Etna magma storage”) strongly increases, reaching the apex on the

6th of the month, about 19 days before the 24 December major unrest (more than 1.5% difference in 5 days). This trend is likely the effect of an overpressurization of the plumbing system

between 1 and 6 km depth, leading to the preferential occurrence of many relatively smaller-than-normal earthquakes shortly before the unrest. During the same time, TS-1b (i.e. within the

shallower Etna plumbing system) shows steady values. An abrupt drop of _b_ started on 22 December only in TS-2 but became sharper on the 24th of the month in both TS-1b and TS-2; it marks

the onset of the major unrest event, being related to stronger earthquakes due to the dykes’ propagation, which at first involved the crustal region between 1 and 6 km depth (TS-2), and then

propagated to the shallower sectors (TS-1b), leading to the eruption. A marked increase in TS-2 on 25 March 2019 is observed, followed by a similar trend in TS-1b on 05 April 2019 (Fig. 5a,

c), which may reflect magmatic fluid-pressure build-up within the “intermediate Etna magma storage”, and subsequent transfer to the shallower sectors of the plumbing system. Soon after, in

May 2019, moderate summit eruptions gradually resumed. At the end of August 2019, almost continuous moderate to strong summit strombolian activity, associated with small lava flows, was

observed, but the time series do not show relevant variations. Finally, the marked increase in TS-2 on 30 October 2019, accompanied by almost steady values in TS-1b, is not associated with

any relevant variation in eruptive activity. DISCUSSION The origin of Etna’s magma in the mantle is still the object of debate, especially considering its collisional geodynamic setting55,

atypical for the emplacement of such a very active basaltic volcano. Nevertheless, from the distribution of earthquake foci and related _b_-values, we found that the rigid downgoing “Hyblean

foreland carbonate slab” establishes a mechanical boundary, limiting to the north-west the “deeper Etna plumbing system”, which is highlighted by an aseismic cusp (Fig. 3a). This

arrangement favours the conveyance of magma from depth toward the present-day central part of the volcano and may have also conditioned, over time, the emplacement of a central-type

stratovolcano. Ascending magma directly pushes on the bounding slab, sometimes triggering earthquakes with _M_L > 3 in specific sectors, such as the “Parmentelli-Ragalna deep cluster”52,

(blue-dashed oval, Fig. 3a). Magmatic intrusions are also responsible for high seismicity characterizing the overlying mid-to-shallow Etna crustal volume, as revealed by the available focal

mechanism solutions, which show a radial distribution of P-axes, indicating a pressure source located beneath the summit caldera4,6,61,62. Nevertheless, earthquake locations are not

distributed radially from the volcano centre, but rather organize into distinct, well-demarcated seismic clusters, which can be related to (i) fracturing processes induced by magmatic

intrusions and/or (ii) reactivation of pre-existing fault zones during magma uprising. In this latter case, intruding magmatic bodies would load the bedrock, triggering seismicity along

faults originally delineated under the effect of tectonic forces, unravelling to some extent the regional structural pattern. Beneath the summit caldera, specifically between about 1 and 6 

km depth, very high _b_-values mark the depth at which gas exsolution and vesiculation of the ascending magma induce further recurrent fracturing processes, which over time produced a

permanent fractured and heated crustal region, characterized by a high fractal dimension of seismicity. Here, a dense network of cracks and conduits overall facilitates the establishment of

a persistent “intermediate Etna magma storage” (Figs. 3a, c and 4b–e). If the foregoing is true, then even more so is its central portion expected to show similar characteristics, since it

stands right beneath the volcano summit, in a region that has been repeatedly injected by magma throughout Etna history. Instead, a limited volume with “relatively low _b_-values” has been

detected there (Figs. 3c, d and 4c, d, h). Concerning the origin of the relatively low _b-_value region, one may think that in an open-conduit condition, like that of Etna volcano, magmatic

volatiles may more easily escape along the central axis of the plumbing system, avoiding the building of such a huge fluid pressure to remarkably increase the _b_-value, as instead observed

in the surroundings. A weak point of this interpretation is the narrowness of this relatively low-_b_ volume since it only occurs between 2.5 and 4.5 km depth. Moreover, despite the local

slight variations in _b_ after December 2018 (region 3, Fig. 3e), it appears to persist to some extent, being laterally bypassed by the effects of the intrusive episode on the central

plumbing system (black-dashed arrow, Fig. 3d), while more marked changes of the _b-_value occurred in all the neighbouring sectors (regions 1a, 1b, 2 and 4; Figs. 3e and 4m–q). The

relatively localized stability of _b_-values in time might decouple its origin from the episodic occurrence of stronger earthquakes during dyke propagation, or variation of fluid pressure,

otherwise, such a big unrest phase would also have noticeably perturbed the local stress conditions and the _b-_value. Another possibility is that local high differential stress may have

contributed to the occurrence of such zone with “relatively low _b-_values” within the “intermediate Etna magma storage”. Albeit in a different volcanic context, some authors27 inferred that

the fast ascent of magma that pressurizes the system without degassing can locally increase the applied stress beneath volcanoes, reducing _b_-values. This can hardly be our case, because

the exsolution and diffusion of magmatic volatiles within the hosting rocks are testified by the surrounding “very high _b_-value zone”, in which _b_ diminishes notably after the December

2018 eruption (region 2, Figs. 3e and 4n–p), as a consequence of reduced fluid pressure due to intense degassing2. Furthermore, the role of thermal gradients18,23 in controlling temporal

variations of _b-_values is very difficult to invoke because, aside from heat transport by magmatic volatiles, whose main effect is that of reducing the applied stress and increasing _b_,

thermal conduction in the crust occurs over much larger times27. As a final issue, a tectonic cause is proposed. The volume with “relatively low _b-_values” is located right where the WSZ

and the ESZ are shifted with respect to Etna’s centre (Figs. 1b, 3c and 4c, d, h); in the hypothesis that they represent two _en-echelon_ fault zones, the regional NW-oriented maximum

horizontal stress field (SHmax, Fig. 1b) affecting the Etna region63,64 may induce dextral strike-slip kinematics on them, in turn producing a highly stressed restraining step-over within

the plumbing system. In this way, a growing compressive tectonic strain may have gradually created a physical obstacle that became critical for magma transfer up to the surface, leading to

the continuous accumulation inside the volcano. It may also explain why no eruptive events capable of efficiently emptying the plumbing system occurred from the second half of 2016 up to the

end of 2019, determining a considerable volumetric gap of erupted magma, accompanied by an increase in the number of earthquakes (Fig. 2b, c) and _b-_values in the time series (Fig. 5a, c).

Such relevant changes in the seismicity pattern were not observed instead during the former 2009–2010 gap (Figs. 2b, c and 5a, b); this observation points to a different condition, which

seems not to include any anomalous accumulation of magma inside the volcano. As a matter of fact, this latter gap starts at the end of a large effusive eruption, which may have considerably

depleted the magma reservoir. Under the action of NW-trending horizontal stress, maximum extension is predicted to be roughly ENE-WSW oriented65, in agreement with the NS to NNW direction of

most Mt. Etna dykes1,4,41,56. This is also the trend of the seismicity characterizing the shallower portion of the plumbing system (within the red-dashed oval, Fig. 1b), suggesting that

magmatic intrusions exploit zones of weakness originally shaped by regional tectonics49. According to this latter hypothesis, even though tectonic processes are related to the movement of

huge crustal plates, they might locally result in small-scale deformation that can be even faster than previously expected, so rapid as to severely perturb the isostatic uprising of magma

within the crust. However, this framework is only conceptual, because, as stated above, most focal mechanisms indicate that seismicity above 10 km depth is mainly controlled by magmatic

intrusions4,6,61,62. For this to hold, one must consider the possibility that, aside from the volcano dynamics, which episodically results in relatively strong earthquakes, regional

tectonics continuously control fault kinematics in the background. Regardless of the causative source, we observe a major change in seismicity patterns and Etna behaviour after mid-2016, and

this condition persisted at least up to the end of the observed period, despite the occurrence of one of the largest known episodes of unrest in December 2018. Increased confining stress

within magma pathways may account for such a condition because it can inhibit bubble expansion and then the transition to foamy magma; in this way, the capacity to drive large volumes up to

the surface and produce paroxysmal explosive eruptions (lava fountains) can be temporarily limited. Upward migration of magmatic volatiles is also reduced, which favours their spreading

within the hosting rocks, increasing _b-_values in the time series from the beginning of 2017 (Fig. 5). The arisen condition would result in almost continuous summit activity, like that

observed from August 2019 until the end of the year, fed by a reduced fresh and foamy magma supply from the deeper plumbing system and the re-mobilization of the liquid residue that fills

the shallower sector. After the December 2018 eruption, it should be noted that at shallow depths a moderate rise of the _b_-value is observed again just northeast of the summit caldera

(region 1b, Fig. 4m, n), suggesting a new localized increase in the fluid pressure beneath the volcano. The time series (Fig. 5) show strongly increased _b-_values within the “intermediate

Etna magma storage”, corresponding to an overpressurization of the system starting on 01 December 2018 and culminating on the 06th of the month, about 19 days before the 24 December major

unrest event. This latter was preceded by an abrupt fall of _b_ that started 2 days in advance within the magma storage and then propagated to the shallower sectors of the plumbing system,

tracking dykes’ propagation for depth to the surface. This observation may represent a compelling case that monitoring the _b_-value could be used not only to investigate magma dynamics, its

eruption, and then to characterize the state of volcanoes—as observed at Mt. Etna before, during, and after the 24 December 2018 strong unrest phase—but also for improving eruption

forecasting. Previous studies at Mt. Etna found evidence for increased _b-_values from 1981 to the first half of 198757,58,59, during a long period of intense eruptive activity but often

supplied by gas-depleted magma. Spatial variations of _b_ have been investigated on this volcano between 1990 and 199728 and also from 1999 to 2001, before the July–August 2001 eruption29.

Finally, time and space variations of _b_ have been investigated from August 1999 to December 200530. Nevertheless, a former analysis comparable to that here presented is not available, and

it cannot be ruled out that similar conditions also occurred in the past, leading to major eruptions (e.g. 1992–1993, 2001, 2002–2003). The study of the _b-_value needs to be framed in

broader correlation with other monitoring measurements to better define patterns linked to volcanic activity. Nonetheless, our results indicate that the systematic spatial and temporal

analysis of this parameter may offer an opportunity to investigate the volcano state and improve the assessment of impending volcanic eruptions, this latter being the biggest concern of the

local population and civil protection authorities. METHODS In this study, we analysed the spatial-temporal distribution and the _b_-value of the seismicity recorded at Etna volcano from 2005

to 201966 by a seismic network of over 30 stations, almost all equipped with broadband three-component seismometers, managed by the “_Istituto Nazionale di Geofisica e Vulcanologia,

Osservatorio Etneo”_ (Catania, Italy). ESTABLISHING SPATIAL AND TEMPORAL VARIATION OF THE _B_-VALUE Reliable computation of the _b_-value must account for possible artifacts, i.e., unreal

changes in seismicity patterns, linked to how the observations are bundled and analysed67. These problems can be reduced in case a large number of earthquakes is included in the FMD because

for a small sampling the log-linear fit of the distribution is degraded by the statistics of small numbers68. Unfortunately, further issues on volcanoes derive from the high background

seismic noise (e.g. volcanic tremor), which may cause decreases in the detection of the smaller events, making the earthquake catalogues incomplete and thus leading to deviation from

linearity at low magnitudes46. Furthermore, non-double-couple seismic signals, related to volcanic dynamics (such as explosion quakes or long-period events), may lead to bimodal FMD16,24 and

must be discarded. Accurate seismic monitoring thus is needed for a proper analysis of the _b_-value, especially in volcanic areas. The Mt. Etna seismic network, operating since 1989, was

greatly enhanced in both the number of stations and hardware at the beginning of 2005, making it one of the best seismically monitored volcanoes worldwide. Network improvement produced

remarkable accuracy in hypocentre locations, mitigating the effect of artifacts linked to systematic mislocation and magnitude estimation, and also enhanced the capability to even detect the

less energetic earthquakes, lowering the Mc. Furthermore, a unique magnitude scale (i.e. _M_L) has been employed to estimate earthquake size; this avoids further systematic problems that

may derive from the translation between different scales since the measurement of the _b_-value strictly rests on the accuracy in the determination of magnitudes22,23,69. Finally, continuity

and homogeneity of recording are assured, since no important network modification or failure occurred, and no changes in analysing procedures were made. _M_L strictly rests below 5 and

source depth is shallower than about 30 km throughout the whole seismic catalogue, likely ensuring no breakdown in the self-similarity of earthquakes’ mechanism9,12,16,60. There is a

question whether it is better to use the entire catalogue or its declustered part for mapping _b_-values in volcanic areas, and several authors31 found that declustering led to only minor

differences accompanied by a lack of information in some volumes. Together, as well as other volcanoes around the world, earthquakes at Mt. Etna mainly occur as “swarms”, whereas

foreshock-mainshock-aftershock sequences are rarely observed (e.g. refs. 70,71). Since earthquakes in volcanic swarms are considered as independent, i.e., triggered by mechanical processes

that are not controlled by previous events (e.g. ref. 72), and there are no sequences in the considered dataset, reliable declustering procedures cannot be applied to the Etna seismic

catalogue. Our analysis comprises the following steps. EARTHQUAKE 3D RELOCATIONS About 13,700 earthquakes, with _M_L in the range 0.2–4.8, were relocated using the ‘tomoDDPS’ software43 and

a 3D velocity model derived by the integration of passive44 and active45 seismic data (the last acquired during the “TOMO-ETNA experiment”, FP7 “MED-SUV”73), revised and optimized for this

work (_Patanè, D_., personal communication) to improve accuracy in hypocentre location. The algorithm has the advantage of using a combination of both absolute and differential arrival-time

readings between couples of close-spaced earthquakes, having the same velocity profile along the ray path between the source and a common receiver. The velocity paths being almost identical,

the difference in arrival times is reasonably attributed only to the spatial offset between the seismic events. Earthquake 3D relocations produced a better clustering of the earthquakes and

reduced the residuals between the observed and theoretical arrival times by about 44%, with an average value of 0.02 s. The final locations are affected by uncertainties of less than 0.7 km

in the epicentral coordinates and less than 1 km in depth. MAPPING THE _B_-VALUE The volume of interest was subdivided into a three-dimensional grid spaced 250 m. For each node, the

_b_-value was calculated taking into account the smallest set of closest earthquakes within a maximum search radius of 3 km, such that the number of events with a magnitude above the local

Mc was greater than a fixed quantity (_N_min). If these conditions were not satisfied, a null value was associated with the local _b_. The more events a data sample includes, the more robust

the estimate of the _b_-value23,35,68. Monte Carlo sampling and bootstrapping of real data, as well as analysis of synthetic catalogues, shows that the error in _b_ is large for small

numbers of events and decreases rapidly to ~15% for 50 events and less than 10% for more than 150 events, reaching about 3% for 1000 events (numbers represent 1σ and are of the same order as

the standard error estimate proposed by _Shi and Bolt_74). To mitigate the effect of an arbitrary a priori _N_min choice on the _b_-value, we considered the mean of the estimated _b_,

derived from four different _N_min (i.e. 50, 100, 150, and 200); for the bins where less than four _b-_values are returned, only the available value/values were used. The associated

uncertainties of final _b_ values were calculated by using the error propagation formula. Information about the radius for the different bins using the two end members of _N_min (i.e. 50 and

200) and the mean of the four _N_min values are shown for the W-E vertical section of Fig. 3 in Supplementary Fig. 1 (see “Supplementary Information”). The estimation of reliable _b_-values

depends critically on the choice of the Mc. In particular, the FMD increasingly deviates from a linear power low for earthquakes with magnitude smaller than Mc, since only a fraction of

them can be detected. It follows that underestimating the Mc may result in unrealistic _b_-values. Conversely, if a too high Mc is considered, a significant number of earthquakes that would

fit the linear power-low will be excluded, severely increasing the uncertainties in the estimation of the _b-_values23,35,68 and, at the same time, decreasing the resolution for mapping23.

_Woessner and Wiemer_47 reported that, when considering subsets of earthquakes, one can find locally higher Mc than that calculated for the whole dataset, and they discuss the need for

identifying such subsets when analysing spatial and temporal variations of seismicity parameters. Accordingly, also in line with several previous studies in volcanic areas (e.g. refs.

21,26,31), we considered the local Mc for every minor subset for _b_-values computing. Mc was calculated using the maximum-curvature method46, which consists of defining the point of the

maximum curvature by computing the maximum value of the first derivative of the FMD curve. Mc for the different bins ranges between about 0.7 and 1.9, being 1.1 for the entire catalogue. The

estimation of the _b_-value and its associated uncertainty was performed by using the _Tinti and Mulargia_75 maximum-likelihood formulation, which is specific for binned magnitudes and

bias-free also for a small dataset. More than 95% of the bins, in both the periods of observation (i.e. from 2005 to 24 December 2018, right before the major unrest event, and from 2005 to

2019) have Mc above the bulk Mc calculated on the entire dataset (i.e. 1.1). For this reason and to also avoid excluding the local contribution of small volcano-tectonic earthquakes that

would fit the FMD, we considered it better not to perform any magnitude-based pre-cut to the catalogue, which could be poorly diagnostic for the analysis of single subsets47. Following our

approach, FMD above Mc for the different bins shows maximum goodness-of-fit with a linear power law (estimated following _Wiemer and Wyss_23) of 92.8% and 94.2%, and average values of 62.4%

and 61.7%, respectively for the two above periods. To validate our approach, as a common practice when dealing with _b_-value computation23, we statistically evaluated the spatial _b_

anomalies for each of the two analysed periods (i.e. 2005 to 24 December 2018 and 2005 to 2019) and the temporal differences between them by using the _p_-test by _Utsu_48, which estimates

the probability that two samples may come from the same population, explicitly accounting for the number of earthquakes contained in each sample. In particular, the test starts from the

condition of two earthquakes’ populations to be the same and verifies how true is this condition. The hypothesis of belonging to the same population is more likely if the two sets share a

considerable number of earthquakes, making the difference between two _b-_values less significant. Differences between _b_-values were considered statistically significant whenever they are

different at the 95% confidence level (i.e. _p_ ≤ 5%) or higher23. CONSTRUCTING _B_-VALUE TIME SERIES The time-variant _b_ was calculated using the same above-reported methodologies, but

averaging the results from selected nodes within two crustal volumes of particular interest, each showing relatively uniform characteristics in the spatial analysis of the _b_-value (region

1b, Fig. 4m, n, TS-1b in Fig. 5 and region 2, Figs. 3e and 4n–p, TS-2 in Fig. 5). We used an event-by-event sliding window of 200 earthquakes, occurring within a maximum radius of 3.3 km

from the selected nodes; the radius was slightly expanded to also include the May 2008 unrest event in the time series. A _b_-value was calculated and associated at the time of the last

event only if the _N_min above the Mc of the set, corrected by a factor of +0.2 (as discussed by _Woessner and Wiemer_47), was greater than 50. The time series of the _b-_values, related Mc

and number of events above Mc + 0.2 are shown in Supplementary Fig. 2 (see “Supplementary Information”). In the time series analysis, to highlight the temporal variations of _b_, we

performed a further processing step that consists in computing the percent differences between each calculated _b_ and the former one. To enhance temporal smoothing, the former _b-_value was

averaged with its preceding values (if present) within a 24-h window; the best result was then achieved by plotting the cumulated series of such variations (Fig. 5). DATA AVAILABILITY

Earthquakes’ initial parameters are available at: https://doi.org/10.13127/ETNASC/ETNARCSC CODE AVAILABILITY The analytical codes for _b_-value computing were written with the MatLab

software, by using formulas reported in the literature. The codes are available from Dr. _Flavio Cannavò_ ([email protected]) upon request. REFERENCES * Bonforte, A., Guglielmino, F.

& Puglisi, G. Large dyke intrusion and small eruption: the December 24, 2018 Mt. Etna eruption imaged by Sentinel-1 data. _Terra Nova_ 31, 405–412 (2019). Article  Google Scholar  *

Calvari, S. et al. The VEI 2 Christmas 2018 Etna eruption: a small but intense eruptive event or the starting phase of a larger one? _Remote Sensing_ 12, 905 (2020). Article  Google Scholar

  * Cannavò, F., Sciotto, M., Cannata, A. & Di Grazia, G. An integrated geophysical approach to track magma intrusion: the 2018 Christmas Eve eruption at Mount Etna. _Geophys. Res.

Lett._ 46, 8009–8017 (2019). Article  Google Scholar  * Alparone, S. et al. Seismological constraints on the 2018 Mt. Etna (Italy) flank eruption and implications for the flank dynamics of

the volcano. _Terra Nova_ 00, 1–11 (2020). Google Scholar  * Bonaccorso, A., Currenti, G., Linde, A., Sacks, S. & Sicali, A. Advances in understanding intrusive, explosive and effusive

processes as revealed by the borehole dilatometer network on Mt. Etna Volcano. _Front. Earth Sci._ 7, 357 (2020). Article  Google Scholar  * Mattia, M. et al. Combined seismic and geodetic

analysis before, during, and after the 2018 Mount Etna Eruption. _Geochem. Geophys. Geosyst._ 21, e2020GC009218 (2020). Article  Google Scholar  * Gutenberg, B. & Richter, C. F.

Frequency of earthquakes in California. _Bull. Seismol. Soc. Am._ 34, 185–188 (1944). Article  Google Scholar  * Scholz, C. H. The frequency-magnitude relation of microfracturing in rock and

its relation to earthquakes. _Bull. Seismol. Soc. Am._ 58, 399–415 (1968). Article  Google Scholar  * Wyss, M. Towards a physical understanding of the earthquake frequency distribution.

_Geophys J. R. Astr. Soc._ 31, 341–359 (1973). Article  Google Scholar  * Main, I. G., Meredith, P. G. & Jones, C. A reinterpretation of the precursory seismic _b_-value anomaly from

fracture mechanics. _Geophys. J._ 96, 131–138 (1989). Article  Google Scholar  * Amitrano, D. Brittle-ductile transition and associated seismicity: experimental and numerical studies and

relationship with the _b_ value. _J. Geophys. Res._ 108, 1–15 (2003). Google Scholar  * Schorlemmer, D., Wiemer, S. & Wyss, M. Variations in earthquake-size distribution across different

stress regimes. _Nature_ 437, 539–542 (2005). Article  CAS  Google Scholar  * Goebel, T. H. W., Schorlemmer, D., Becker, T. W., Dresen, G. & Sammis, C. G. Acoustic emissions document

stress changes over many seismic cycles in stick-slip experiments. _Geophys. Res. Lett._ 40, 2049–2054 (2013). Article  Google Scholar  * Gulia, L. & Wiemer, S. Real-time discrimination

of earthquake foreshocks and aftershocks. _Nature_ 574, 193–199 (2019). Article  CAS  Google Scholar  * Mogi, K. Magnitude-frequency relation for elastic shocks accompanying fracture of

various materials and some related problems in earthquakes. _Bull. Earthquake Res. Inst. Univ. Tokyo_ 40, 831–853 (1962). Google Scholar  * Main, I. G. A characteristic earthquake model of

the seismicity preceding the eruption of Mount St. Helens on 18 May 1980. _Phys. Earth Planet. Interiors_ 49, 283–293 (1987). Article  Google Scholar  * Vinciguerra, S., Gresta, S., Barbano,

M. S. & Distefano, G. The two behaviours of Mt. Etna Volcano before and after a large intrusive episode: evidences from _b_ value and fractal dimension of seismicity. _Geophys. Res.

Lett._ 28, 2257–2260 (2001). Article  Google Scholar  * Warren, N. W. & Latham, G. V. An experimental study of thermally induced microfracturing and its relation to volcanic seismicity.

_J. Geophys. Res._ 75, 4455–4464 (1970). Article  Google Scholar  * Wiemer, S. & McNutt, S. R. Variations in the frequency-magnitude distribution with depth in two volcanic areas: Mount

St. Helens, Washington, and Mt. Spurr, Alaska. _Geophys. Res. Lett._ 24, 189–192 (1997). Article  Google Scholar  * Kanamori, H. & Anderson, D. L. Theoretical basis of some empirical

relations in seismology. _Bull. Seismol. Soc. Amer._ 65, 1073–1095 (1975). Google Scholar  * Wyss, M., Shimazaki, K. & Wiemer, S. Mapping active magma chambers by _b_ values beneath the

off-Ito volcano, Japan. _J. Geophys. Res._ 102, 20413–20422 (1997). Article  Google Scholar  * Frohlich, C. & Davis, S. D. Teleseismic _b_ values: or, much ado about 1.0. _J. Geophys.

Res._ 98, 631–644 (1993). Article  Google Scholar  * Wiemer, S. & Wyss, M. Mapping spatial variability of the frequency-magnitude distribution of earthquakes. _Adv. Geophys._ 45, 259–302

(2002). Article  Google Scholar  * McNutt, S. R. Volcanic seismology. _Annu. Rev. Earth Planet. Sci._ 32, 461–491 (2005). Article  Google Scholar  * Power, J. A., Wyss, M. & Latchman,

J. L. Spatial variations in the frequency-magnitude distribution of earthquakes at Soufriere Hills Volcano, Montserrat, West Indies. _Geophys. Res. Lett._ 25, 3653–3656 (1998). Article 

Google Scholar  * Wiemer, S., McNutt, S. R. & Wyss, M. Temporal and three-dimensional spatial analysis of the frequency-magnitude distribution near Long Valley caldera, California.

_Geophys. J. Int._ 134, 409–421 (1998). Article  Google Scholar  * Jolly, A. D. & McNutt, S. R. Seismicity at the volcanoes of Katmai National Park, Alaska; July 1995-December 1997. _J.

Volcanol. Geothermal Res._ 93, 173–190 (1999). Article  CAS  Google Scholar  * Murru, M., Montuori, C., Wyss, M. & Privitera, E. The locations of magma chambers at Mt. Etna, Italy,

mapped by _b_-values. _Geophys. Res. Lett._ 26, 2553–2556 (1999). Article  Google Scholar  * Murru, M., Montuori, C., Console, R. & Lisi, A. Mapping of the _b_ value anomalies beneath

Mt. Etna, Italy, during July-August 2001 lateral eruption. _Geophys. Res. Lett._ 32, L05309 (2005). Article  Google Scholar  * Murru, M., Console, R., Falcone, G., Montuori, C. & Sgroi,

T. Spatial mapping of the b value at Mount Etna, Italy, using earthquake data recorded from 1999 to 2005. _J. Geophys. Res._ 112, B12303 (2007). Article  Google Scholar  * Wyss, M., Klein,

F., Nagamine, K. & Wiemer, S. Anomalously high b-values in the South Flank of Kilauea volcano, Hawaii: evidence for the distribution of magma below Kilauea’s East rift zone. _J.

Volcanol. Geothermal Res._ 106, 23–37 (2001). Article  CAS  Google Scholar  * Sanchez, J. J., McNutt, S. R., Power, J. A. & Wyss, M. Spatial variations in the frequency-magnitude

distribution of earthquakes at Mount Pinatubo volcano. _Bull. Seismol. Soc. Am._ 94, 430–438 (2004). Article  Google Scholar  * Bridges, D. L. & Gao, S. S. Spatial variation of seismic

_b_-values beneath Makushin Volcano, Unalaska Island, Alaska. _Earth Planet. Sci. Lett._ 245, 408–415 (2006). Article  CAS  Google Scholar  * Kato, A. et al. Preparatory and precursory

processes leading up to the 2014 phreatic eruption of Mount Ontake, Japan. _Earth, Planets Space_ 67, 111 (2015). Article  Google Scholar  * Wiemer, S. & Wyss, M. Mapping the

frequency-magnitude distribution in asperities: an improved technique to calculate recurrence times? _J. Geophys. Res._ 102, 15115–15128 (1997). Article  Google Scholar  * Wiemer, S. &

Katsumata, K. Spatial variability of seismicity parameters in aftershock zones. _J. Geophys. Res._ 104, 13135–13151 (1999). Article  Google Scholar  * Proietti, C. et al. Towards a

satellite-based approach to measure eruptive volumes at Mt. Etna using Pleiades datasets. _Bull. Volcanol._ 82, 35 (2020). Article  Google Scholar  * Behncke, B. et al. Lidar surveys reveal

eruptive volumes and rates at Etna, 2007-2010. _Geophys. Res. Lett._ 43, 4270–4278 (2016). Article  Google Scholar  * Richter, C. F. _Elementary Seismology_. (W. H. Freeman & Company,

1958). * Aloisi, M. et al. Imaging the multi‐level magma reservoir at Mt. Etna volcano (Italy). _Geophys. Res. Lett._ 38, L16306 (2011). Article  Google Scholar  * De Novellis, V. et al.

DInSAR analysis and analytical modeling of Mount Etna displacements: the December 2018 volcano-tectonic crisis. _Geophys. Res. Lett._ 46, 5817–5827 (2019). Article  Google Scholar  *

Paonita, P. et al. Intense overpressurization at basaltic open-conduit volcanoes as inferred by geochemical signals: the case of the Mt. Etna December 2018 eruption. _Sci. Adv._ 7, 36

(2021). Article  Google Scholar  * Zhang, H., Thurber, C. & Bedrosian, P. Joint inversion for Vp, Vs, and Vp/Vs at SAFOD, Parkfield, California. _Geochem. Geophys. Geosyst._ 10, Q11002

(2009). Article  Google Scholar  * Scarfi, L. et al. Slab narrowing in the Central Mediterranean: The Calabro‐Ionian subduction zone as imaged by high resolution seismic tomography. _Sci.

Rep._ 8, 5178 (2018). Article  CAS  Google Scholar  * Dìaz-Moreno, A. et al. New insights on Mt. Etna’s crust and relationship with the regional tectonic framework from joint active and

passive P-wave seismic tomography. _Surv. Geophys_. https://doi.org/10.1007/s10712-017-9425-3 (2017). * Wiemer, S. & Wyss, M. Minimum magnitude of completeness in earthquake catalogs:

examples from Alaska, the Western United States, and Japan. _Bull. Seismol. Soc. Am._ 90, 859–869 (2000). Article  Google Scholar  * Woessner, J. & Wiemer, S. Assessing the quality of

earthquake catalogues: estimating the magnitude of completeness and its uncertainty. _Bull. Seismol. Soc. Am._ 95, 684–698 (2005). Article  Google Scholar  * Utsu, T. On seismicity. _Report

of the Joint Research Institute for Statistical Mathematics. Inst. Stat. Math. Tokyo_ 34, 139–157 (1992). Google Scholar  * Firetto Carlino, M. et al. Time and space scattered volcanism of

Mt. Etna driven by strike‐slip tectonics. _Sci. Rep._ 9, 12125 (2019). Article  CAS  Google Scholar  * Azzaro, R., Branca, S., Gwinner, K. & Coltelli, M. The volcano-tectonic map of Etna

volcano, 1:100.000 scale: an integrated approach based on a morphotectonic analysis from high-resolution DEM constrained by geologic, active faulting and seismotectonic data. _Ital. J.

Geosci._ 131, 153–170 (2012). Google Scholar  * Wallace, P. J., Plank, T., Edmonds, M. & Hauri, E. H. Volatiles in Magmas. In _Encyclopedia of Volcanoes_, Second Edition, (ed.

Sigurdsson, H.) 163–183 https://doi.org/10.1016/B978-0-12-385938-9.00007-9 (Elsevier Inc, 2015). * Patanè, D. et al. Insights into magma and fluid transfer at Mount Etna by a multiparametric

approach: a model of the events leading to the 2011 eruptive cycle. _J. Geophys. Res.: Solid Earth_ 118, 3519–3539 (2013). Article  Google Scholar  * Gerstenberger, M., Wiemer, S. &

Giardini, D. A systematic test of the hypothesis that the _b_ value varies with depth in California. _Geophys. Res. Lett._ 28, 57–60 (2001). Article  Google Scholar  * Spada, M., Tormann,

T., Wiemer, S. & Enescu, B. Generic dependence of the frequency-size distribution of earthquakes on depth and its relation to the strength profile of the crust. _Geophys. Res. Lett._ 40,

709–714 (2013). Article  Google Scholar  * Lentini, F., Carbone, S. & Guarnieri, P. Collisional and postcollisional tectonics of the Apenninic-Maghrebian orogen (southern Italy).

_Geological Soc. Am. Special Paper_ 409, 57–81 (2006). Google Scholar  * Patanè, D. et al. Interplay between Tectonics and Mount Etna’s volcanism: insights into the geometry of the plumbing

system. New frontiers in tectonic research at the midst of plate convergence. (ed. Schattner, U.) 73–104 https://doi.org/10.5772/1766 (IntechOpen Limited, 2011). * Gresta, S. & Patanè,

G. Variation of _b_ values before the Etnean eruption of March 1981. _Pure Appl. Geophys._ 121, 287–295 (1983). Article  Google Scholar  * Gresta, S. & Patanè, G. Changes in _b_ values

before the Etnean eruption of March-August 1983. _Pure Appl. Geophys._ 121, 903–912 (1983). Article  Google Scholar  * Patanè, D., Caltabiano, T., Del Pezzo, E. & Gresta, S. Time

variation of _b_ and _Qc_ at Mt. Etna (1981–1987). _Phys. Earth Planet. Int._ 71, 137–140 (1992). Article  Google Scholar  * Urbancic, T. I., Trifu, C. I., Long, J. M. & Toung, R. P.

Space-time correlations of _b_-values with stress release. _Pure Appl. Geophys._ 139, 449–462 (1992). Article  Google Scholar  * Patanè, D. & Privitera, E. Seismicity related to 1989 and

1991-93 Mt. Etna (Italy) eruptions: kinematic constraints by fault plane solution analysis. _J. Volcanol. Geothermal Res._ 109, 77–98 (2001). Article  Google Scholar  * Scarfì, L., Messina,

A. & Cassisi, C. Sicily and southern Calabria focal mechanism database: a valuable tool for local and regional stress-field determination. _Ann. Geophys._ 56, 1–17 (2013). D0109. Google

Scholar  * Palano, M. et al. GPS velocity and strain fields in Sicily and southern Calabria, Italy: updated geodetic constraints on tectonic block interaction in the central Mediterranean.

_J. Geophys. Res._ 117, B07401 (2012). Google Scholar  * Scarfì, L., Langer, H., Messina, A. & Musumeci, C. Tectonic regimes inferred from clustering of focal mechanisms and their

distribution in space: application to the Central Mediterranean Area. _J. Geophys. Res.: Solid Earth_ 126, e2020JB020519 (2021). Article  Google Scholar  * Sylvester, A. G. Strike-slip

faults. _Geol. Soc. Am. Bull._ 100, 1666–1703 (1988). Article  Google Scholar  * Alparone, S. C. et al. Mt. Etna Revised and Concise Seismic Catalog from 1999 (EtnaRCSC) [Data set].

_Istituto Nazionale di Geofisica e Vulcanologia (INGV)_. https://doi.org/10.13127/ETNASC/ETNARCSC (2020). * Habermann, R. E. Man-made changes of seismicity rates. _Bull. Seism. Soc. Am._ 77,

141–159 (1987). Google Scholar  * Howell, B. F. Jr. On the effect of too small a data base on earthquake frequency diagrams. _Bull. Seism. Soc. Am._ 75, 1205–1207 (1985). Article  Google

Scholar  * Zuniga, F. R. & Wyss, M. Inadvertent changes in magnitude reported in earthquake catalogs: their evaluation through _b_-value estimates. _Bull. Seism. Soc. Am._ 85, 1858–1866

(1995). Article  Google Scholar  * Gambino, S., Mostaccio, A., Patanè, D., Scarfì, L. & Ursino, A. High-precision locations of the microseismicity preceding the 2002–2003 Mt. Etna

eruption. _Geophys. Res. Lett._ 31, L18604 (2004). Article  Google Scholar  * Alparone, S. C. et al. Instrumental seismic catalogue of Mt. Etna earthquakes (Sicily, Italy): ten years

(2000-2010) of instrumental recordings. _Ann. Geophys._ 58, 4 (2015). Google Scholar  * Van Stiphout, T., Zhuang, J., & Marsan, D. Seismicity declustering. _Community Online Resource for

Statistical Seismicity Analysis_, https://doi.org/10.5078/corssa-52382934 (2012). * Coltelli, M. et al. The marine activities performed within the TOMO-ETNA experiment. _Ann. Geophys._ 59,

4 (2016). Article  Google Scholar  * Shi, Y. & Bolt, B. A. The standard error of the magnitude-frequency _b_ value. _Bull. Seismol. Soc. Am._ 72, 1677–1687 (1982). Article  Google

Scholar  * Tinti, S. & Mulargia, F. Confidence intervals of _b_ values for grouped magnitudes. _Bull. Seismol. Soc. Am._ 77, 2125–2134 (1987). Google Scholar  Download references

ACKNOWLEDGEMENTS We are grateful to the “gruppo analisi dati sismici” of the “Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo”, for their unceasing work in earthquakes

localization. We are indebted to Boris Behncke, for his precious help in accurately reconstructing the different phases of Etna volcano activity. We thank Steve Conway for revising the

English. AUTHOR INFORMATION AUTHORS AND AFFILIATIONS * Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo, Catania, Italy Marco Firetto Carlino, Luciano Scarfì, Flavio

Cannavò, Graziella Barberi, Domenico Patanè & Mauro Coltelli Authors * Marco Firetto Carlino View author publications You can also search for this author inPubMed Google Scholar *

Luciano Scarfì View author publications You can also search for this author inPubMed Google Scholar * Flavio Cannavò View author publications You can also search for this author inPubMed 

Google Scholar * Graziella Barberi View author publications You can also search for this author inPubMed Google Scholar * Domenico Patanè View author publications You can also search for

this author inPubMed Google Scholar * Mauro Coltelli View author publications You can also search for this author inPubMed Google Scholar CONTRIBUTIONS M.F.C., M.C., L.S. and D.P. conceived

and designed the research. D.P., L.S. and G.B. derived the crustal 3D velocity model and localized the earthquakes. F.C. wrote the codes for the time-space calculation of the _b-_values and

for plotting the results (using the MatLab software). M.F.C., L.S., F.C. and D.P. determined the methodology for _b-_value computing. M.F.C. wrote the text and made the figures (using the

CorelDraw software). M.C. and M.F.C. inferred the influence of tectonics on magma dynamics. CORRESPONDING AUTHOR Correspondence to Marco Firetto Carlino. ETHICS DECLARATIONS COMPETING

INTERESTS The authors declare no competing interests. PEER REVIEW PEER REVIEW INFORMATION _Communications Earth & Environment_ thanks Jamie Farrell, Katie Jacobs and the other,

anonymous, reviewer(s) for their contribution to the peer review of this work. Primary Handling Editors: Luca Dal Zilio, Joe Aslin. ADDITIONAL INFORMATION PUBLISHER’S NOTE Springer Nature

remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. SUPPLEMENTARY INFORMATION SUPPLEMENTARY INFORMATION RIGHTS AND PERMISSIONS OPEN ACCESS

This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as

long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third

party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the

article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright

holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. Reprints and permissions ABOUT THIS ARTICLE CITE THIS ARTICLE Firetto Carlino, M., Scarfì, L.,

Cannavò, F. _et al._ Frequency-magnitude distribution of earthquakes at Etna volcano unravels critical stress changes along magma pathways. _Commun Earth Environ_ 3, 68 (2022).

https://doi.org/10.1038/s43247-022-00398-6 Download citation * Received: 18 May 2021 * Accepted: 17 February 2022 * Published: 23 March 2022 * DOI: https://doi.org/10.1038/s43247-022-00398-6

SHARE THIS ARTICLE Anyone you share the following link with will be able to read this content: Get shareable link Sorry, a shareable link is not currently available for this article. Copy

to clipboard Provided by the Springer Nature SharedIt content-sharing initiative