Scattering Attenuation Images of the Control of Thrusts and Fluid Overpressure on the 2016–2017 Central Italy Seismic Sequence

Deep fluid circulation likely triggered the large extensional events of the 2016–2017 Central Italy seismic sequence. Nevertheless, the connection between fault mechanisms, main crustal‐scale thrusts, and the circulation and interaction of fluids with tectonic structures controlling the sequence is still debated. Here, we show that the 3D temporal and spatial mapping of peak delays, proxy of scattering attenuation, detects thrusts and sedimentary structures and their control on fluid overpressure and release. After the mainshocks, scattering attenuation drastically increases across the hanging wall of the Monti Sibillini and Acquasanta thrusts, revealing fracturing and fluid migration. Before the sequence, low‐scattering volumes within Triassic formations highlight regions of fluid overpressure, which enhances rock compaction. Our results highlight the control of thrusts and paleogeography on the sequence and hint at the monitoring potential of the technique for the seismic hazard assessment of the Central Apennines and other tectonic regions.

• Scattering attenuation detects the control of thrusts and lithology on post-seismic fracturing and fluid migration during the AVN sequence • Overpressurized fluids compact low-scattering rocks at thrusts' roots before earthquakes • Detecting fluid overpressure and fracturing suggests an unexploited monitoring potential of scattering attenuation

Supporting Information:
Supporting Information may be found in the online version of this article.
Correspondence to: S. Gabrielli, simona.gabrielli@ingv.it Citation: Gabrielli, S., Akinci , 1979). During the AVN, attenuation fluctuations detect permeability changes (dilatation and crack closure), pore pressure rise, and fluid diffusion , implying rupture driven by damage (pore pressure drop) mechanisms, fault weakening, and connectivity (Convertito et al., 2020). As seismic attenuation is equally sensitive to fracturing and fluid migration, measuring and modeling in space its attenuation mechanisms (scattering attenuation and absorption) provide valuable and additional information on the presence of medium heterogeneities and structural irregularities (e.g., cracks), either fluid-filled or dry (Di Martino, De Siena, Serlenga, et al., 2022;Napolitano et al., 2020;Reiss et al., 2022). Gabrielli et al. (2022) performed a 2D mapping of peak delay time and coda attenuation tomography before and after the three mainshocks of the AVN. The peak delay parameter, used as a proxy of seismic scattering, showed that structural elements (e.g., the Monti Sibillini thrust) and lithological differences (e.g., the Umbria-Marche and Lazio-Abruzzi domains) control scattering losses. Coda-attenuation mapping (as a proxy of seismic absorption) showed a difference in fluid circulation between the pre-sequence and the AVN, with high-absorption anomalies developing across the seismogenic zone during the deep migration of CO 2 -bearing fluids across the fault network. Fluid migration and changes in pore pressure during the AVN sequence were likely to cause the opening and closure of cracks in the medium, as observed from velocity variations (Soldati et al., 2019). Peak delay is therefore a likely candidate to detect fluid-related structural changes. Furthermore, peak delay imaging has been calibrated in the laboratory during deformation experiments leading to shear. King et al. (2022) demonstrate that peak delays can detect time-dependent fracturing in the sample, imaging both a post-deformation high-scattering fault zone and low-scattering zones where aftershocks are most likely. Being sensitive to the medium's structural features and fluid content, 3D mapping of peak delay has successfully detected fluid-filled faults (Di Martino, De Siena, Serlenga, et al., 2022), but it has never been applied to tectonic settings or performed as a function of time. A 4D tomographic peak-delay model focusing on different phases of the AVN could constrain the mechanisms of fracturing and fluid migration that are considered the most likely trigger of seismic sequences in the Northern-Central Apennines (Miller et al., 2004;Walters et al., 2018).
In this study, we investigate the variation in peak delay parameters both in 3D and through time, to understand the structural evolution of pore pressure, fluid migrations, and associated fracturing across the shallow crust. The results show the control of the main tectonic, geological, and paleogeographic features on the AVN mainshocks and the related fluid migration and storage.

Geological Setting
The area is characterized by an east-verging stack of thrust sheets involving the Meso-Cenozoic sedimentary successions of two main paleogeographic domains ( Figure 1): the Umbria-Marche pelagic basin (Late Triassic-Late Miocene), and the Lazio-Abruzzi carbonate platform (Late Triassic-Miocene) (Cosentino et al., 2010). These two domains are overlain by the Messinian foredeep siliciclastic deposits of the Laga Formation (LF). At depth, the Umbria-Marche and the Lazio-Abruzzi units rest on a Triassic sequence, mainly composed of evaporites and dolomites (Anidriti di Burano, Dolomia Principale).
The principal N-S and NE-SW-oriented thrusts within the study area are the Monti Sibillini thrust system (MST), bounding the southern sector of the Umbrian Arc, the Gran Sasso thrust (GST), bounded by the Monti Sibillini thrust to the west and the eastern Acquasanta thrust system (ACQ) (Billi & Tiberti, 2009;Scisciani et al., 2006; Figure 1). The structural setting also comprises high-angle normal faults deriving from Mesozoic, Miocene and Quaternary tectonic phases . The slip tendency (i.e., the ratio of the resolved shear and normal stress on a given fault surface, Morris et al., 1996) of NW-SE-oriented, high-angle normal faults, is expected to be relatively high, contrary to thrusts. Nonetheless, various authors documented the involvement and, possibly, inversion of thrusts during normal fault earthquake sequences in the Central and Northern Apennines (Bonini et al., 2019;Buttinelli et al., 2018). The cause would be the rheological conditions along some portions of the faults, either with or without the involvement of pressurized fluids Scognamiglio et al., 2018).
In Figure 1, we plot the tectonic features (thrusts and normal faults) obtained within the framework of the RETRACE 3D project Di Bucci et al., 2021). The five cross-sections described in the same figure with their structural setting are used to present the scattering attenuation results. This area has been struck by repeated seismic sequences in the past few decades, such as those in Colfiorito (1997 Mw 6.0) and L'Aquila (2009, Mw 6.3, Rovida et al., 2022), both associated with the migration of fluids and variations in pore fluid pressure (Malagnini et al., 2012). The 4D seismic tomography of the 2009 L'Aquila event, presented by Chiarabba et al. (2022), shows a post-failure migration of fluids, from the hypocenter to the surface, with an increase of water at the springs of the Gran Sasso aquifer. The seismicity of the AVN filled the seismic gap between the two Colfiorito and L'Aquila sequences. Its mainshocks are generated on a portion of the Mt. Gorzano normal fault systems and by the reactivation of an inherited major thrust ramp (the MST) within an extensional regime (Bonini et al., 2019). The relation between highly-fractured medium and overpressurized CO 2 -rich fluid, sourced from the melting of subducted carbonates migrated upward to the shallow crust (Chiodini et al., 2020) is supported during the AVN by seismic anisotropy (Pastori et al., 2019) and tomographic models (Chiarabba, Buttinelli, et al., 2020;.

Dataset and Method
We considered the combined weak-and strong-motion dataset recorded before and during the seismic sequence, already used in previous studies (Akinci et al., 2020;Gabrielli et al., 2022). The pre-sequence and AVN datasets comprise 4,731 waveforms (362 events recorded at 17 seismic stations- Figure S1 in the Supporting Information S1) and 5,373 waveforms (655 events recorded at 64 seismic stations- Figure S2 in the Supporting Information S1), respectively. The stability of the images to the changes in network configuration has been tested by Gabrielli et al. (2022). We divided the seismic sequence into three time periods to further evaluate the changes in scattering attenuation through time:  Peak delay times measure seismic envelope broadening, providing quantitative information on the strength of multiple scattering in the Earth (Calvet et al., 2013;Takahashi et al., 2007). In the near field, the parameter is also affected by multiple reflections/reverberations generated by structures with strong impedance contrasts (Gabrielli et al., 2020;King et al., 2022). To calculate and map peak delay times, we used the most recent version of the open-access code MuRAT, applied both in volcanic (De Siena et al., 2016;Di Martino, De Siena, Serlenga, et al., 2022;Reiss et al., 2022) and tectonic settings (Borleanu et al., 2017;Gabrielli et al., 2022;Napolitano et al., 2020). We define peak delay time as the lag time between the S-wave arrival and the maximum amplitude of the envelope for each filtered seismogram (Calvet et al., 2013;Takahashi et al., 2007). The seismograms are filtered with a band-pass Butterworth filter (fourth order) in the 1-2 Hz frequency band (central frequency 1.5 Hz), following the previous analysis in Gabrielli et al. (2022). The peak delay variation is computed as the difference between the measured peak delay time of each filtered ith waveform ( PD ( ) ) and the theoretical peak delay at the related hypocentral distance: Δlog 10 ( ) = log 10 PD ( ) − log 10 PD ( ) where PD ( ) = ( ) + ( ) ⋅ log 10 . ( ) and ( ) are the regression coefficients, and R is the hypocentral distance.
We mapped the spatial variations of envelope broadening after extracting its log-linear hypocentral dependence ( Figure S3 in the Supporting Information S1) with a regionalization approach. The study area is divided into 2.5 × 2.5 km horizontally and 1 km vertically-sized blocks, and we removed the cells crossed by less than 5% of the total of the rays. Then, we allocated to each block the averaged value of the logarithmic deviation of S-wave peak delay time. Positive Δlog 10 ( ) values indicate high-scattering (red) zones with heterogeneous and mostly fractured crust, while negative values are low-scattering (blue) rigid or compacted zones. More details about the processing and mapping are in Text S1 in the Supporting Information S1.

Results and Discussions
Figure 2 shows 3D views of the peak delay models obtained during the pre-sequence (top) and sequence (bottom) together with the location of the main normal faults and thrusts. Figure 2a represents the Norcia and Visso areas from the south, while Figure 2b shows the Amatrice zone from the north. We also present the spatial variation of peak delay along the five cross-sections marked in Figure 1 for the pre-sequence (Figure 3a) and the  Figures S4 and S5 in the Supporting Information S1, together with the mainshock (colored stars), thrusts (black lines and yellow triangles), and fault planes (red boxes) (Scognamiglio et al., 2018). The same vertical and horizontal sections are used to show changes in peak delay following the Amatrice (Figures 4a and S6 in the Supporting Information S1), Visso (Figures 4b and S7 in the 10.1029/2023GL103132 6 of 10 Supporting Information S1) and Norcia (Figures 4c and S8 in the Supporting Information S1) mainshocks. We added Capitignano mainshock's location in the figures to aid interpretation; this event is not part of our analysis but occurred 13 days after our dataset's last event.

Control of Tectonic and Paleogeographic Structures on Seismic Scattering
The comparison between the pre-sequence and the sequence models shows an increment of scattering attenuation after the mainshocks across the hanging wall of the MST, west of Norcia (Figure 2a). This thrust represents the contact between paleogeographic domains with different lithological features: the carbonates of the Umbria-Marche domain and the siliciclastic sediments of the Laga formation. During both sequences, the areas of low scattering between MST and ACQ remain spatially consistent, whereas the high scattering follows the shape and curvature of the MST in the sequence (Figure 2a) and is confined to the Norcia and Castelluccio basins (Figure 1). Moving south, in the area of the Amatrice earthquake (Figures 2b, S4, and S5 in the Supporting Information S1), the leading structure becomes ACQ and its splay, ACQs (Figure 1). These controllers focus most of the high pre-sequence scattering between them (Figure 3a, CC′) and at the base of the Mt. Gorzano Fault that will produce the earthquake (CC′-DD′). During the sequence, this high scattering distributes across the whole thrust (Figure 3b, CC′), with the splay bounding the bottom of the high-scattering volumes below Amatrice (Figure 2b). Substantial variations appear between 4 and 6 km depths ( Figures S4 and S5 in the Supporting Information S1), where the MST and ACQ thrusts act as rheological barriers.
Average scattering clearly defines the shallow siliciclastic deposits topping the Triassic formation before the sequence (Figure 3a, CC′-DD′, white color). During the sequence, high scattering anomalies spread across these deposits (Figure 3b, CC′-DD′), filling the hanging wall of the Mt. Gorzano fault, which becomes a structural barrier between western high and eastern low scattering. The volumes east of this fault and those below the ACQ further south (Figure 3b, EE′) are the only ones showing a drastic decrease in scattering attenuation during the sequence in our model. The Capitignano Mw5.0 event and subsequence, comprising four M > 5 events , occurred on 18 January 2017 just in this low scattering area. We note that the Capitignano subsequence events were not included in our dataset, occurring 13 days after the last earthquake of our dataset.

The Role of Thrusts and Pore Pressure in Controlling Seismic Fracturing, Strain and Fluid Migration
The increase in scattering attenuation near the locations of the mainshocks and the distribution of the anomaly across thrusts suggests a close relationship between this parameter, rheological boundaries, and the regions where fracturing and strain focus. This relationship is supported by deformation experiments performed on sandstone rock samples, where peak-delay tomography analog to the one applied here has been performed using data from acoustic emissions-the rock-physics equivalent of earthquakes (King et al., 2022). These authors demonstrated that increased scattering attenuation was associated with higher strain and fracturing, while their tomographic images detected the shear zone formed at the end of the experiment. Their results also identify stages of fracture nucleation and post-deformation fracturing with increases in peak delay and are supported by comparing peak delay variations with historical seismicity (Napolitano et al., 2020). This association would explain why the L'Aquila and Colfiorito areas, struck by similar seismic sequences in 2009 and 1997, appear as high scattering across the entire period of study. It supports the conclusion that stress is redistributed across the thrusts through fluids migrating within fractures produced by the mainshocks, visible with high-scattering.
Above and surrounding each mainshock of the AVN during the pre-sequence, significant increments in scattering are detected in regions that initially exhibited low and very-low scattering (Figures 2 and 3, AA′-CC′). King et al. (2022) observe that post-deformation low-scattering zones indeed correspond to locked regions of strain shadow and speculate that future earthquakes could nucleate across these zones in tectonic settings. However, their study analyses dry sandstone rocks, while in our region the role of fluid overpressure and diffusion is central Chiarabba, Buttinelli, et al., 2020;Convertito et al., 2020;Malagnini et al., 2022). Low scattering anomalies correspond to compact and un-fractured regions, but also to the fluid-filled rocks that can also compact when fluid pore pressure increases. This interpretation is also supported by the experimental works of Di Martino et al. (2021) and Di Martino, De Siena, and Tisato (2022), who demonstrated that the pore space topology (pore dimensions and distribution) is a primary control of peak delay. Therefore, the fluid-filled rocks generally appear as low-scattering anomalies and can only be detected by the contemporary modeling of seismic absorption (Borleanu et al., 2017).
We infer that the low scattering observed within the Triassic units across the hanging wall of the MST before the seismic sequence (e.g., section BB′) identifies volumes compacted by high pore pressure. The volumes change to high scattering due to fluid overpressure and diffusion, fracturing the overlaying rocks and distributing fluids after the mainshocks through shallow sediments. These results indicate that the stress conditions within the Triassic units are significantly affected by pore overpressure and that static friction along the MST fault surface could be relatively low, supporting the active role of this thrust in the nucleation and propagation of the Norcia mainshock (Bonini et al., 2019;Buttinelli et al., 2018Buttinelli et al., , 2021Curzi et al., 2023;Scognamiglio et al., 2018) despite its expected low slip tendency (average dip angle of ∼35°). According to the distribution of aftershocks along the deep ramp of ACQ, the same observation can be done for the ACQ and GF (Figure 3, CC′-EE′). Fluid overpressure spreads high-scattering anomalies across the thrust from the hypocenter of the Amatrice mainshock, nucleated in an area of low scattering during the pre-sequence (Figure 3, CC′). Even more relevant appears the change from low to high scattering at the root of the ACQ and the western compartmentalization of fracturing to the hanging wall of the GF (Figure 3, DD′).

Evolution of Fluid Migrations and Stress Release After Each Mainshock
The changes in scattering after each mainshock confirm the control of the paleogeographic and structural settings on fracturing following stress release and fluid migrations, marking their spatial redistribution. After the Amatrice mainshock, high scattering anomalies previously distributed within the MST footwall and the ACQ hanging progressively rise to fill the Norcia and Castelluccio basins (compare Figures 3a and 4a, BB′). This progression passes through an increase of scattering at the bottom tip of the MST thrust, around the Visso hypocenter (Figure 4b, AA′). The Norcia sequence fully compartmentalized high scattering within the MST hanging-wall ( Figure 4c, AA′-BB′).
Similar compartmentalization across the ACQ thrust and Triassic formations is visible by the end of the Norcia sequence around the Amatrice earthquake (Figure 4c, CC′). While during the pre-sequence high scattering is focused between ACQ-ACQs and the top of the Triassic formation, after Amatrice mainshock it had spread through the medium (compare Figures 3a and 4a, CC′). South of the mainshocks, high scattering progressively migrates west, permeating the low-scattering anomaly at the root of ACQ (compare Figures 3a and 4a, DD′). By the end of the Norcia sequence, the GF fault represents the primary structural boundary that confines the western fractured high-scattering volumes, with the eastern Triassic units and carbonates showing a decrease in scattering attenuation (Figure 4c, DD′).
At the south end of the studied area, we observe the progress of western compartmentalization of high scattering across the hanging wall of the GF, below the Campotosto lake (Figures 4a and 4c, EE′). A low-scattering anomaly developed during the Amatrice sequence, partially replaced by high scattering during the Norcia sequence. By the end of the sequences, ACQ has become the interface between the high-scattering shallow sediments and deeper, compact, low-scattering sediments and basement, where the 18 January 2017 Capitignano earthquake will nucleate. Our detailed temporal analysis shows that changes from low to high scattering around the mainshocks are drastic, likely caused by the sudden release of stress through fracturing and fluid migrations. The compartmentalization of seismic scattering associated with fracturing following stress-release is controlled by the presence of the thrust-bounded Triassic deposits, which have the capability of trapping crustal fluids (Curzi et al., 2023;Trippetta et al., 2010), acting as rheological asperities that are bounded and barred by the main thrusts.
Low-scattering anomalies marking the volumes where the mainshocks nucleate before the sequences, turning to high scattering afterwards, appear related to overpressured fluids compacting crustal rocks, then released by each sequence. This interpretation is supported by observing that the only regions where scattering consistently decreases by the end of the sequences are the deep portion of the pre-Triassic deposits across the hanging wall of the GF and, especially below the ACQ near Campotosto Lake, Figures 4a and 4c, DD′-EE′, where the Capitignano earthquake occurred. The comparison with section EE′ in Figure 3a, during the pre-sequence, suggests that pore pressure accumulation started soon after the Amatrice mainshock across these regions. This inference is supported by the spatial correlation through time of the low-scattering anomalies with increases in absorption, linked to fluid migrations, and following the sequence evolution (Gabrielli et al., 2022). We conclude that the scattering maps detect fluid overpressure build-up where the Capitignano earthquake will nucleate. This inference suggests an unexploited seismic monitoring and hazard assessment potential for peak delays.

Conclusions
By imaging peak delay variations, used as a proxy of scattering attenuation, we detect increases in pore pressure, fluid flow, and fracture density produced by stress changes, as previously hypothesized by rock-physics experiments. This sensitivity allows the detection of the control of regions with different rheological and lithological properties, like thrusts and Triassic formations, on the redistribution of strain and storage of fluids after the mainshocks of the AVN. The main thrusts act as rheological boundaries, where low scattering, detected in the hanging-wall of MST in the pre-sequence and associated with rock compaction, is then substituted by high scattering during the sequence phase, caused by fracturing. Our results support the active role of the main thrusts in the nucleation and propagation of the stress associated with mainshocks, even if their slip tendency is considered relatively low. The buildup of fluids in the Triassic succession prior to an earthquake causes an increase in pore pressure and a decrease in scattering attenuation. This presence of fluids changes the stress conditions and frictional properties of the thrusts supporting their role in the evaluation of AVN. During earthquakes, the release of stress through the upward movement of fluids causes fracturing in different areas of the fault and thrust network. Imaging across the AVN confirms the change from low to high scattering around the hypocenter of each mainshock. The opposite happens before the nucleation of the Capitignano earthquake, with progressively decreasing scattering in its hypocentral area. We conclude that peak delay has considerable diagnostic potential, particularly for the study region's seismic monitoring and hazard assessment.