Sea-Level Rise and Shoreline Changes Along an Open Sandy Coast: Case Study of Gulf of Taranto, Italy

: The dynamics of the sandy coast between Castellaneta and Taranto (Southern Italy) has been influenced by many natural and anthropogenic factors, resulting in significant changes in the coastal system over the last century. The interactions between vertical components of sea-level changes and horizontal components of the sedimentary budget, in combination with anthropogenic impact, have resulted in different erosion and accretion phases in the past years. Local isostatic, eustatic, and vertical tectonic movements, together with sedimentary budget changes, must be considered in order to predict the shoreline evolution and future marine submersion. In this study, all morpho-topographic data available for the Gulf of Taranto, in combination with Vertical Land Movements and sea-level rise trends, were considered by assessing the local evolution of the coastal trend as well as the future marine submersion. Based on the predicted spatial and temporal coastal changes, a new predictive model of submersion was developed to support coastal management in sea-level rise conditions over the next decades. After that, a multi-temporal mathematical model of coastal submersion was implemented in a Matlab environment. Finally, the effects of the relative sea-level rise on the coastal surface prone to submersion, according to the Intergovernmental Panel on Climate Change Assessment Reports (AR) 5 Representative Concentration Pathways (RCP) 2.6 and RCP 8.5 scenarios, were evaluated up to 2100.

On a global scale, satellite altimetry revealed the absolute sea-level variations relative to the reference ellipsoid, showing a sea-level rise of 3.3 mm/year (NASA Goddard Space Flight Center data). In the Mediterranean basin, however, tide gauge records revealed the sea-level height to the land upon which benchmarks are grounded, and a sea-level rise of 1.8 mm/year in the last century [28,29]. Sea level is expected to be 500-1400 mm higher than present by 2100 AD [21,22,27,30]. According to Church et al. [23,24] and the Intergovernmental Panel on Climate Change (IPCC) scenarios [25][26][27] during the 21st century, it is very likely that the global mean sea-level rise will exceed values observed in the last three decades, increasing up to [8][9][10][11][12][13][14][15][16] mm/year during 2081-2100.

Local Sea-Level Changes
Sea-level records in the Mediterranean over the last century show a sea-level variability depending on the location of the individual stations [28,29,31,32]. Rates vary in accordance with the length of the records, the local Vertical Land Movements (VLM), as well as the changes in atmospheric dynamics. In Antonioli et al. [33], the tide gauge data analyzed for the time span 2000-2013 at given localities along the Italian coasts, show sea-level rates ranging from 3.9 ± 0.2 mm/year (Ancona, Adriatic Sea) to 6.8 ± 0.1 mm/year (Cagliari, Tyrrhenian Sea), up to as much as 10.7 mm/year (Venice, Adriatic Sea). This is due to the decadal sea-level variability. Vecchio et al. [31] analyzed a set of long sea-level recordings (>60 years) for the northern Adriatic Sea, and estimated rates of around 1.2 mm/yr. However, at the subsiding tidal station of Venice Punta della Salute (Venice), rates reached up to 2.4 mm/year. In the Italian peninsula, as on other coasts of the Mediterranean and beyond, sea-level data show the occurrence of a continuous sea-level rise over the last decades compared to the previous centuries. Based on these rates, together with the contribution of VLM and the climatic scenarios reported in the IPCC Assessment Reports (AR) 5, and Representative Concentration Pathway (RCP) 8.5, a relative sea-level rise of approximately 0.8 m is expected by 2100 in the northern Adriatic Sea, and about 0.6 m for the Gulfs of Cagliari and Taranto. However, several authors [21,22,30,34] predict a more severe global sea-level rise of up to about 2 m by 2100, applying regional climate models [32] and a variety of statistical approaches in combination with VLM [31]. The latter value would cause a dramatic submersion of global coastal areas, and would affect human settlements, industrial and commercial facilities, archaeological and cultural heritage sites, and natural areas. Sea-level rise warnings and coastal hazards have been issued globally [35], as well as in Italy [2,[16][17][18][19][20]33,36]. In addition, human activities may produce land subsidence due to compaction or extraction of fluids from the alluvial sequences, determining local subsidence up to several mm/yr, thereby increasing local sea-level rise [2,[37][38][39][40].

Coastal Behavior in Function of Sea-Level Rise
While many studies have focused on the sea-level rise, few have described its effects along the coasts. This may be due to the difficult approach to sedimentary coastal budget and its dynamics, which only the common Bruun rule model [41] can help overcome. Currently, the Italian coast is influenced by a decrease in the sedimentary input owing to natural and anthropogenic conditions, e.g., dams and hydrological adjustments, changes in vegetation cover and climate, erosional/depositional capacity of the rivers, defense interventions along the coastal stretches. These conditions are making the coastal system changes, caused by atmosphere dynamics and rising sea level [31,33,38], difficult to assess. Thus, only a fallacious answer may explain how coastal landscapes respond to the landform re-organization during an accelerated sea-level rise. This is especially the case for those tracts, such as sedimentary coasts, that may undergo rapid morphological changes in river drainage basins in order to attain new equilibrium profiles [9,10,42,43]. Predicting a future scenario for every coast set within Valentin's milestone diagram [44] ( Figure 1) is not an easy task. Hence, many areas along the Mediterranean coasts will be susceptible to submersion by future sea-level rise and severe flooding [20,29,31].  [44] which defines the changes in the shoreline in combination with relative sea-level rise.

Objective
Projected sea-level rises could determine significant future surface loss on the low-lying coastal regions in the Mediterranean basin.
Furthermore, sea-level rise could cause permanent surface submersions in many Mediterranean coastal areas. Thus, it is important to assess the surfaces that could be submerged at any given moment of the future [2,19,33,45].
This paper proposes a new mathematical model for predicting possible submersion scenarios as the one developed for the area of the northern coast of the Gulf of Taranto in Southern Italy. It will consider the changes in both the vertical sea-level components and the horizontal shoreline components, beginning from the past articulation and behavior of the mobile coastal systems. To highlight the coastal response in function of sea-level rise, knowledge of coastal changes and decadal sea-level changes will be reviewed, while, to assess the coastal response linked to the sea-level rise, various existing models will be considered e.g., [9,41,43,46,47].
The aim of this paper is to develop a submersion model of the relative sea-level rise scenario estimated for the coastal plain on the Gulf of Taranto by 2100.

Geomorphological Setting
The Metaponto coastal plain covers most of the northern coast of the Gulf of Taranto. The plain takes its name from the ancient Magna Grecia colony of Metapontum, whose remains are still standing today. It stretches NNE to SSW from the Taranto area, and runs along the foredeep of the Bradanic Trough, up to the area of Roseto Capo Spulico, at the eastern border of the Apennine Chain ( Figure 2a).
The coastal plain is the result of the sedimentary body accumulated during the last 7 kyear (10 3 years) on the local basement. It was shaped by marine processes which produced a ravinement surface at the base of an old cliff. This coastal plain was conditioned by the interplay between the eustatic sea-level rise during the post Last Glacial Maximum (LGM) transgression and the regional uplift [48,49]. The coastal plain has been heavily supplied by siliciclastic sediments discharged by the main rivers crossing the Apennine chain (Agri, Sinni, Cavone, Bradano, and Basento Rivers). It is composed of quartz, polymineralic lithics and bioclasts [48,49] which are distributed SSW to NNE by the longshore drift.
In particular, these main phases can be recognized: • Marine Isotope Stage (MIS) 2-Last Glacial Maximum (LGM): sea-level stand at about 120 m below the present sea level and deep engravings on the coastal wedge up to the shelf break • 20-6/7 kyear: sea-level transgression with a rate of 8 mm/year and re-shaping of incised valleys along the Gulf of Taranto • from 6/7 kyear to the second half of twentieth century: sea-level rise allowing the growth of the coastal plain and beach-dune system progradation conditioned by NE longshore drift.
The Metaponto plain is bordered landward by a progradated mobile coastal system [55,56] consisting of a series of polyphasic, quasi-parallels dune belts, extending more than 1000 meters inland and characterized by an altitude ranging between 8 and 17 meters [56,57] (Figure 2b). The dune belts are characterized by an accentuated lateral continuity, and are largely colonized by the "Macchia Mediterranea", e.g., vegetation assemblages stabilizing the innermost dunes. The back-dune areas are characterized by low-lands that, in the past, hosted a lagoon, reclaimed at the beginning of the twentieth century. A primary dune belt is bordered and supplied by a medium/fine sandy beach with a gently sloping profile, currently undergoing strong erosion. The entire coastal area was in progradation up to the late 1950s. In many areas, dunes show strong evidence of significant shoreline retreat, highlighted by cliffs shaped in the dune deposits ( Figure 3). This retreat has been induced by a negative sedimentary input as a result of: (і) the anthropogenic modification of the river basin catchments with the realization of large dams and hydrological works during the second half of the twentieth century and, (іі) more recently, the conditioning of the longshore drifting due to the building of some touristic harbors and coastal defenses ( [58] and references therein). This mobile coastal system is characterized by an intermediate surf zone domain between reflective and dissipative hydrodynamic regimes, producing marked changes of nearshore sea-floor features [4,49]. Originally, the studied mobile coastal system was characterized by an intermediate morpho-hydrodynamic state, with rhythmic bars located where wave interactions cancel net bottom stress.
At present, these bars act as natural breakwaters, reducing wave energy in the proximity of the shoreline. This influences the currents that allow sediment movement along the entire littoral zone up to a seaward annual closure depth equal to 7.3 m [59].
The study area corresponds to the easternmost part of the plain, stretching between Castellaneta and Taranto. As in the rest of the plain, the mobile coastal system is characterized by well-defined beach-dune belt-back dune areas; several orders of parallel dune belts with elevations of up to about 15 m above sea level (a.s.l.), crowned by a backshore with widths ranging from 25 to 40 m and a foreshore width of up to 15 m. Four short rivers, the Lato, Lenne, Patemisco, and Tara, are mainly supplied by karst aquifers from East to West across the coastal plain without significant sediment supply to the mobile coastal system.

Tectonics and Isostasy
The Taranto area is characterized by a weak vertical tectonic uplift of the Late Pleistocene and the Holocene. This can be inferred by the elevation of the inner margin of the Marine Isotope Stage (MIS) 5.5 deposits in the Taranto area at 23 m a.s.l. and up to 40 m between Metaponto and Policoro (e.g., [52][53][54][60][61][62]).
The different altitudes of the MIS 5.5 terrace indicate decreasing values of long-term uplift rates, from West to East, along the Gulf of Taranto, from 0.39 mm/year close to the Apennine Chain, to 0.26 mm/year in the Taranto area [40,52].
To evaluate the isostatic components, the Lambeck et al. models [11,13] were applied. They consider three-layers, an elastic lithospheric thickness of 110 km, an upper mantle viscosity of 3 × 10 20 Pa·s and a lower mantle viscosity of 3 × 10 22 Pa·s, providing a rate of the Vertical Land Movement at about 0.45 mm/year in the Gulf of Taranto.
An Interferometric Synthetic Aperture Radar (InSAR) analysis [63] did not show evidence of anthropogenic-induced displacement in the Taranto area. This reveals the general tectonic stability in the NE end of the gulf and a very slight subsidence trend in the northern part of the City of Taranto. To assess the current rates of VLM along the coastal zones of the Mediterranean region, a continuous GPS (cGPS) dataset, analyzed in the frame of the SAVEMEDCOASTS Project (www.savemedcoasts.eu [64]), was used. This study focuses on the uplift velocity estimated at the cGPS station of MMET (Matera METaponto) which falls within the investigated area located in Metaponto (Matera, Southern Italy). This station is managed by the Italian Space Agency (ASI/CGS) and is integrated in the Rete Integrata Nazionale GNSS (RING) network (http://ring.gm.ingv.it/ [65]). Vertical GPS velocity has been determined in the IGb08 realization of the global ITRF08 absolute geocentric reference frame [31,64,66] providing an uplift of 0.29 ± 0.14 mm/year.

Tide and Sea-Level Trend
The current sea-level trend between Castellaneta and Taranto was obtained using sea-level data recorded by the Istituto Superiore per la Protezione e la Ricerca Ambientale (ISPRA), tidal station located at the Molo Sant'Eligio, in the Taranto harbor. Data were continuously recorded with a sampling rate of 10 min using an ultrasound transducer SIAP ± MICROS ID0710 (from 1999 to 2010, ISPRA, Taranto, Italy), and a radar sensor SIAP + MICROS TLR (from 2010 to 2020, ISPRA, Taranto, Italy). During the analysis, gaps and outliers were recognized and corrected to avoid artificial signals in the time series. The analysis highlighted a semidiurnal tidal range of 0.4 ± 0.1 m and a mean sea-level rise of 1.04 ± 0.5 mm/year in the time span 1999-2020 (20 years) ( Figure 4). The tide gauges located along the Apulian coasts, managed by ISPRA and the Autorità di Bacino della Puglia (AdBP), have shown a continuous sea-level rise in the last two decades, with similar trends for Taranto ( Figure 5).   (Table S1,  supplementary material).
From each image, the shoreline was extracted, the weather conditions at the time of acquisition were considered, and the barometric effects reduced, using the Theiler et al. method [69].

LIDAR and TLS
Laser Imaging Detection and Ranging (LIDAR) data were derived from the 2008-2009 surveys performed by the Ministero dell'Ambiente, at a resolution 2 × 2 m. In addition, Terrestrial Laser Scanner (TLS) surveys were performed at different times during the winter of 2018 by means of the Faro Focus X130 (Polo Scientifico Tecnologico Magna Grecia, Taranto, Italy), at a resolution of 2 cm. LIDAR and TLS data and Digital Elevation Models (DEMs) were combined to extract the shoreline configuration and the shore morpho-topography ( Figure 6). For the considered coastal stretches, DEMs reveal a mean shoreface slope between 2 and 5 degrees which, in response to the tide amplitude, determines a daily horizontal shoreline displacement of 12.5 ± 3 meters under a low hydrodynamic condition. A reference value to evaluate the shoreline rate changes was obtained by integrating long-term changes, detected by orthophotos and satellite images, with short-term changes detected by TLS data.

Data Processing
An analysis of the spatial and temporal changes in the Ionic littoral zone between Castellaneta and Taranto was performed by digitalizing all the shorelines via aerial photographs, orthophotos, satellite images, and LIDAR data available for that stretch of coast. They were collected for a medium-term analysis and integrated with TLS data acquired at different hours on different days. From all the collected images, the shorelines were digitalized in ArcGIS for each year. Digital Shoreline Analysis System (DSAS) tools were also used [69] in order to obtain the shoreline 1947 to 2018 changes. During the processing, uncertainty about the position of the shoreline had to be defined. The following factors were considered: For the entire stretch of the sandy coast, 345 transects (50 m spacing) were extracted perpendicularly to a reference baseline placed 50 meters seaward. The distance between the baseline and each shoreline intersection point made it possible to obtain metric measurements and all shoreline change rates, expressed as meters of change along the transects per year. To determine the different behaviors of the sandy coastal stretches, these estimates were obtained for three different portions of the sandy coast, each in lateral continuity from West to East (Figure 7). For each time span, the minimum, maximum, weighted average, and error range were evaluated. The choice of the weighted average is related to the density data for a given time range (Figures 8 and 9).

Predictive Model of Submersion Surfaces Connected to Sea-Level Rise for 2100
Predicted coastal changes in response to sea-level rise may reflect a different behavior in function of the geomorphological elements of a sandy coast. These may determine three different scenarios: і) a barrier erosion with dispersed sediments on the shoreface just above the wave base, as described by Bruun [41]; іі) a rolling-over of the barrier, migrating onshore through mass relocation [43]; and ііі) an overstepping or in-place drowning in response to a fast sea-level rise with a complete submersion of the barrier. A mathematical submersion model was implemented according to the concept that the 3D geometric shoreline displacement and subsequent migration of the intertidal zone are caused by a sea-level rise and horizontal shoreline movements. The development of a submersion model along a sandy coast requires knowledge of the following parameters: і) sea-level trend; іі) VLM rates, including isostatic adjustments; and ііі) shoreline erosion/accretion. The model was implemented mathematically in a Matlab environment by considering the components conditioning both the vertical and horizontal coastal displacements for the sea level and shoreline changes: where: ∆z-relative sea-level rise (m) v RSLR -relative sea-level rate (mm/year) ∆t-prediction time span (year) ∆z tide -tide amplitude (m) ∆x-longitude shoreline displacement (m) ∆y-latitude shoreline displacement (m) v LRR -shoreline rate changes (m/year) β-normal shoreline angle (degrees).
The output provides all the points corresponding to the submersion surface predicted up to 2100, presented in a GIS-layer format with the xyz coordinates in world metric reference.

Shoreline Changes
The multi-temporal analysis of the shoreline changes highlighted a marked erosion between 1947 and 1967, followed by a significant shoreline accretion up to 1998 (already shown by [48,49]) and, finally, another erosion period up to 2018. Important changes in the entire sandy coast were revealed, together with river mouth changes ( Figure 10).
The first and most important event is connected to the shoreline accretion close to the Tara River. The 120 meter accretion which occurred in 1987 was significant, and can be attributed to the construction of the Taranto Harbor ( Figure 11).
The Patemisco River mouth showed a substantial modification in width during 1960, subsequent to the armored anthropogenic construction devised to deviate the longshore drift. In this area, a marked past change on the shoreline and dune ( Figure 12) can be noted. In particular, the nearshore bar movements caused by wave diffraction may be observed (Figure 12e).
The Lenne River revealed different full-mouth events, marked by the building of dykes. The latter hinder transversal sediment transports and amplify sand accumulation near the river mouth during storm events ( Figure 13). Compared to other rivers, the Lato River is less affected by anthropogenic contributions, and shows little accretion or variable mouth shape, connected to the sea-state regime and longshore currents ( Figure 14).     A variable dynamic is recognizable along the Castellaneta beach, close to the Lato River. This is directly connected to the sea-state influence together with the behavior of the Lato River mouth in the past. Indeed, the hot spot values of the shoreline migration for the transects, corresponding to the Lato River mouth (Figure 10), may be identified. During 1947 and up to 1961, a shoreline erosion trend of −0.38 ± 0.17 m/year occurred, with a maximum net movement of about 7.49 m. This was followed, in the subsequent years , by a shoreline accretion having the same maximum net movement.
From the Lenne River mouth, crossing the Pino di Lenne beach and Chiatona beach (Figures 2-6), a significant shoreline erosion rate of −0.42 ± 0.19 m/year occurred up to the 1970s. This was followed by a general accretion rate of 1.26 ± 0.66 m/year, occurring, in particular, in proximity to the Chiatona railway station, and to the numerous holiday resorts built along the coast. Subsequent sediment redistribution took place between the 1970s and the 2010s. Then, during the period 1990-2010, an erosion rate of −0.36 ± 0.18 m/year occurred.
During the mid-twentieth century , the shoreline change rate between the Marina di Ferrara and the Lido Azzurro highlighted the retreat of some portions of the dune systems, at a rate of −0.81 ± 0.36 m/year. Furthermore, with the construction of the Lido Azzurro holiday resort, an acceleration of this retreat occurred, showing that the current position of the shoreline corresponds to the position of the dune scarp in the 1970s. As of 1989, there has been an erosion rate of −0.23 ± 0.11 m/year, and, as of 1998, an erosion rate of −0.40 ± 0.13 m/year. On the contrary, the building of the breakwater structures in the Taranto harbor, together with the supply trend of the Tara River, changed the shoreline configuration dramatically during the 1970s-1990s.

Geomorphological Scenario for the Future
An analysis of the dune-beach system topography indicates geomorphological scenarios similar to those described by Davidson-Arnott [47], e.g., a sea-level rise corresponds to a shoreline retreat and beach/dune erosion as a result of the reduction in the sedimentary budget connected to a river transport decrease. A conceptual geomorphological model for the sea-level rise along the coast stretching between Castellaneta and Taranto was deduced by analyzing a sequence of aerial photos from 1947 to present, as well as by using direct surveys ( Figure 15). In the first half of the twentieth Century, a progradation of the mobile coastal system was observed [48,49,56,70] with a well-extended dune system up to 1 km landward stabilized by vegetation cover.
Between the 1950s and 2000s, two phases of coastal erosion were observed throughout the entire littoral area: sea-level rise and increase in the number of storm events, and a decrease in the amount of natural sedimentary nourishment coming from the land. A foredune erosion was observed, and the foredunes along the coastal stretches near the Lato River were completely destroyed. Furthermore, erosion of the secondary dunes was also observed. The building of holiday resorts in the area of the Chiatona and Lido Azzurro coastal stretches caused foredune leveling and increased vulnerability to effects of the sea-level rise.
Within the framework of an increase in the speed of the sea-level rise up to 2100 [27,31], together with an increase in storm intensity [32], the coasts around the Gulf of Taranto could be subjected to a shoreline regression, coupled with a partial foredune erosion and subsequent marine ingressions through river mouths, not to mention submersion of low-lands behind the dune ridge.

Simulating the Impact of a Rising Sea Level-A Difficult Methodology
The sandy coast of the Gulf of Taranto is characterized by a complex dynamics system, making the exact shoreline position very difficult to assess. First of all, it can be assumed that if the coast is considered to be in a steady-state, with no significant sediment movements, the shoreline migration should be subjected to only vertical displacements of the relative sea level. By considering the sea-level rise recorded at the Taranto-ISPRA station, different horizontal displacement rates, in function of the coastal slopes (Figure 16), can be derived. Second of all, geometric horizontal shoreline movements due to sea-level rise were subtracted from the observed shoreline movements in order to obtain effective shoreline change (Table 1) and to apply a correct value of horizontal displacement for the three coastal slopes considered in Figure 7.  The observed shoreline movements show very different results from the geometric horizontal migrations induced by a sea-level rise (Table 1 and Figure 10). According to this methodology, the sea-level contribution to the shoreline retreat is lower than 10%. Factors contributing to the shoreline retreat, other than the sea-level rise, such as the gradient of longshore sediment transport, negative sediment input due to dam and hydrological adjustments or construction of holiday resorts along the coast, were all considered. Albeit the investigated area is part of a bay, the beach dynamics have most probably been wave-dominated up to now, but increasing sea levels could modify this behavior.

Sea-Level Rise and Predicted Shoreline Positions for Two IPCC Scenarios
In order to estimate both the sea level for 2050 and 2100 AD and the subsequent sea-level rise scenario along the coast of Taranto, two regional sea-level (SL) projections, discussed in the Fifth Assessment Report of the IPCC-AR5 [23,24] and made available by the Integrated Climate data Center-ICDC of the University of Hamburg, rescaled for the Mediterranean Sea (http://icdc.cen.uni-hamburg.de/1/daten/ocean/ar5-slr.html [71]), were used. These data consist of mean values at the upper 95% and lower 5% confidence bounds of the SL, obtained by adding the contributions from geophysical sources driving long-term sea-level changes. The SL projections are based on two different Representative Concentration Pathways (RCP 2.6 and RCP 8.5), and estimated by including the thermosteric/dynamic contribution, obtained from the 21 Coupled Model Intercomparison Project Phase 5 (CMIP5). In addition, the atmosphere-ocean general circulation models (AOGCMs), the surface mass balance and dynamic ice-sheet contributions from Greenland and Antarctica, the glacier and land water storage contributions, the glacial isostatic adjustment, and the inverse barometer effect [23][24][25]31] were also employed. The predicted sea-level rise along the coast of Taranto was finally estimated as the sum of the above contributions and the tectonic uplift inferred from the elevation of the marine terraces of MIS 5.5 (placed at an elevation of 43 m), and the nearest Global Navigation Satellite System (GNSS) MMET stations, following the method previously applied in other areas of the Mediterranean [17,19,27,31,33,64]. The results are reported in Table 2 and Figure 17. It can be observed that both long-term and instrumental VLM rates are in agreement and provide similar mean uplift velocity of about 0.3 mm/year (MMET is uplifting at 0.29 ± 0.14 mm/year [29]; www.savemedcoasts.eu [64]), while the MIS 5.5 marine terrace provides values between 0.39 and 0.26 mm/year [40,52], supporting detailed relative sea-level (RSL) projections for 2050 and 2100. Over the last decades, the sea-level rise in the Gulf of Taranto has been causing a generalized landward shoreline displacement and the foredune erosion, while secondary dune belts have only been partially eroded owing to the stability provided by the vegetation cover (Figure 15). With a future rise in the sea-level, the sandy coast could become extremely vulnerable due to the submersion of the entire shoreface and the marine ingression through the inlets (Figure 18).
In view of a sea-level rise for 2100, the model developed in the present study describes a framework where inlets, constituted by rivers, determine preferential pathways to marine ingression. According to the Brunel and Sabatier model [9], shoreline retreats determine the destruction of the dunes during storms as some of the sand is projected behind the dunes (e.g., washover). Therefore, a shoreline retreat is accompanied by an increase in the elevation of the land behind the dunes which conditions the future position of the shoreline, and compensates, at least partially, for the sea-level rise [72]. To date, no methodologies which include this process have been proposed. Hence, it is necessary to interpret the present results attentively. In any case, river flow regimes will certainly be affected by the increase in the sea-level and, consequently, their ability to carry sediment onto beaches will be altered.

Conclusions
Over the past years, the Apulian coastal systems have been subjected to relevant erosion, in part, caused by meteo-marine processes and climate changes, and, in part, by anthropogenic activities. This coastal erosion, together with sea-level rise and storm flooding, is increasing the vulnerability of the coasts.
The authors' analysis, based on available data and predictive models, enabled the reconstruction of the past dynamics of the coastal area between Castellaneta and Taranto in function of the present-day trend, as well as the development of a submersion model for future sea-level rise scenarios.
Since the mid-twentieth century, the coast has been undergoing various changes in response to sea-level rise: • First half of the twentieth Century-progradation of the coastal system stretching from Castellaneta to Taranto, plantations of Pino d'Aleppo along the entire coast stabilized the dune system; • 1950s-1970s-shoreline erosion highlighted along the entire littoral stretch, the building of holiday resorts in proximity to railway stations (Castellaneta Marina, Chiatona, Lido Azzurro); • 1970s-construction of the Taranto Harbor outside the Mar Grande inlet with a substantial modification of the shoreline, and significant changes caused by the Tara, Patemisco, Lato and Lenne river mouths; • 1970s-1990s-shoreline accretion for most of the littoral area, mainly for the beaches in proximity to the Tara river mouth, except for the Lido Azzurro and Marina di Ferrara coastal stretches, where a significant shoreline erosion and dune retreat were observed; • 1999-2016-sea-level rise, recorded at the tide gauge in Taranto, at a rate of 1.04 ± 0.5 mm/year, and shoreline erosion recorded throughout the entire littoral area, with the exception of the area near the Tara mouth where accretion had continued; • 2016-2018-shoreline stabilization and slow accretion in some coastal stretches, such as the Pino di Lenne beach, Chiatona beach, and Tara river mouth and an exception was made for Castellaneta, Marina di Ferrara and Lido Azzurro where a slow erosion was observed.
Considering the RCP 8.5 of sea-level projections for 2050 and 2100, together with VLM, the model output shows how open beaches will be subjected to a significant shoreline retreat with a foredune erosion and subsequent marine ingressions ( Figure 18). This behavior will accelerate locally in coastal stretches already undergoing erosion, such as at Lido Azzurro and Metaponto. This will cause a complete foredune destruction with part of the sediments transported nearshore. On the other hand, adjacent to the coastal areas of Chiatona and Pino di Lenne, the model shows a significant foreshore submersion up to the dune scarp and marine ingressions through inlets of the river mouths.