RSL reconstruction

The harbour basin is the outcome of an environmental history in which long-term sea-level rise has played a key role. A total of 31 RSL index points has been used to frame the Holocene sea-level evolution of the eastern Ligurian Sea (Fig. 1). At ~8550 BC, the oldest index points place the RSL at ~35 m below current Mean Sea Level (MSL). Younger index points indicate that RSL rose rapidly until ~5050 BC followed by a slowdown in the rising rates. At ~4050 BC, multiple RSL index points constrain the RSL to ~−5 m MSL. Since this period, index points delineate a significant reduction in the rising rates that become minimal during the last ~4000 years, when the total RSL variation was within 1.5 m MSL. The reconstructed RSL history reflects the pattern observed in the mid to northern portion of the western Mediterranean36,37. The significant slowdown in rising rates after ~5550 BC is consistent with the final phase of North American deglaciation38 while the further decrease in rising rates is related to the progressive reduction in glacial meltwater inputs that were minimal during the last ~4000 years39,40.

Figure 1 Relative sea-level reconstruction of the eastern Ligurian Sea. The RSL history is based on 31 index points deriving from lagoonal sediment archives of the Arno and Versilia coastal plains and fossil Lithophyllum byssoides rims from Northern Corsica. The blue boxes represent index points from lagoons and salt marshes. The black boxes are L. bissoides-derived index points. The dimensions of the boxes denote the 2 s altitudinal and chronological errors associated with each index point. The map is an original document drawn using Adobe Illustrator CS5 (http://www.adobe.com/fr/products/illustrator.html). Full size image

Drilling the harbour basin

The terrestrial, marine and freshwater biological indicators used to reconstruct palaeoenvironmental evolution of Portus Pisanus were extracted from an 890-cm sedimentary core (PP3, 43°35′55.33″N, 10°21′41.71″E; +2 m MSL; Fig. 2) drilled ~5 km from the sea on the southern portion of the modern Arno River alluvial-coastal plain, close to Livorno and the Pisa Hills. During ancient times, this part of the alluvial plain, located ~10 km south of Pisa, was fed by a former branch of the Arno River41 here reported as the Calambrone River. According to previous studies15,18,20, the sedimentary core was taken from the former harbour area, active from the Archaic Etruscan (6th–5th centuries BC) to the Late Roman periods (5th century AD). In the Middle Ages, the port shifted westwards, as testified by the building of the fortified harbour basin of Leghorn including medieval towers dated from 1300 and 1400 AD (Fig. 2), and then again further westwards in modern times, consistent with the progradation dynamics of the delta. The core was recovered behind the innermost beach-ridge of the Arno Delta complex, where a meter-thick back-barrier succession occurs, recording the development and prolonged persistence of a wide lagoon basin during the Late Holocene. Back-barrier deposits, mainly represented by fossiliferous clay-silt interbedded with sandy layers, are overlain by a progradational suite of coastal-alluvial facies, deposited during the recent phase of decelerating sea-level rise.

Figure 2 Study area and location of the archaeological site. (A) Geomorphological map of the study area (modified from CARG Regione Toscana). 1 current beach; 2 shallow swale; 3 wetland; 4 beach ridge, superimposed dunes; 5 alluvial plain; 6 residual relief; 7 Livorno urban/industrial area; 8 mountains and hills; dashed line: 17th century AD coastline; dotted line: 12th century AD coastline; arrow: current drift. (B) Photograph of the archaeological site. Full size image

The chronology of the core PP3 is based on nine accelerator mass spectrometry radiocarbon (14C) dates (Fig. 3). Dated samples (short-lived terrestrial samples: seeds, small leaves of annual plants) were calibrated [2-sigma (σ) calibrations, 95% of probability] using Calib-Rev 7.1 with IntCal13. According to the 14C chronology, the core covers the last 8000 years (Fig. 3). The age model (Fig. 3) was calculated using Xl-Stat2017 and Calib-Rev 7.1. The dates obtained in between each 14C dating are modeled and therefore are liable to mask some of the temporal variability in the depositional patterns. The calculated model displays an average 2-sigma range of 50 years (P < 0.001) for the whole sequence. All the calibrated ages are shown/discussed as BC/AD to fit with the archaeological-historical data and are presented at the 2-sigma range (95% probability). The two scales, BC-AD and BP, are both displayed on each figure. While the average chronological resolution of the core stratigraphy is 9 years per cm−1 (1.11 mm per yr−1), a homogeneity test (Monte Carlo simulation, standard test, P value < 0.001) suggests two abrupt changes in the sedimentation rate at 740 cm depth (4300 ± 70 BC) and 290 cm depth (200 ± 30 BC).

Figure 3 Details of Portus Pisanus basin in North Tuscany, Italy. The lithology of the core in the basin area, with the influence of marine components, is reported according to depth. The main sedimentary environments are plotted on a linear depth-scale. The radiocarbon dates are depicted as 2σ calibrations (95% of probability). The age model is superimposed on the 2σ calibration curve, and a linear-model was also added showing a theoretical continuous sedimentation rate. The timescale is shown as BP and BC-AD. Full size image

Biological proxies in the harbour basin

Terrestrial data retrieved from Portus Pisanus’ harbour basin were analysed using a cluster analysis (descending type). Each cluster was summed to generate pollen-derived vegetation patterns and assigned to a potential location, from the intertidal zone to the hinterland (Fig. 4), referring to modern patches of vegetation along the Ligurian coast (local data and Vegetation Prodrome42). To ascertain the ordination of terrestrial data according to the “sea” factor, a second cluster analysis (descending type; Fig. 5A) was calculated using the vegetation patterns and the marine proxies [dinoflagellate cysts and marine components (foraminifera, marine bivalves, debris of Posidonia oceanica)]. Three vegetation communities (backshore scrubs, shrubland, and coastal pine-oak woodland) are linked to a marine influence (from the supratidal zone to the lower coastal zone; Figs 4–5A) whereas the other communities (mixed oak forest, wet meadow, fen trees, freshwater plants) are related to fluvial inputs from the Calambrone River (coastal alluvial zone and marsh-swamp zone). Cross-correlations (vegetation patterns versus marine proxies; Fig. 5B) also indicate that seawater has influenced proximal vegetation patterns (positive correlations on the null lag score: Lag 0 = 0.665, Lag 0 = 0.547 and Lag 0 = 0.307, P = 0.05) around the basin. The vegetation group “warm woodland” is set apart as this cluster is mainly related to a third influence, agro-pastoral activities (Fig. 5A) that are observed in the area after 3350 ± 90 BC. Agriculture is composed of cereals (Poaceae cerealia), olive trees (Olea europaea), common grape vines (Vitis vinifera), and other trees (Prunus). The associated anthropogenic indicators are common weeds (Centaurea, Plantago and Rumex).

Figure 4 Pollen-based ecological clusters from the harbour basin for the last 8000 years. A cluster analysis (paired group as algorithm, Rho as similarity measure) was used to define the ecological assemblages. Each cluster was summed to create pollen-derived vegetation patterns. The potential location of each cluster, from the intertidal zone to the hinterland, is indicated. Full size image

Figure 5 Environmental-based clusters from the harbour basin for the last 8000 years. (A) A cluster analysis (paired group as algorithm, Correlation as similarity measure) was used to define the environmental assemblages (marine versus fluvial influence). (B) The two cross-correlograms (P = 0.05) depict the marine influence on ecosystems around the basin. Full size image

Pre-harbour facies (sensu Marriner and Morhange)

Marine versus fluvial influence since 8000 years is represented by the importance of marine indicators (dinoflagellate cysts and marine components) and supratidal-intertidal scrubs in and around the pre-harbour (Fig. 6). A principal components analysis (PCA) was run to test the ordination of ecosystems by assessing major changes in the area including vegetation patterns, dinoflagellate cysts and marine components (Fig. 7). Environmental dynamics (marine versus fluvial influence) in the basin is indicated by the axis-1 of a principal components analysis (PCA-Axis1). The PCA-Axis1 (61% of total variance) is positively loaded by vegetation patterns indicative of a saline-xeric environment, dinoflagellate cysts and marine components. The negative scores correspond to freshwater vegetation types (Fig. 7). The PCA-Axis1 reflects the ecological erosion of wetlands by the intrusion of seawater into the freshwater-fed plains, raising salinity in the hinterland, with land fragmentation and salt-water intrusion into the groundwater table, in and around the basin. The PCA-Axis1 can be considered as a proxy for marine ingression in/around the pre-harbour (with a main physical impact and several secondary influences such as salt spray and salinization)43.

Figure 6 Reconstructed marine influence in Portus Pisanus during the last 8000 years. The marine influence (components per cm−3) and the backshore scrubs (%) are displayed as a LOESS smoothing (with bootstrap and smoothing 0.05) plotted on a linear timescale (BP and BC-AD). The Posidonia oceanica debris (presence/absence) are indicated by green marks along the marine curve. The ostracods (presence/absence) are displayed as dots. Fire activity in the area is shown as charcoal concentrations (fragments per cm−3) plotted on a linear timescale. The relative sea level, displayed as MSL, is depicted for the last 8000 years. The shipwreck curve from the Mediterranean region is also plotted on a linear timescale76. The brown-shaded horizontal section indicates the harbour development and the blue-shaded horizontal section shows the period when the sea-level stabilized. Full size image

Figure 7 Reconstructed environmental dynamics in Portus Pisanus during the last 8000 years. The PCA-Axis1, plotted on a linear age scale (BP and BC-AD), is displayed as a LOESS smoothing (with bootstrap and smoothing 0.05) and a matrix plot. A boxplot was added to mark the extreme scores. The loading of each cluster is indicated below the PCA-Axis1 curve. The agro-pastoral activities (agriculture and weeds, %) are plotted on a linear age scale, and also displayed as a matrix plot. The brown-shaded horizontal section indicates the harbour phase and the blue-shaded horizontal section shows the period when sea-level stabilized. Full size image

While permanent inputs of seawater and fluvial freshwater were recorded in the pre-harbour, occurrences of higher marine influence are locally underlined by a combination of peaks in marine indicators, the occurrence of different ostracod ecological groups (shallow marine; brackish-marine and euryhaline), important increases in Posidonia oceanica debris (Fig. 6), and strong positive deviations in the PCA-Axis1 scores (Fig. 7). Our reconstruction shows two early periods of seawater influence (at 5800 ± 40–5425 ± 55 BC and 4750 ± 60–4500 ± 60 BC) when the site was a marine invagination, followed by an unstable phase when sea-level stabilized along the Tuscany coastline (Fig. 1), between 4250 ± 60 BC and 2000 ± 45 BC (Fig. 8). Potential discontinuities in environmental dynamics were assessed using a homogeneity test on the PCA-Axis1 (Monte Carlo simulation, Pettitt and Buishand tests). The outcome indicates that the environmental dynamics are not uniform, underlining a major break around 2000 ± 45 BC (P value < 0.001). A second homogeneity test, only applied to the period 3350–6050 BC, highlights a second important break around 4250 ± 60 BC (P value < 0.001, Monte Carlo simulation, Pettitt and Buishand tests), when the pre-harbour evolved into a leaky lagoon. During this period (4250 ± 60–2000 ± 45 BC), the last main peak in the PCA-Axis1 corresponds to a wave-dominated delta and occurred within the chronological interval of the 4.2 ka BP event44, suggesting the potential role of climate in influencing the pre-harbour’s evolution. A later phase was recorded at 1250 ± 40–850 ± 40 BC when the pre-harbour evolved into a delta plain, corresponding to the 3.2 ka BP event45,46.

Figure 8 Geographic maps showing three evolutionary stages of the lagoon, in relation to Portus Pisanus as documented by historical sources and archaeological data. Maps were produced by integrating stratigraphic (cores and trenches) and geomorphological data (Pranzini, 2007) with historical cartography. Dots indicate cores used to draw the maps (key cores are highlighted by red dots). The modern shoreline is depicted on each map for reference. The grey arrow indicates the direction of the predominant wind (Libeccio). (A) Roman period - a wide lagoon basin, hosting Portus Pisanus as mentioned in literary sources; (B) late Middle Ages - the accretion of arcuate beach ridges, belonging to the Arno Delta strandplain, led to an increase in the degree of confinement of the lagoon basin. Construction of the maritime harbour of Livorno in a seaward position with respect to the lagoon; (C) 17th century AD - the rapid accretion of strongly arcuate sets of beach ridges led to the siltation of the lagoon that was transformed into a wetland, physically detached from the Ligurian Sea. Portus Pisanus abandonment and expansion of the fortified maritime harbour of Livorno. Full size image

Among the main periods of marine influence in the harbour basin, before the establishment of Portus Pisanus, the period 4250 ± 60–2050 ± 45 BC was the hinge phase, culminating in a wave-dominated delta. These periods, characterized by inland salt intrusion, significantly affected the agricultural productivity of this coastal area (Fig. 7). Prolonged marine inundation appears to have led to the salinization of agriculturally productive soils resulting in diminished output for long periods of time (Fig. 7).

Portus Pisanus

Because the area has been frequented by archaic ships since at least the 6th–5th centuries BC, the first phases of navigation had to be carried out on a delta plain characterized by wetlands34,35. According to our reconstruction (Figs 6–7), marine influence increased in the harbour basin after 200 ± 30 BC, when a naturally protected lagoon developed and hosted Portus Pisanus up to the 5th century AD (according to archaeological evidence14,34,35). During this period, a first peak in charcoal fragments is recorded at 180 ± 30 BC. A second inflection in charcoal fragments occurred at 550 ± 25 AD and is synchronous with a major fall in agro-pastoral activities (Fig. 7). From 1000 ± 20 to 1200 ± 20 AD, a first decrease in marine influence was documented before the last major marine phase at 1300 ± 20 AD. The decline of the protected lagoon started at 1350 ± 15 AD and ended at 1500 ± 10 AD, when the basin evolved into a coastal lake, concomitant with the development of agriculture, then to a floodplain (1700 ± 10 AD) and finally a soil atop a fluvial plain (1850 ± 5 AD).

Stratigraphic data from several cores (Fig. 8) and archaeological trenches (Fig. 9) undertaken at Santo Stefano ai Lupi (Fig. 2), along with prominent geomorphological features (e.g. outcropping beach ridges and residual reliefs), were used to produce three maps, which illustrate the landscape evolution of the southern portion of the Arno alluvial-coastal plain, between the Roman period and the Modern age (Fig. 8). The extension and the environmental characteristics of Portus Pisanus are based on facies correlations (chronological framework based on 14C dates and archaeological data). The core PP3 (and also PP1) formed the type stratigraphy of the basin (Fig. 8).