Integrated on‐land‐offshore stratigraphy of the Campi Flegrei caldera: New insights into the volcano‐tectonic evolution in the last 15 kyr

Silicic calderas are volcanic systems whose unrest evolution is more unpredictable than other volcano types because they often do not culminate in an eruption. Their complex structure strongly influences the post‐collapse volcano‐tectonic evolution, usually coupling volcanism and ground deformation. Among such volcanoes, the Campi Flegrei caldera (southern Italy) is one of the most studied. Significant long‐ and short‐term ground deformations characterise this restless volcano. Several studies performed on the marine‐continental succession exposed in the central sector of the Campi Flegrei caldera provided a reconstruction of ground deformation during the last 15 kyr. However, considering that over one‐third of the caldera is presently submerged beneath the Pozzuoli Gulf, a comprehensive stratigraphic on‐land‐offshore framework is still lacking. This study aims at reconstructing the offshore succession through analysis of high‐resolution single and multichannel reflection seismic profiles and correlates the resulting seismic stratigraphic framework with the stratigraphy reconstructed on‐land. Results provide new clues on the causative relations between the intra‐caldera marine and volcaniclastic sedimentation and the alternating phases of marine transgressions and regressions originated by the interplay between ground deformation and sea‐level rise. The volcano‐tectonic reconstruction, provided in this work, connects the major caldera floor movements to the large Plinian eruptions of Pomici Principali (12 ka) and Agnano Monte Spina (4.55 ka), with the onset of the first post‐caldera doming at ca. 10.5 ka. We emphasise that ground deformation is usually coupled with volcanic activity, which shows a self‐similar pattern, regardless of its scale. Thus, characterising the long‐term deformation history becomes of particular interest and relevance for hazard assessment and definition of future unrest scenarios.


| INTRODUCTION
Caldera systems generally show complex ground deformation phenomena mainly associated with variations of magma production, uprise, storage and migration in the volcano plumbing system through different timescales (e.g., Cashman & Giordano, 2014;Kennedy et al., 2018). Overall, long-term ground deformation represents the surface expression of the incremental addition of magma within complex reservoirs (e.g., Acocella, 2019;Marsh, 1984;Townsend et al., 2019). Usually, the ground movements are coupled with volcanic activity and may either precede or follow eruptions. However, they can also occur with no associated eruptive activity, especially at silicic calderas (Acocella et al., 2015;Newhall & Dzurisin, 1988). Reconstruction of the complex interplay between volcanism and ground deformation is thus one of the main challenges in studying caldera systems, especially in active volcanic areas. Thus, an accurate tephrostratigraphic framework becomes the primary element in defining and quantifying the area affected by volcano-tectonic processes and reconstructing the main structures governing these phenomena and their evolution (e.g., Lucchi et al., 2007;Martì et al., 2018). The stratigraphic model, obtained by coupling the facies analysis of continental and marine sedimentary records with tephrochronology, allows reconstructing the timing of the constructive and destructive volcanic phases characterising calderas, which may include large submerged portions (e.g., Santorini, Greece; Rabaul, Papua New Guinea; Deception Island, Antarctica; Campi Flegrei, Italy). Nonetheless, accurate reconstructions of the volcano-tectonic evolution are fundamental for any attempt to provide eruptive forecasting scenarios.
The analysis of the marine sedimentary record provides a powerful tool for reconstructing both the stratigraphic framework and volcano-tectonic evolution of CFc.
Notably, the study of the marine and volcanic succession exposed at La Starza coastal cliff nearby Pozzuoli Town (Figure 1a) allowed the reconstruction of the long-term caldera ground vertical movements over the last 15 kyr (Cinque et al., 1985;Giudicepietro, 1993;Isaia et al., 2019;Rodriquez, 1964;Vitale et al., 2019). However, because more than one-third of the caldera is currently submerged, marine geophysical data can significantly help to achieve a more detailed stratigraphic reconstruction.
With this scope, several studies have assessed the stratigraphy of the offshore caldera basin-fill using seismic reflection datasets of different vintages and resolution (Aiello et al., 2012Fusi et al., 1991;Milia, 1998;Natale et al., 2020;Pescatore et al., 1984;Sacchi et al., 2014Sacchi et al., , 2019Steinmann et al., 2016Steinmann et al., , 2018. In this work, we benefit from a blend of high-and very-high-resolution seismic profiles, and provide a detailed stratigraphic framework for the submerged part of the caldera and integrate the marine record of the Pozzuoli Gulf with both the La Starza marine succession on-land and the whole continental tephrostratigraphic records of the caldera. Thus, this study significantly updates the stratigraphic framework of the Pozzuoli Gulf, giving rise to a renewed understanding of the long-term interplay between caldera floor motion and sea-level rise as a proxy for the volcano-tectonic events that accompanied the eruptive activity within the caldera.

| Campi Flegrei caldera
The CFc is located along the Tyrrhenian Sea margin of southern Italy (Figure 1; Vitale & Ciarcia, 2018 and references

Highlights
• A detailed analysis of high-resolution seismic profiles in the offshore of Campi Flegrei caldera is presented • We provide a reconstruction of the offshore stratigraphy of the last 15 kyr • A correlation between seismic units with the exposed La Starza succession is provided • The offshore caldera sector recorded the interplay between volcanism, doming, and sea-level rise • Stratigraphic reconstruction allowed us better to define the relative age of some coastal eruptions therein). Volcanism in the broader CFc area started from at least 80 ka and was characterised by a scattered explosive activity associated with monogenetic cones, lava domes and maars Sbrana et al., 2015;Scarpati et al., 2013;Vitale & Isaia, 2014), including some volcanic edifices in the urban area of Naples. The present 12 km sized caldera ( Figure 1a) is the complex result of at least three very large explosive eruptions: the Campanian Ignimbrite (39.8 ka, up to about 250 km 3 DRE, Giaccio et al., 2017Giaccio et al., , 2008Costa et al., 2012;Scarpati & Perrotta, 2016, Scarpati et al., 2020Silleni et al., 2020); Masseria del Monte Tuff (29.3 ka, 17 km 3 DRE, Albert et al., 2019) and Neapolitan Yellow Tuff (NYT, 14.9 ka, 20-40 km 3 DRE, Deino et al., 2004;Orsi et al., 1992;Scarpati et al., 1993). The current morphological layout chiefly derives from the NYT caldera collapse and subsequent intracaldera activity (e.g., Orsi et al., 1996;Rosi & Sbrana, 1987;Vitale & Isaia, 2014). Volcanism in the last 15 kyr mainly took place in distinct periods alternating with quiescence time lapses. Three main eruptive epochs ranging between 15 and 10.6 ka (Epoch 1), 9.6 and 9.1 ka (Epoch 2) and 5.5 and 3.5 ka (Epoch 3) have been set out in literature Isaia et al., 2015;Smith et al., 2011). Most eruptions were distributed within the collapsed area and subordinately along its borders, varying in style and magnitude with a volume of erupted magma exceeding 1 km 3 during the major event of Pomici Principali (12.3 ka) and Agnano-Monte Spina (4.55 ka) Romano et al., 2020). Eruptions in the last 5 kyr mainly occurred in the eastern caldera sector (Bevilacqua et al., , 2016Neri et al., 2015), clustered along structural discontinuities (Isaia et al., , 2009Vitale & Isaia, 2014) and, in at least one documented case, simultaneously occurred in different caldera sectors (Pistolesi et al., 2016). The estimated cumulative volume of erupted magma (DRE) varies during Epochs 1, 2 and 3, with values of 4.2 ± 0.7, 0.5 ± 0.1 and 2.6 ± 0.5 km 3 , respectively (Bevilacqua et al., 2016). F I G U R E 1 (a) Digital Elevation Model of the Campi Flegrei caldera, showing on-land caldera rims (modified after Steinmann et al., 2018;Vitale & Isaia, 2014), offshore ring fault zones (this work), dome surface area (Bevilacqua et al., 2020); dome faults (Natale et al., 2020), crater rims, scoria cones and lava domes. Dashed boxes indicate the location of time-structure maps of Figures 10 and 11. (b) Locations of the seismic reflection lines used in this work Significant ground deformation and seismo-tectonic activity accompanied volcanic activity . The ground deformation of the caldera floor is wellrecorded in its central sector by the marine-transitional sedimentary sequence of the La Starza unit Vitale et al., 2019) that includes well-dated tephra layers. Presently, this succession crops out along the Pozzuoli coast, extending from Monte Nuovo volcano to the west and Bagnoli to the east (Figure 1a), and it has been drilled by numerous wells Isaia et al., 2019;Sbrana et al., 2015). The caldera floor motion usually preceded and accompanied volcanic activity . Following the Agnano-Monte Spina Plinian eruption (AMS, 4.55 ka, De Vita et al., 1999;Smith et al., 2011), ground subsidence occurred rapidly, followed by an uplift phase that raised marine sediments to 20 m a.s.l. (Isaia et al., 2009). The last eruption, Monte Nuovo, occurred to the west of Pozzuoli town ( Figure 1a) in 1538 CE, following a few centuries of ground uplift (Guidoboni & Ciuccarelli, 2011).
In the early 1950s and the 1970-1972 and 1982-1984 periods, CFc experienced significant unrest (up to 1 m/ year of uplift and peaks of thousands of earthquakes per month), which is still ongoing today at slower rates (0.1 m/ year of uplift; hundreds of earthquakes per month; INGV, 2021). Nevertheless, the normalised ground deformation shape (vertical and horizontal components) is remarkably constant regardless of the total displacement during uplift and subsidence phases. Furthermore, several works point out that the normalised deformation shows an axisymmetrical bell-shape, with a maximum at its centre and radially decreasing from the caldera centre to its periphery during uplift and subsidence phases (e.g., Bevilacqua et al., 2020). For this reason and in analogy with the historical monitoring data, we consider the uplift and subsidence phases as having this deformation shape.

Pozzuoli Gulf
The Pozzuoli Gulf encompasses an area of ca. 40 km 2 , covering a considerable part of the CFc (Figure 1a). In the last 15 kyr, the interaction between sea-level rise, ground deformation, erosion and sedimentation strongly affected its morphology (e.g., Sacchi et al., 2014;Somma et al., 2016;Versino, 1972). The offshore hosts the complementary portion of the subaerial caldera structures, including in the central sector a 5 km radius dome, roughly centred in Pozzuoli town, with the shallower sector (culmination) ranging between 20 and 40 m b.s.l. and bounded by a shelf break at variable depths ( Figure 1a). The peripheral area (caldera collar) encompasses the deepest part of the gulf, down to 120 m b.s.l, with the Epitaffio and Bagnoli valleys to the west and east, respectively. The caldera collar is an annular morpho-structural depression between the dome and the ring fault zone. This sector also comprises some morpho-structural highs, made of remnants of monogenetic volcanic edifices that bound the Pozzuoli Gulf from Capo Miseno to Punta Epitaffio to the west and from Nisida Island and Coroglio cliff to the east. The sector from Porto Miseno to Nisida localises in correspondence of the caldera ring fault zone, which is confined by two main boundaries (Figure 1a). The inner caldera boundary is defined by segmented normal faults with tens of metres of displacement associated with the NYT caldera collapse (Steinmann et al., 2016). The depth and the amount of the NYT down throwing are comparable to that estimated from on-shore drill holes Rosi & Sbrana, 1987). Southward, the outer caldera boundary is associated with the CI caldera collapse, with its footwall corresponding to the morpho-structural high made up of pre-CI volcanic banks (Steinmann et al., 2018). These fault strands join with the on-land lineaments associated with the CI caldera faults from Monte di Procida and Posillipo (e.g., Vitale & Isaia, 2014). Intense seafloor degassing and local intrusions localise in the area included by these two main fault zones, such as the Punta Pennata structure, along the coast of Bacoli, and the Bagnoli intrusion (Steinmann et al., 2018).
The most recent studies on the CFc offshore sector on stratigraphy and the interplay between caldera activity and sea-level fluctuations have been proposed by Steinmann et al. (2016Steinmann et al. ( , 2018. The authors used a 3D approach extending the identified bounding unconformities to the whole dataset to improve the previous 2D conceptual reconstruction (Steinmann et al., 2016). Their work also assessed an evolutionary model of the caldera system more thoroughly as recorded by their sedimentary units. In the present study, the original three units above the NYT defined by previous works have been deeply revised and detailed up to 11 units, with a different age interpretation, giving rise to a new understanding of the interplay between caldera floor motion, volcanic activity and sea-level rise.

| Seismic profiles
Our reconstruction of the Pozzuoli Gulf stratigraphy was mainly based on a set of high-resolution seismic profiles (Seistec_2013) acquired using the uniboom IKB-Seistec profiler (Mosher & Simpkin, 1999;Sacchi et al., 2019;Simpkin & Davis, 1993) by ISMAR-CNR of Naples in 2013 (Figure 1b). With a source frequency range between 1 and 20 kHz and a central frequency of 2 kHz, the theoretical vertical resolution is 0.2 m, with an average velocity of 1500 m/s. The horizontal resolution at central frequency is ca. 0.5 m at 150 ms TWT (Figure 2a). The seismic signal penetration is strongly influenced by the high signal attenuation of the volcaniclastic rocks, especially in tuffs and fluid-saturated sediments, and it typically reaches Two Way Travel (TWT) depths between 100 and 150 ms.
To gather information on the deeper post-NYT reflectors, which are not adequately imaged in the Seistec dataset due to loss of signal, and to identify the base of the analysed succession (top of NYT), we relied on the multichannel seismic lines from the GeoB08 dataset ( Figure 1b, acquired during the CAFE-7/3 expedition in 2008; Steinmann et al., 2016). These profiles have a vertical resolution of ca. 2 m with a central source frequency of 250 Hz, out of 0.3-1 kHz range, with an average velocity of 1500 m/s. The horizontal resolution at central frequency is ca. 1 m at 150 ms, with a signal penetration exceeding 350 ms (Steinmann et al., 2016). The processing procedure of Steinmann et al. (2016) improved overall imaging depth compared with previous works in literature, allowing to follow seismic units beneath the shallow multiple reflections and the fluid-bearing layers. Seismic lines belonging to the two datasets were imported and analysed in the GeoSuite AllWorks (Geo Marine Survey System, 2012) and Eliis Paleoscan software packages (Eliis, 2018). To estimate the thickness of seismic units, we used a TWT time velocity of 1650 m/s, which is considered to be representative for shallow (<200 m), mostly unconsolidated volcaniclastic and marine deposits (Sacchi et al., 2014;Steinmann et al., 2018).

| Stratigraphic approach
The stratigraphic reconstruction proposed in this work is based on integrating on-land outcrop and borehole data retrieved from the literature (i.a., Isaia et al., 2019) with multiscale seismic reflection data (Figure 1b). Concerning the offshore, we started our reconstruction on three reference units, which include the unit NYT (Steinmann et al., 2016), the well-dispersed transparent unit marker B (Natale et al., 2020) and the shallowest layers sampled by gravity corings (Sacchi et al., 2014). The methodological approach used to increase the stratigraphic resolution includes identifying the Unconformity Bounded Stratigraphic Units (UBSUs), particularly indicated for the volcanic environment (Lucchi, 2019). The UBSUs Interpreted N-S oriented seismic profile GeoB08-050 highlighting the relationships between lower sequences and caldera ring faults. V.E.: Vertical Exaggeration. For locations, see bottom-right insets reconstructed offshore have been correlated with the corresponding marine successions (named 'intervals') and their bounding unconformities defined on the La Starza on-land by Isaia et al. (2019). We emphasise that the application of UBSU correlations might suffer limitations due to different factors, including undersampling of onland data due to incomplete records, hiatuses and erosion. Despite these limitations, this methodologic approach represents a valid alternative to overcome the lack of direct stratigraphic information (e.g., well-logs).
Using the correlation between gravity core data (Sacchi et al., 2014) and seismic stratigraphic units, reflectors with high-amplitude and medium-frequency reflections were interpreted as coarse or poorly sorted volcaniclastic material intercalated within fine-grained, well-sorted marine sediments (Natale et al., 2020). The layered succession of marine deposits with different textures (e.g., lithology, grain size, matrix/particle ration, sorting) causes a significant variation of the acoustic impedance that produces strong reflection amplitudes. On the other hand, reflectionfree bodies correlate to successions with prevailing homogeneous medium-fine grain-sized marine sediments (e.g., marker B of Natale et al., 2020). These deposits occasionally contain laterally discontinuous coarse beds likely represented by reworked volcaniclastic and marine deposits.

| Sedimentation rates
In order to estimate the sedimentation rates for the main seismo-stratigraphic units identified here, we measured the thickness of each unit and sub-unit in a fixed area of the caldera collar on each seismic profile (Supplementary table). Then, we averaged the thickness values and calculated the sedimentation rates by dividing the average thickness by the duration of the corresponding time interval. We performed the computation at the caldera collar insofar this area is less affected by volcano-tectonic motions, and thus potentially, all the seismostratigraphic units are recorded. Thus, the established rates are to be considered a maximum estimation. More details of this procedure are reported in the supplementary material.

| Definition of seismic units
To correctly display the recognised unconformities that bound the seismic units in different morphological settings, we show some full-length seismic profiles from the Seistec dataset, oriented about N-S and E-W ( Figures   2-4). We complemented the high-resolution profiles with three overlapping GeoB08 lines to highlight seismic interpretation consistency and further extrapolate the reconstructed stratigraphic framework at depth.
The first two complementary lines ( Figure 2) extend from the inner continental shelf of the central sector to the south of Pozzuoli bay for ca. 5 km. They illustrate the main seismic units, with greater resolution at shallow depth for Seistec-402 ( Figure 2a) and higher penetration for the GeoB08-050 (Figure 2b), where the main structural features of the dome and the ring fault zones are evident. The caldera faults in these lines do not cut units younger than the NYT.
Conversely, the 5 km long seismic profiles Seistec-401 ( Figure 3a) and GeoB08-121 ( Figure 3b) extend from Pozzuoli to Capo Miseno, crossing the Punta Pennata structure. The latter is an uplifted zone located south of the southwestern segments of the inner caldera ring (Steinmann et al., 2018). Furthermore, two ca. W-E-oriented seismic lines, 6 km long Seistec-801 ( Figure 4a) and 7 km long GeoB08-065 ( Figure  4b), broadly stretch from Baia, across the central dome, to Nisida Island. Finally, the submerged dome, the offshore extension of the continental dome (Bevilacqua et al., 2020), is strongly affected by high-angle normal faults (Natale et al., 2020), which lowered the dome culmination.
We identified different discontinuity surfaces (indicated, from bottom to top, H0 to H10) that, together with stratal architecture and seismic features, allowed us to distinguish eleven seismic units overlying the unit NYT (Table 1), which is the relative base of our high-resolution reconstruction. We labelled them as S0 to S10, with some of them including up to three subunits (Table 1). We used different criteria to identify these units, including, in order of priority: (i) boundaries geometry; (ii) lateral continuity; (iii) seismic facies and (iv) internal structure. The chronostratigraphy of S10 and the upper part of S9 has been calibrated with gravity cores (Natale et al., 2020;Sacchi et al., 2014). Unlike younger units, NYT, S0, and S1 are evident only on GeoB08 lines.
In the following paragraphs, we illustrate the main features of each seismic unit, as identified in the available seismic profiles and shown in Figure 5 and Table 1. To improve the readability of the seismic features, we indicate the unit-bounding horizon (Hn) with the same color code as the overlying unit.

| Unit NYT (Neapolitan Yellow Tuff)
Seismic unit NYT is the base of the post-caldera basin fill, as defined offshore by Steinmann et al. (2016). NYT generally exhibits chaotic seismic facies with mainly discontinuous, low-to-high amplitude reflections. To the west, the seismic facies evolves laterally to more layered medium-amplitude reflectors ( Figure 4b). NYT overlays continuous high-amplitude reflectors laterally and is capped by paraconformity H0 (U2 in Steinmann et al., 2016Steinmann et al., , 2018. The NYT unit is displaced by normal faults across the caldera ring zone, showing significant lateral thickness variations (Figures 2b, 3b, and 4b).

| Units S0-S1
An overall sub-transparent facies characterises unit S0 with the local occurrence of some medium-amplitude, more or less continuous reflectors. It has an average thickness of ca. 20 m, which keeps almost constant moving from the present shelf to the basin until the inner border of the caldera, where it is confined. Unit S0 is bounded upward by a sharp seismic amplitude and continuity variation, which defines unconformity H1.
Compared with S0, the overlying unit S1 is defined by increased amplitude and continuity and a slight decrease of the reflection frequency. The sequence has an average thickness of ca. 10 m with no significant lateral variations (Figures 3b and 4b).

| Unit S2
A sharp variation in reflection attributes characterises the transition between S1 and S2 units across H2. Unlike the underlying units, S2 shows high-amplitude, high-frequency continuous reflectors. It includes three subunits, from bottom to top: S2a, S2b and S2c (Table 1, Figure 5a,b). The lowermost S2a subunit exhibits constant thickness (ca. 8 m) with high-amplitude reflectors and upward frequency decrease (Figures 2 and 3). It is distributed throughout the Pozzuoli Gulf and represents a continuous seismostratigraphic marker for the two seismic datasets.
Seismic facies and geometric features sharply change within overlying subunit S2b, showing marked transparent facies with poor lateral continuity and, locally, some parallel reflectors. Its thickness is lesser than that of S2a and noticeably increases towards the northern part of the gulf (Figures 2a, 3a, 5a,b).
The uppermost S2c subunit lies concordantly on top of S2b, showing more reflective facies, with continuous medium-amplitude and medium-frequency reflectors. The subunit is topped by paraconformity H3 (Figure 5a,b). This contact is represented by a sharp, high-amplitude, undulated to plane reflector.

| Unit S4
S4 is composed of three subunits (S4a-c, Figure 5c). S4a is located at the base and is characterised by stratified medium-amplitude, medium-frequency reflectors, with higher amplitude upward, alternating with transparent reflectors mainly in its lower part with increasing thickness towards the collar area. Subunit S4b is characterised by sub-transparent facies with locally poorly continuous reflectors mainly towards its top and shows a gentle decrease in thickness from the basin to the dome. Finally, a sharp contact characterises the transition to the overlying S4c (Figure 5c), represented by high-amplitude continuous reflectors with significant thickening passing from the dome apical zone to the collar.

| Unit S5
Sequence S5 features two sub-units (S5a-b). Subunit S5a is a well-recognisable stratigraphic marker in both Seistec and GeoB08 lines (Marker B in Natale et al., 2020) that exhibit transparent facies with discontinuous low-amplitude reflectors towards the top (Figures 2-5c). It is characterised by marked thickness variations recording up to ca. 6 m in the collar zone. Notwithstanding the occurrence of erosive features in the dome area, the thickness of the unit is broadly preserved and is about 2 m. Subunit S5b includes medium-to-high amplitude parallel to sub-parallel reflectors, with a slightly variable thickness (averagely around 0.5 m; Figure 5d), generally preserved in the dome area.

| Unit S6
Unit S6 overlies with a paraconformity horizon H6 and is characterised by high and very high-amplitude, continuous to chaotic reflectors in the basin. The basal contact evolves laterally to an angular unconformity towards the dome apical zone.
This unit encompasses subunits S6a-c, separated by internal unconformities and representing discrete depositional events (Table 1, Figure 5e). Although they display similar facies and lobate geometry, the uppermost subunit shows larger thicknesses (up to 2.5 m) along the slope and the collar. All three subunits are very thin or absent in the apical dome area, except in structural lows, thus predating the main abrasion surface above. Onlap terminations are present along the slope of the dome (Figure 5e).

| Unit S7
S7 shows a lobe-like, transgressive pattern in paleobathymetric lows and along N-S-oriented lines (Figure 5e,f). This unit comprises the basal subunit S7a, defined mainly by transparent facies with interspersed high-amplitude continuous reflectors, displaying unconformable contact with the underlying sequence. On the other hand, a sharp change in seismic facies, which shifts towards higher amplitude and frequencies, is present within overlying subunit S7b, which displays continuous to chaotic reflectors. Moreover, in the eastern sector, in the Bagnoli valley, blocky to chaotic features characterise this subunit, reaching a maximum thickness of ca. 10 m (Figure 5f).

| Unit S8
In the basin, S8 overlies paraconformity H8, which evolves to an angular unconformity towards the dome area, producing a pronounced abrasion surface (Figure 5g, area enclosed by the dashed line in the map inset of Figure 5). We divided this unit is into three subunits, as significant facies and geometric variations are recognised within the F I G U R E 5 (a, b) Landward thickening of unit S2, bounded on top by horizon H3, passes landward from a paraconformity to an unconformity, as seen on Seistec GeoB08 lines, respectively. (c) Main characteristics of units S3 and S4. (d) Unit S5 as observed in domeradial line Seistec-223, with the alternating amplitude reflectors. (e) Unit S6 and its three subunits and unit S9 seismic expression. (f) Anticlinal structure capped by chaotic blocky and boulder rich deposit. (g) Major abrasion surface (H8) topping the apical dome area and associated with S7 unit. (h) Two generations of abrasion surfaces and prograding wedges on the eastern part of the gulf. (i) Regressive stacking pattern marked by shifted onlap terminations and subdivision of S8 unit. (j) Stratigraphic position of Capo Miseno tephra located in the upper part of S8b subunit. (k) Transparent unit S9 with interbedded 79 CE tephra, and unit S10, with projected cored interval by Sacchi et al. (2014). (l) Map inset showing the location of the seismic lines, dashed line encloses the extension of the abrasion surface sequence. On N-S-oriented lines (Figures 2 and 3), subunits S8a and S8b show prograding onlap terminations and lobate features (Figure 5i).
The lowermost subunit S8a is marked by chaotic midamplitude and high-frequency reflectors showing a slight thickness increase from east to west. S8b is mainly composed of laterally discontinuous high-amplitude reflectors. This unit has a flat-lying top surface, and on the western side, it pinches out on the Punta Pennata structure ( Figure  3a), whereas on the flank of the dome, it passes laterally to chaotic and irregular reflectors.
Finally, S8c, which represents a stratigraphic marker on both datasets, is composed of parallel stratified high frequency and amplitude reflectors that mantle the underlying sequences. Furthermore, this horizon was cored and calibrated by three gravity cores (Sacchi et al., 2014). To the west and east, its overall thickness increases towards Capo Miseno ( Figure 5j) and Nisida (Figure 3), respectively. 4.1.10 | Units S9-S10 Unit S9 features low-amplitude, locally transparent reflectors with a high-amplitude, continuous reflector halfway its thickness and includes the Vesuvius 79 CE and Ischia 60 CE tephras, as recognised by Sacchi et al. (2014). The overall thickness in the collar zone ranges between 4 and 7 m, reaching its maximum value in the Epitaffio and Bagnoli valleys (Figure 4).
Laterally discontinuous high amplitude reflectors characterise S10 with internal geometries and valley-ponding features in the Epitaffio valley, where it reaches its maximum thickness. This sequence has been calibrated with gravity cores and correlated to the 1538 CE Monte Nuovo tephra (Figure 5k).

| Deformation structures
Seismic images reveal that faults and intrusion-related structures thoroughly affect the caldera offshore. Two main fault systems are present: ring faults and dome high-angle faults (Figure 1a). The ring fault zones, where numerous fault strands are present, are the place of an intense fluid uprise (Figures 2 and 3). The inner ring faults displace the NYT unit and generally bound unit S0, whereas they are topped by S1 (Figure 2). Nevertheless, in the western sector, these faults deform sediments up to unit S5 (Figure 3). In concomitance with the ring faults, the signal is often blanked due to degassing (Figures 2 and  4). Nearby Porto Miseno (Figure 1a), a deformation structure occurs (Punta Pennata structure; . The latter is a complex intrusion in the footwall of the inner ring fault, likely related to a shallow intrusion (Steinmann et al., 2018). In fact, the caldera ring faults represent preferential pathways for the ascent of magma, promoting the emplacement of intrusions (Magee et al., 2013;Saunders, 2004;Smith & Bailey, 1968). Furthermore, in the eastern part of ca. E-W profiles (Figure 4), a welldefined anticlinal structure occurs, already interpreted by Sacchi et al. (2014) and Steinmann et al. (2016) as related to a sub-surficial intrusion.
The central dome, broadly included in the dashed circular area in Figure 1, is affected by near-vertical faults that displace the seismostratigraphic succession up to unit S6. It is worth noting as the main abrasion surface associated with S7-S8 units cuts this fault system. These structures generally show smaller displacements than ring faults (Natale et al., 2020)

| Correlation of the seismostratigraphic units with the La Starza succession
We correlated the reconstructed stratigraphy of seismic units with the coeval La Starza succession exposed on-land. The latter represents the depocentre of the depositional basin ensuing the NYT collapse, whereas the offshore seismic units were formed in progressively more distal sectors. Figure 6 shows the proposed correlation between the eleven units identified in seismic profiles with the main stratigraphic intervals (a to h3 in Figure 6b) defined for the La Starza unit by Isaia et al. (2019) and Vitale et al. (2019). Figure 7 shows a detailed visual correlation between seismic sequences ( Figure 5) and their exposed counterparts observed at La Starza cliff and in a tunnel excavation (after Isaia et al., 2019).
The analysis of the La Starza succession indicates that, in the last 15 kyr, the proximal sector was characterised by alternating marine and subaerial environments . In contrast, the Gulf of Pozzuoli has always been submerged during this time interval, including phases of prevailing net uplift of the seafloor when a sharp marine surface was formed in the inner shelf area. Therefore, to compare the seismic units with La Starza deposits, we use the interval names (a-h) of Isaia et al. (2019), as indicated in Figure 6b and reported in Table 1.
Based on the seismic stratigraphic analysis, we interpret the first intracaldera sedimentary unit (S0) above the NYT as a result of reworking of NYT volcaniclastics, the deposition of marine and volcaniclastic sediments, corresponding to tephras likely generated during the period of the early Epoch 1 (14.9-14.3 ka), when explosive eruptions mainly resulted in the building of monogenetic tuff edifices, such as Mofete, Gauro, Bellavista, and La Pietra (Figures 1a, 2,  3, and 6a, b, d). The Unit S0 becomes significantly thicker landward and shows some distinct reflectors corresponding to the La Pietra tuff cone (Figure 2).
Likewise, in GeoB08 images, unit S1 is defined by the repeated alternation of transparent and reflective layers, suggesting a recurrent tephra input into the basin. We associate the reflective horizons of this unit with volcanic products generated mainly in the central northern sector of the caldera , such as La Pigna, Archiaverno and Minopoli (14.3-12.3 ka;Smith et al., 2011). Because this unit shows no significant lateral thickness variations, in this time frame, the Pozzuoli Gulf behaved as a relatively more distal sedimentary sink, where volcaniclastic inputs dominantly filled the accommodation space created by the sea-level rise. The highly reflective subunit S2a correlates with the basal interval a of the La Starza succession, with a constant thickness of ca. 8 m. We suggest that the stronger reflectors in the central upper part of the subunit (Figure 5a-c) correlate to Soccavo 1 tephra and the high magnitude eruption of F I G U R E 6 (a) Chronostratigraphic log of the post-NYT volcanism, (b) La Starza stratigraphic intervals (modified after Isaia et al., 2019), (c) Seistec seismic type section with also indicated subhorizons, (d) GeoB08 seismic type section, (e) seismic sequences and corresponding color code adopted in the figures Agnano-Pomici Principali (PP; 1 km 3 DRE) that occurred at 12.1 ka, in the central part of the caldera (Figure 6; Bevilacqua et al., 2015;Di Vito et al., 1999;Smith et al., 2011). At this stage, the succession nowadays exposed at La Starza was formed in a coastal-marine transitional environment , suggesting a balance between the volcaniclastic input and sea-level rise. As for the underlying unit S1, lateral thickness variations are not evident (Figure 5a-c), possibly highlighting the absence and/or very low rates of ground deformation, both uplift or subsidence.
Weakly reflective subunits S2b and S2c display a landward (northward) increase in thickness and are correlated to the La Starza interval b, where a few intercalated tephra layers occur within the foreshore to shoreface marine sediments. This feature indicates the general subsidence of the central sector fed by mainly marine sediments and, to a lesser extent by tephra deriving from minor eruptions in the northern sector, unable to balance the contemporaneous sea-level rise.
The overlying H3 surface is a paraconformity in the basin, but evolves to an angular unconformity moving landward, truncating older reflectors (Figure 5a-c). This surface correlates to the lowermost and well-developed angular unconformity visible on-land located at the base of interval c, which truncates the underlying marine F I G U R E 7 Excavation front within the La Starza unit, with angular unconformity H3 between units S2 (intervals a-b) and S3 (interval c), inset shows outcrop picture of interval a; (b) Starza interval c-d; (c) interval d, e, f at La Starza tunnel; (d-e) Detail of interval g1 at La Starza outcrops and (f) at the tunnel excavation front; (g) Epoch 3a pyroclastic deposits in contact with La Starza unit; (h) Epoch 3a tephra overlying La Starza Unit g1 and topped by a paleosol separating it from Epoch 3b tephra; (i-j) Chaotic and coarse transitional deposits located at 30 and 20 m a.s.l. at La Starza and La Pietra, respectively; (k) Abrasion surface (H8) observed on-land truncating AMS tephra and covered by marine sands and Accademia tephra; (l) Epoch 3 tephra exposed on top of La Starza terrace sediments (Figure 7a), whose top is dated at 10.46 ka (Isaia et al., 2019 and reference therein). This unconformity and the overlying lacustrine sediments indicate shallower sedimentary conditions compared with the underlying strata. This feature highlights the first clear evidence of ground uplift, likely related to doming. Consequently, we suggest that the onset of the inner caldera doming and the shallowing of the central sector, highlighted by the occurrence of lacustrine sediments at La Starza, occurred between 10.46 and 9.6 ka (interval c, Isaia et al., 2019, Figures 6 and 7a-c). It is not trivial to disclose whether this uplift phase lasted at least until 9.6 ka, or if the uplift pulse was shorter in time but large enough to preserve lacustrine conditions in the central sector despite the rapid sea level in this time interval (see Figure 9a for the sea level reference curve).
On the other hand, the occurrence of seismites in the same lacustrine beds suggests coeval seismo-tectonic activity, typical of doming-related unrests, as reported by Vitale et al. (2019). Furthermore, the presence of gravity flows in subunit S3 strengthens the hypothesis of protracted uplift. In this scenario, unit S3 is thus interpreted as the offshore counterpart of the lacustrine and coastal sediments of interval c, broadly corresponding to the short-lived absence of volcanic activity between Epochs 1 and 2. In some continental locations, this brief hiatus allowed the formation of the paleosol A (Di . The overlying S4 corresponds to interval d, characterized by foreshore marine sand intercalated with tephra beds related to Epoch 2 volcanic activity. In particular, the highly reflective reflectors alternating with subtransparent seismic facies and westward increasing thickness suggest the correlation of unit S4a with the eruptive products of Baia-Fondi di Baia on the western side of the bay (Pistolesi et al., 2017). The volcaniclastic contribution from younger phreatomagmatic eruptions of Epoch 2, such as Costa San Domenico and Pigna San Nicola, might be responsible for the facies characteristics of subunits b and c. Moreover, the S4c highly reflective facies and its stacking pattern are compatible with the culmination of uplift of Epoch 2, which resulted in increased coarser sedimentary input likely due to remobilisation of volcaniclastic products testified by the emersion of the La Starza unit (interval e). This uplift has occurred around 9.1 ka  and is constrained to at least before 8.59 ka, the age of a Cardium fossil shell at La Starza cliff located immediately above the interval e (Giudicepietro, 1993;Isaia et al., 2019).
The full transparent facies of thicker unit S5a suggests a close correlation with the interval f, characterised by finemedium massive sands, which marks the reestablishment of progressively distal and deeper conditions between 9.1 and 5.86 ka, linked to a new submersion of La Starza. Unit S5b is correlated to interval g1, marked by a recurrent alternation of layers rich in Posidonia rizholiths and medium-size marine sands interpreted as having formed in a shoreface and offshore environment, respectively, within an overall shallowing pattern. In seismic images, this feature is reflected in alternating high amplitude and transparent reflectors (Figures 5d and 7d-f). This setting testifies a renewed doming with a cyclic shallowing of the paleoenvironment, suggesting discontinuous increments of ground uplift, possibly linked to volcano-tectonic crises that predate the eruptive activity of Epoch 3 (Figure 7b).
The strong amplitude reflectors within overlying S6a are interpreted as tephra layers associated with the eruptions that occurred at the onset of Epoch 3a (Averno 1, Agnano 1-2; Smith et al., 2011;Isaia et al., 2019) and included in interval g2. We suggest that the following subunits S6b and S6c, similarly characterised by very strong-amplitude reflectors, include pyroclastic/volcaniclastic products of the main eruptive units of Epoch 3a The overlying unit S7, which includes a well-visible transparent layer, represents a stratigraphic marker within the offshore area. S7 attains a significantly larger thickness within the collar zone depocentres than in the dome area and on-land. This marker correlates to shallow marine to coastal marine deposits of the Pozzuoli unit (interval h2), which testifies to a short-lived phase of ground subsidence that followed the Agnano-Monte Spina Plinian eruption and the associated caldera collapse of the Agnano Plain (Isaia et al., 2009). This subsidence allowed the coastline to shift northward up to the San Vito Plain (Sbrana et al., 2015) (Figure 1a) and the deposition of coastal sediments that now are exposed, uplifted, 20 m a.s.l. on La Pietra tuff (Isaia et al., 2009) (Figure 7i,j). In the Bagnoli offshore, unit S7 displays chaotic and convoluted/slumped reflectors (Figure 5f). This feature may be associated with gravity flows and landslides triggered by intense volcano-and seismotectonic activity, resulting in faulting, fracturing, and paleoliquefaction phenomena .
A large (over 22 km 2 ) abrasion surface developed on the shelf above S7 and older sequences (map inset in Figure  5) and has been correlated to a truncation contact exposed on-land (Figure 7k) that represents the proximal segment of H8 due to later subaerial exposure. It is worth noting that this large and preserved abrasion surface almost entirely wiped away the abrasion surface (H6) associated with the previous uplift phase, which is only preserved on the eastern part of the shelf (Figure 5h), with its corresponding infralittoral prograding wedge (Hernandez-Molina et al., 2000).
Seismic facies of overlying unit S8, including strongly reflective layers, are interpreted as marine sediments embedding large amounts of eruptive products of Epoch 3b (Figures 5i and 6c), emplaced during the Averno 2-Solfatara (4.28 ka), Astroni (4.2 ka), Capo Miseno and Nisida eruptions (3.9-3.7 ka; Figures 5j,k and 7g,h,l). Pyroclastic products associated with volcanism younger than Solfatara eruption overlie abrasion surface H8 (Figure 7k). In contrast, the sub-transparent facies of unit S9, calibrated by shallow coring by Sacchi et al. (2014), correlate to the centuries-long quiescence period between Epoch 3b and the Monte Nuovo eruption (1538 CE). Finally, the uppermost unit S10 reaches a maximum thickness of ca. 1.5 m in the Epitaffio valley (Figure 1a), but elsewhere it is only represented by a thin layer, correlated to Monte Nuovo tephra ( Figures  5k and 7l).

| Seismic attributes and grain size comparison
To strengthen the proposed correlation, we compared the RMS (Root Mean Square) amplitude attribute of seismic units S2 to S5 to the grain size data of correlative sediments at La Starza provided by Isaia et al. (2019). The RMS amplitude attribute analysis is the post-stack evaluation of signal strength regardless of its polarity, and its visualisation helps to identify the acoustic impedance contrast. In general, several factors such as cementation, change in bed thickness, porosity, and alteration might control the seismic amplitude. However, considering the relatively more undisturbed depositional environment of the caldera collar and the characteristics of marine sediments at La Starza, generally consisting of loose sediments, we assume that the grain size and the sorting play the primary control.
Overall, a good fit exists between the RMS amplitude and sorting and grain size for transparent and reflective facies (Figure 8). The two highest amplitude horizons correspond to poorly sorted, generally medium to coarse layers, whereas the transparent facies correspond to well-sorted, fine sands.
From bottom to top, the marine sediments of intervals b and c are constituted by fine to medium sand separated by angular unconformity H3, corresponding to a distinct emergence phase on-land (Figure 8a,b). Sediments of interval b show a coarsening upward trend associated with a sorting decrease in correspondence of channelised beds, which peaks at H3 (Figure 8c-e). By contrast, interval c and overlying marine sands d show a fining upward trend and a good sorting until interval e, corresponding to a pedogenised horizon (paleosol B, marker level #2, Figure  8a). The development of paleosol B indicates a subaerial exposure. It is represented in seismic images by the strongly reflective subunit S4c with a high RMS value, which is interpreted as the result of the uplift phase of the central sector and consequent sediment erosion and transport of poorly sorted material offshore.
The well-sorted fine sands (interval f at La Starza) corresponds to the transparent, low RMS amplitude of subunit S5a. The fine sands pass upward to less well-sorted and progressively coarser sands (interval g1), deposited in foreshore settings and topped by disconformity H5ab. Above this surface, coarse volcaniclastic deposits associated with a peak of poor sorting occur (Figure 8d). The coarsening in interval g1 is reflected in seismic images by high amplitude reflections in subunit S5b.

| Sedimentation rates
We estimated the time evolution of sedimentation rates in the CFc offshore sector based on the proposed stratigraphic reconstruction. This evaluation has to be considered a first approximation of the accumulation rate and does not include other factors such as sediment compaction. Subsequently, we compared our results with the sealevel change curve calculated by Lambeck et al. (2011) to qualitatively evaluate the contribution to sedimentation of Sea Level Rise (SLR; Figure 9a) and with the erupted magma rates (Figure 9c, Bevilacqua et al., 2016). A rapid SLR has occurred between ca. 15 and 7 ka but with significant variations in SLR rates, with peaks in velocity close to 18 mm/yr (Figure 9a; Lambeck et al., 2011). On the contrary, during the last ca. 7 ka, the sea level rose at a constantly decreasing rate (Figure 9a).
We calculated sedimentation rates for the offshore (mm/yr; green dashed line Figure 9b) and La Starza sections where data are available (blue dashed line in Figure  9b). Age uncertainty of data points representing individual rate estimates is indicated by horizontal error bars (details in the supplementary material). In addition, the curve of sedimentation rates for the last 15 kyr was drawn, joining the points with a curve instead of straight lines to enhance the interpretation. Data points are placed at the end of the reference interval.
Following the peak of the NYT eruption (and early Epoch 1 tuff eruptions), the curve shows a millennialscale decreasing trend inverted at ca. 12 ka. This peak (ca. 15 mm/yr) is related to the deposition of PP tephra, emplaced during the first large eruption following the NYT eruption (Figure 9b,c), probably with the contribution of additional eruptions (e.g., the Soccavo 1 eruption).
After that, a sudden drop in both erupted volume and sedimentation rate is observed, with sedimentation rates falling to ca. 5 mm/yr, a value broadly maintained until the start of Epoch 2 at 9.6 ka (Figure 9b,c). The gentle increase observed at around 10.5 ka, and the peaks at 9.1 ka both on-land and offshore are interpreted as related to post-PP doming and Epoch 2 volcanism/ doming contributions. The peaks in sedimentation rates offshore occur at ca. 9.4 ka and largely (ca. 26 mm/yr) at ca. 9.1 ka (Figure 9b). Although minor volumes have been erupted in Epoch 2 compared with the ca. 12 ka peak, these peaks reasonably reflect the occurrence of eruptions more proximal to the gulf, such as Baia-Fondi di Baia, Costa San Domenico and Pigna San Nicola, and also to an excess of sediment supply related to the doming that culminated at the end of Epoch 2 (Figure 9b,c). In this time frame, a possible contribution of SLR-related erosion (Schwartz, 1967) is not excluded because of the ca. 18mm/yrs peak in SLR (Figure 9a). Instead, a low (1.4 mm/yr) sedimentation rate is estimated at La Starza for the same time interval when this sector emerged and a paleosol developed (Figure 9b). This discrepancy is probably related to the erosion of coastal areas and the excess of sediment input in the offshore sector. On the other hand, a good correspondence exists between the rates of interval f (La Starza) and S5a (9.1-5.86 ka), with overlapping values around 1.6 mm/yr (Figure 9b). With the Epoch 2 uplift phase, the gulf became a relatively more proximal basin, but the ensuing long quiescence, joint to the subsidence, and slowing SLR rates led to lower sedimentation rates.
After a 2.7 ka-long quiescence, volcanic activity resumed with the onset of Epoch 3, resulting in increased sediment production (Figure 9b,c), with a constantly declining SLR rate (Figure 9a). Sedimentation rates, offshore and on-land, and erupted volume rates follow an increasing trend peaking at the end of Epoch 3a when the Agnano-Monte Spina eruption occurred (Figure 9b). This event is the largest-volume eruption in the last 10 ka and is associated with a sedimentation peak at near 30 mm/yr in the caldera collar. The following rate point is related to the sedimentation during the short-lived volcanic quiescence and rapid subsidence-uplift related to the Pozzuoli Unit. The surplus of sediment production (Figure 9b) is related to the intense ground deformation , resulting in the deposition of material in the caldera depocentres (Figures 4 and 5f). The following sedimentation peak (ca. 25 mm/yr) couples with Epoch 3b volcanism at ca. 4.3 ka, and specifically, it may be associated with products of the Averno-Solfatara and Astroni eruptive sequence (Figure 9c).
The lack of volcanic input after Epoch 3b is mirrored in a minimal sedimentation rate (ca. 1.5 mm/yr), which is only interrupted by the peak in correspondence of the Monte Nuovo eruption. We remark that the calculated sedimentary rate for unit S9, ca. 1.5 mm/yr is comparable with the one calculated for other quiescence phases (i.e., interval f and subunit S5a), as shown in Figure 9b. In summary, an excellent agreement exists between erupted volumes, sedimentation rates and resulting collar thickness, highlighting the coupling between sedimentation and volcanic activity (Figure 9b,d).

| Improving the stratigraphic position of coastal eruptions
The presented study is an up to date stratigraphic reconstruction which may allow us to improve the stratigraphic position of some eruptions along the coastline, which mainly produced tuff cones and rings such as Nisida, Capo Miseno, Nisida bank, Porto Miseno and Bacoli (Di Renzo et al., 2011;Fedele et al., 2011;Insinga et al., 2006;Smith et al., 2011). This feature has been debated in the literature, and analysis of seismic reflection profiles recently allowed Steinmann et al. (2018) better to define the relative chronology attribution of several vents and also to provide an updated erupted volume estimation that in some cases is significantly increased compared with previous estimations (Bevilacqua et al., , 2016. Here, we present a further improvement to constrain the stratigraphic position of volcanic events using our fine-tuned offshore-onland stratigraphic framework. On the eastern side of the gulf, along the Bagnoli valley, we suggest that the anticlinal structure ( Figure  4), located offshore of Bagnoli, represents a locally deformed zone related to a small volume intrusion, in agreement with Steinmann et al. (2016). A bell-like shape is well visible on ca. E-W seismic profiles ( Figures  4 and 10c,e). On the other end, the deformation can also be observed on ca. N-S lines (Figure 10f-g), including the GeoB08-035 (Figure 10g), which turns into a ca. E-W direction in its northern end (Figure 10a), allowing us to see the last clear evidence of deformation. With this in mind, considering the at least 800 m long NNE-SSW alignment of the crests and ca. 600 m width, as visible in the time-structure map (Figure 10a), we interpret this structure as a forced fold sensu Jackson et al. (2013), related to a shallow magmatic intrusion (Figure 10b). The maximum amount of displacement is ca. 12-15 m (using 1650 m/s average velocity; Figure 10c-g). It is interesting to note the occurrence of fluid uprise on the crest of the forced fold (Figure 10c,f). Such as already suggested by some authors (Natale et al., 2018;Steinmann et al., 2016), the depth of this intrusion might be below the data penetration limit (300-350 ms). Approximating the volume underlying the forced fold as a triaxial hemiellipsoid, the minimum displaced volume is ca. 3 mm 3 . Moreover, the downlap terminations (Figure 5f) of the S7 unit constrain the intrusion timing at around 4.4-4.3 ka, corresponding to the beginning of Epoch 3b similarly to other lava domes (Accademia dome 4.3 ka, Olibano dome 4.3 ka, Solfatara dome 4.4-4.28 ka) emplaced close to the coastline (Bevilacqua et al., 2016;Isaia et al., 2015;Smith et al., 2011;Figure 1a). In support of this comparison, the estimated volume is in the same order of magnitude as these exposed domes (Bevilacqua et al., 2016). Another hypothesis considers such structure as a fold associated with a reverse fault. However, the vicinity of other lava domes of the same age, and the lack of evidence of such faults in the whole seismic dataset, lead us to consider the first interpretation as the most likely.
The two tuff cones of Nisida and Banco di Nisida similarly occur on the eastern sector within the caldera ring fault zone (Figure 1a). The Island of Nisida is the subaerial emergence of a larger tuff cone complex which was already included by Steinmann et al. (2018) in the volcanic Epoch 3 in agreement with the available literature age determinations (i.a., Di Renzo et al., 2011;3.98 ± 0.53 ka for Nisida Island). The arrangement of the different volcanic mounds constituting the underwater complex suggests that volcanic activity occurred broadly aligned N-S, with the northward migration of the eruptive activity ( Figure  10h). Furthermore, our reconstruction allows us to define more precisely the stratigraphic position of the Nisida eruption within Epoch 3b since Nisida tephra is above S7 unit (Pozzuoli Unit, 4.4 ka) and within the upper part of unit S8, covered by the top of subunit S8c (Figures 4a and  5j), resulting in a supposed age of 3.7 ka.
Banco di Nisida is a large tuff cone emplaced above the erosive head of the Canyon Dhorn, which drained sediments up to unit S2. Sediment drainage was halted from the volcanic bank eruption (Figure 10h). The stratigraphic attribution of the cone to the upper part of unit S2 allowed to constrain its age to after the PP eruption (12 ka) and before unit S3 (10.5 ka), collocating it in the late-stage Epoch 1 (Figure 10h) (e.g., Di Vito et al., 1999).
The western coast of the gulf is scattered with remnants of tuff cones and tuff rings. The stratigraphic position of these eruptions is non-trivial due to the lack of clear stratigraphic ties and the distance from the inner caldera basin. In correspondence with Porto Miseno, unit S2 has an anomalous thickness increase that we associate to the presence of the submerged flanks of the Bacoli tuff cone and Porto Miseno tuff ring, separated by an unconformable contact (Figure 10i). This local increase in thickness is also observed by Steinmann et al. (2018) for their unit M2, which comprises our unit S2. The eruptions that built these two monogenetic edifices occurred short after PP (12 ka), placing these lithosomes at the end of Epoch 1, as hypothesised by Smith et al. (2011).

F I G U R E 1 0 Legend on next page
On the same line, we observe the main architectural elements (Bischoff et al., 2021) of the Capo Miseno tuff cone (Figure 10i), which is, however, partially reworked moving towards SE. This eruption probably occurred in the late stage of Epoch 3b, almost at the same time as Nisida to the eastern part of the bay. The mutual relationships between these two lithosomes are not easy to disclose, but seismic horizons correlation suggests Nisida tephra be slightly younger than Capo Miseno. However, given the vertical resolution of the seismic dataset, the contemporaneity cannot be ruled out, as is a common feature in several caldera settings (e.g., Blong, 1994;Pistolesi et al., 2016).

| Insights into caldera deformation in the last 15 kyr
The integration of the reconstructed on-land stratigraphy with the offshore seismic stratigraphic architecture provides insights into the relationships between ground deformation, sea-level rise and sedimentary evolution. This paragraph qualitatively summarises the ground deformation history in the last 15 kyr. With the occurrence of the first post-caldera uplift at 10.5 ka, the deformation broadly followed a cyclic pattern between uplift and subsidence phases, broadly centred near Pozzuoli, likewise today. This cyclicity, among other features, is recorded offshore by alternating regressive and transgressive stacking patterns, respectively. However, in the long term, the uplift prevails, as testified by exhumed marine sequences.
To better highlight the impact of the ground deformation on the depositional system, we used lateral thickness changes, reflector terminations, and geometry (convergence or divergence) as indicators of volcano-tectonics perturbations at the basin scale. Finally, we summarise the ground evolution in a stratigraphic sketch showing periods of quiescence, doming or stability with a different colour referenced to individual seismic units/subunits and corresponding La Starza intervals (Figure 11a).
Because seismic units S0, S1 and S2a do not show significant lateral thickness variations, they form a predoming sequence (Figure 11a). Thus, they are considered sea-level-controlled, aggrading sequences deposited without significant/detectable volcano-tectonic deformation. On the contrary, in younger seismic units (S2b-S10), the depositional system adapts to deformation, resulting in various stratal termination and depositional features, such as onlaps and angular unconformities, erosive truncations and gravity flows (Figure 11b).
Following the PP eruption (12 ka), the system records a subsidence localised on-land. The landward thickness increase of subunits S2b and S2c suggests a subsidence stage that deepened the central caldera sector, as much as 15-20 m, between ca. 11.9 and 10.5 ka when interval b formed (Figures 8a and 11a,b). Evidence of marine ingression following the PP eruption is reported at different on-land locations, including La Starza, San Vito Plain and several drill holes and excavations (i.a., Di Isaia et al., 2019;Sbrana et al., 2015).
The subsidence of the caldera central sector started at ca. 10.5 ka has been inverted by an uplift event, testified by the unconformity H3, which evolved to an abrasion surface at La Starza. This uplift event, amounting up to 30-35 m, is also corroborated by the occurrence of gravity flows and a shift of the depocentres, with the overlying unit S3 displaying an increasing thickness from the dome towards the collar area (Figures 2, 3, 5a,b, and 11b).
Between 9.6 and 9.1 ka, marine sediments of interval d were deposited at La Starza above the lacustrine beds of sequence c (Figure 6a), suggesting that the rate of SLR (Figure 9a) exceeded previous ground uplift of S4a, corroborated by ca. 10 m of subsidence (Figure 9d, Isaia et al., 2019), while the transparent subunit S4b, with probably indicating a stasis in volcaniclastic input. The following strongly reflective S4c. is associated with the forced marine regression, concurrently with the emersion of La Starza unit at around ca. 9.1 ka and no later than 8.6 ka . The emergence of La Starza was likely accompanied by the formation of an abrasion surface that isolated the cliff, with not less than 20-25 m of uplift at La Starza (Figures 9d and 11a).
Nonetheless, between ca. 9.1 and 5.9 ka, when the sea level rose from −30 to −5 m (Figure 9a), the sedimentary facies at La Starza suggests a significant paleobathymetric depth increase (down to 60/80 m b.s.l. according to Isaia et al., 2019; Figures 8a,b and 9d). This sequence event at F I G U R E 1 0 (a) Time-structure map of H6 horizon above the Bagnoli bulge. Thin grey lines indicate 2D seismic lines used to reconstruct the time-structure map. Survey map location in Figure 1a. (b) Sketch of an intrusion-related forced fold (modified after Jackson et al., 2013). (c-e) ca. E-W line Seistec-801, GeoB08-065, Seistec-905 showing the geometry of the forced fold close across sections orthogonal to the fold axis. (f-g) Interpreted ca. N-S line GeoB08-037 and GeoB08-035 pointing out the deformation parallel to the fold axis. (h) Interpreted N-S trending line GeoB08-036 running along the coast of Nisida. Several features are present, like erosive truncations, inner ring fault, the two tuff complexes, and their relationships with the stratigraphic frame. (i) Interpreted N-S line GeoB08-094 running along the western coast sharpen the stratigraphic relationships between Bacoli and Porto Miseno volcanoes of the Epoch 2 and Capo Miseno of Epoch 3b. BAC, Bacoli; BN, Banco di Nisida; CM, Capo Miseno; NIS, Nisida; PM, Porto Miseno La Starza corresponds to the full transparent facies of subunit S5a, which highlights the lack of primary pyroclastic input (Figure 8e). Thus, with the development of an uplifted dome, the Pozzuoli Gulf has become a more proximal basin than the pre-doming conditions and is rather constant in the last 6 kyr. Hence, in this timeframe, we use the position of the on-lap terminations (labelled from 1 to 6 in Figure 11b Cruden et al., 2017) deformation reconstructed at La Starza by Isaia et al. (2019) (Figure 9d).
Horizons H5a and H5b are truncated in the shelf area (point 1 in Figure 11b); however, the same marine unit is exposed on-land at La Starza (interval f, g1; Figure  5d-h). On the other hand, unit S6 onlaps above its basal unconformity within the present slope (point 2, Figure  11b) due to marine regression caused by ground uplift. The La Starza unit remained in subaerial conditions from 5.25 to 4.4 ka, following the uplift phase that continued during deposition of unit S6 (Epoch 3a), which further developed an abrasion surface preserved from subsequent erosion in the eastern sector (Figure 5h). The onlap at the base of unit S7 (point 3, Figure 11b) moves over 500 m landward, accompanying the marine transgression of the Pozzuoli Unit (Isaia et al., 2009), which deposited marine-transitional sediments on-land (Figure 7i,j).
A large abrasion surface is associated with the subsequent regression phase. This structure is a distinctive morphological feature in the Gulf, which eroded the previous units and was associated with a well-developed prograding wedge (Sacchi et al., 2014; Figure 7k). High-angle faults, which imply negligible horizontal extension, do not dislocate the abrasion surface, thus are not to be linked to the ground resurgence but likely to volcano-tectonic collapses (e.g., Ruch et al., 2012). The occurrence of gravity flows also testifies the regressive trend flows within unit S8 ( Figure  11b). This regression is associated with the uplift phase that accompanies the Epoch 3b activity , with the corresponding onlap at the base of subunit S8b (point 4, Figure 11b) located 400 m seaward, recording the Epoch 3b uplift, which recovered and exceeded the previous shortlived subsidence (Isaia et al., 2009. The relative stability of point 5, and later landward shift of point 6, indicating a subsiding trend, characterise the end of Epoch 3b (Figure 11a,b). In summary, the aggregate vertical displacement since 10.5 ka, considering subunit S2a as the last predoming unit, is at least 100 m between the caldera collar and the shelf. Most of this aggregate value is due to Epoch 3 uplifts (Figure 9d, Isaia et al., 2019).
Nonetheless, not all the ground deformation visible in the infill of the last 15 kyr is related to the resurgent dome or the ring faults. In fact, in addition to the Bagnoli intrusion, a significant ground deformation occurred in the area east of Porto Miseno (Figures 1a and 3). The Punta Pennata structure is characterised by an asymmetric ca. 700 m N-S length with a tapered southern end and a subvertical, fault-bound northern edge (Figure 11c). Here a ca. 20 m localised uplift occurs in the footwall of the inner ring fault, for over 1500 m along the E-W direction (Figure 11d; bysmalith, sensu Cruden et al., 2017, and references therein). This structure was active independently from the corresponding ring fault, as the depth of reference horizons in the hangingwall along strike is not downthrown, as visible from the time-structure map of Figure  11d. Using our detailed stratigraphic framework, we track the temporal evolution of the deformation occurring here (Figure 11c). A first intrusion occurred around 10.5 ka, with the sharp displacement from subunit S2a to unit S3. Erosive unconformity probably formed concurrently with the uplift during the 9.6-9.1 ka time interval as testified by growth strata in unit S3. Weak continuous activity and differential compaction might explain thickness variation within units S5a and S6, whereas on-lap terminations within S8 unit suggest a renewed activation during Epoch 3b (Figure 11c). In general, volcanic activity in the last 15 ka has often clustered along the ring faults (Bevilacqua et al., 2016) and is already being considered in the recent work about hazard assessment (Bevilacqua et al., 2017). Moreover, continental deep drill holes (e.g., AGIP-ENEL joint-venture) have also encountered peripheral subvolcanic bodies at different depths (Rosi & Sbrana, 1983).

| CONCLUSIONS
This paper provides an integrated on-land-offshore stratigraphy of the last 15 kyr for the Campi Flegrei caldera. The interpretation of 2D reflection seismic profiles allowed us to reconstruct a high-resolution stratigraphic framework of the offshore post-NYT caldera sedimentary record and better understand the complex interplay between the caldera floor deformation (uplift and subsidence), sea-level change and volcanic activity. We correlated the main UBSUs with the marine-transitional sequence exposed on-land (La Starza unit). An excellent agreement exists between the seismic signature of the individual seismostratigraphic units and the sedimentological features of the exposed counterpart.
We identified eleven units above the NYT volcanic unit tested and fine-tuned by calculating the sedimentation rates for each unit/subunit and comparing them with the sedimentation rate of the corresponding intervals of La Starza sequence, calibrated with the sealevel rise. Results indicate that the sedimentation rate increases during the syn-eruptive volcaniclastic supply and uplift as more sediments are eroded from the dome apical region and transported seaward. This feature indicates that the Campi Flegrei volcanic activity strongly controlled the sedimentation within the Pozzuoli Gulf, as supported by comparing sedimentation rates and estimations of erupted magma volumes. Furthermore, it is worth noting that the sedimentary system reacts and adapts to the deformation of the caldera floor, which systematically produces significant variations of the accommodation space and the localisation of the depocentre over time.
Following the large-magnitude Pomici Principali eruption (12 ka), the sedimentation recorded the deformation affecting the caldera floor, resulting in the development of unconformities, lateral thickness variations and variable onlap stacking patterns. Centuries-long subsidence was interrupted with the first uplift at around 10.5 ka. The volcanic activity was generally coupled with the uplift, even though minor unrest phenomena not followed by eruptions cannot be excluded during the CF history of the last 15 kyr. The configuration of seismic stacking patterns and onlap terminations highlight the considerable uplift of the central part of the caldera that accompanied the Epoch 3 volcanism. This ground deformation stage was characterised by alternating uplift and subsidence phases, such as the transgressive event of the Pozzuoli Unit, following the large-magnitude Agnano-Monte Spina eruption (4.55 ka). The subsequent volcanic activity has been accompanied by important ground uplift well-recorded by erosive surfaces. In addition, the seismic analysis allowed us to frame the relative stratigraphic position of coastal volcanic eruptions and offshore intrusions, which contributes to refine an integrated postcaldera stratigraphic evolution.
This study underlines the benefit of a joint on-landoffshore approach to reconstruct the whole caldera stratigraphic framework and more comprehensively assess the long-term volcano-tectonic evolution, which is of the utmost importance in a highly inhabited, restless volcanic area such as Campi Flegrei caldera.