Sedimentary archives of climate and sea-level changes during the Holocene in the Rhône prodelta (NW Mediterranean Sea)

Abstract

A 7.38 m long sediment core was collected from the eastern section of the Rhône prodelta (NW Mediterranean) at 67 m water depth. A multi-proxy study (including sedimentary facies, benthic foraminifera, ostracods, and clay mineralogy) provides a multi-decadal to century-scale record of climate and sea-level changes during the Holocene. The early Holocene is marked by alternative silt and clay layers interpreted as distal tempestites deposited in a context of rising sea level. This interval contains shallow infra-littoral benthic meiofauna (e.g., Pontocythere elongata,Elphidium spp., Quinqueloculina lata) and formed between ca. 20 and 50 m water depth. The middle Holocene (ca. 8.3 to 4.5 ka cal. BP) is characterized, at the core site, by a period of sediment starvation (accumulation rate of ca. 0.01 cm yr${}^{-\mathrm{1}})$ resulting from the maximum landward shift of the shoreline and the Rhône outlet(s). From a sequence stratigraphic point of view, this condensed section, about 35 cm thick, can be identified on seismic profiles as a maximum flooding surface that marks the transition between delta retrogradation and delta progradation. The transition between the early Holocene deposits and the middle Holocene condensed section is marked by a gradual change in all proxy records. Following the stabilization of sea level at a global scale, the late Holocene is marked by the establishment of prodeltaic conditions at the core site, as shown by the lithofacies and by the presence of benthic meiofauna typical of the modern Rhône prodelta (e.g., Valvulineria bradyana,Cassidulina carinata, Bulimina marginata). Several periods of increased fluvial discharge are also emphasized by the presence of species commonly found in brackish and shallow-water environments (e.g., Leptocythere spp.). Some of these periods correspond to the multi-decadal to centennial late Holocene humid periods recognized in Europe (i.e., the 2.8 ka event and the Little Ice Age). Two other periods of increased runoffs at ca. 1.3 and 1.1 ka cal. BP are recognized, which are likely to reflect periods of regional climate deterioration that are observed in the Rhône watershed. Conversely, the Migration Period Cooling (ca. 1.4 ka cal. BP) and the Medieval Climate Anomaly (ca. AD 950–1250) correspond locally to periods of increased dryness.

Introduction

Deltas comprise a subaerial delta plain, where river processes dominate; a coarser-grained delta front, where river and basinal processes interact; and a muddy submarine prodelta dominated by oceanic processes (Bhattacharya and Giosan, 2003; Galloway, 1975; Postma, 1995). Most of the world's deltas were initiated during the early Holocene, between ca. 9.5 and 6 ka cal. BP, owing to a deceleration of global sea-level rise (Stanley and Warne, 1994). They constitute key element of the continental margin system as they represent the first sink of sediments delivered by rivers (Trincardi et al., 2004).

Over the last decennia, numerous investigations have documented the land–sea evolution of these systems including the Amazon delta (Nittrouer et al., 1986), the Mekong delta (Ta et al., 2002; Xue et al., 2010), the Yellow River delta (J. P. Liu et al., 2004; Liu et al., 2007), and the delta of the Po River (Amorosi et al., 2008; Cattaneo et al., 2003). The Rhône delta is one of the most important of the Mediterranean Sea and has also been widely investigated combining seismic, sedimentological, and micropaleontological approaches (e.g., Boyer et al., 2005; Fanget et al., 2013a, b, 2014; Gensous et al., 1993; Labaune et al., 2005).

In this paper, we study the evolution of the Rhône prodelta in terms of sedimentary environments, during the last 10.5 kyr, as marked by changes in lithofacies and benthic meiofaunal assemblages (i.e., foraminifera and ostracods) and in relation to the Holocene sea-level rise and climate changes. This study shows that (1) major phases of sea-level rise and delta evolution can be clearly identified based on several independent proxy records, and that (2) changes in fluvial discharge inferred, in particular, from ostracod assemblages in the upper part of the core are linked to the last major periods of rapid climate changes of the Holocene (Mayewski et al., 2004; Wanner et al., 2014).

Regional geological and climatic setting

Geological evolution of the Rhône subaqueous delta

In the Gulf of Lion (NW Mediterranean), the Rhône delta (Fig. 1) occupies a deep valley incised during the Messinian and infilled with ca. 2 km of Plio-Quaternary sediments (Lofi et al., 2003), mainly delivered by the Rhône River (Aloïsi et al., 1977). For the last ca. 500 kyr, borehole data demonstrated that shelf deposits are primarily made up of forced-regressed sequences formed in response to 100 kyr glacio-eustatic cycles (Bassetti et al., 2008; Frigola et al., 2012; Sierro et al., 2009). These authors also demonstrated that higher frequency cycles, as well as suborbital climate changes, were nicely recorded within paleo-prodeltaic sedimentary archives. Following the Last Glacial Maximum (LGM, ca. 21 ka cal. BP; Mix et al., 2001), rapid sea-level rise led to the retrogradation of Rhône delta, and formation of a wedge of transgressive (backstepping) deposits thickening landward. The most prominent feature is an elongated paleo-deltaic complex, named the Early Rhône Deltaic Complex (Berné et al., 2007), which formed during the Younger Dryas and the Preboreal (Fig. 1). After ca. 7 ka, stabilization of sea level allowed the progradation of a series of regressive deltaic lobes (Fanget et al., 2014), corresponding to the overall eastward migration of the Rhône distributaries. In total, these transgressive and regressive deposits form the Rhône subaqueous delta that reaches, along the modern delta front, up to 50 m in thickness and pinches out at a present water depth of ca. 90 m (Gensous and Tesson, 1997).

Offshore and onshore morphology of the Rhône deltaic system. This map is based on the compilation of Berné et al. (2007) and Vella et al. (2008). The successive shifting of the Rhône distributaries under natural and/or anthropic influence during the middle and late Holocene led to the formation of several deltaic lobes. Steps in the present-day morphology correspond to the position of the paleo-delta fronts that can be linked to known paleo-distributaries of the Rhône: (1) early Saint Ferréol, (2) Peccais, (3) Bras de Fer, (4) Pégoulier, and (5) modern Roustan distributary. The map of relict morpho-sedimentary units in the Rhône delta plain (paleo-shorelines, beach ridges and onshore deltaic lobes) is based on L'Homer et al. (1981), Arnaud-Fassetta (1998), Vella (1999), and Provansal et al. (2003). Thick line and black dot correspond respectively to chirp seismic profile (Fig. 2) and sediment core RHKS-55 presented in this study.

[Figure omitted. See PDF]

The early Holocene deposits (called seismic unit U500), which rest on a wave ravinement surface (called D500), formed a transgressive parasequence (Labaune et al., 2005) made of tempestite deposits (Fanget et al., 2014). They are separated from middle and late Holocene deposits by a condensed section which corresponds on seismic profile to a maximum flooding surface (MFS, called D600). The age of the MFS varies along-strike between ca. 8 and 3 ka cal. BP, and reflects, at a given site, the duration of condensation and/or erosion (Fanget et al., 2014). After the stabilization of global sea level (ca. 7 ka cal. BP), the middle and late Holocene Rhône outlets shifted progressively eastward, under natural and/or anthropic influence. As a result, several deltaic lobes are formed (Fig. 1) (Arnaud-Fassetta, 1998; L'Homer et al., 1981; Provansal et al., 2003; Rey et al., 2005; Vella et al., 2005, 2008). Saint Ferréol, which is related to the “Rhône de Saint Ferréol” Channel, is the first and largest paleo-deltaic lobe. It started to prograde around 7 ka cal. BP (L'Homer et al., 1981). The Ulmet lobe, located to the east and linked to the “Rhône d'Ulmet” Channel formed simultaneously to the Saint Ferréol lobe. To the west of the Saint Ferréol lobe, the Peccais lobe, related to the “Rhône de Peccais” Channel, appeared to be posterior to the erosion of the St Ferréol lobe (Rey et al., 2005; Vella et al., 2005). During the Little Ice Age, the Bras de Fer lobe, linked to the “Rhône de Bras de Fer” Channel, formed between AD 1587 and 1711 (Arnaud-Fassetta, 1998). Until AD 1650, the “Rhône de Bras de Fer” Channel is considered as synchronous to the “Rhône du Grand Passon” Channel (Arnaud-Fassetta, 1998). The “Rhône de Bras de Fer” Channel shifted to the east up to the present-day position of the Grand Rhône River after several severe floods that occurred in AD 1709–1711. The progradation of these lobes is primarily influenced by changes in sediment fluxes and thus by climate (Arnaud-Fassetta, 2002; Bruneton et al., 2001; Provansal et al., 2003).

Holocene climate and its regional characteristics

High fluctuations in rainfall and low-amplitude temperature variations are observed during the Holocene (Davis et al., 2003; Mayewski et al., 2004; Seppä et al., 2009; Wanner et al., 2008, 2011). Examination of globally distributed paleoclimate records led to identification of six multi-decadal to century-scale cooling events interrupting periods of relatively stable and warmer climate (Mayewski et al., 2004; Wanner et al., 2008, 2011, 2014). These periods are known as rapid climate changes (RCCs; Mayewski et al., 2004) or cold relapses (CRs; Bassetti et al., 2016; Jalali et al., 2016; Wanner et al., 2014). The 8.2 ka cal. BP cold relapse (CR0, Table 1) took place at the beginning of the Holocene (Barber et al., 1999), a period of progressive warming that induced ice cap melting and freshwater outbursts to the oceans from North American glacial lakes.

Chronology of Holocene cold relapses (CRs) based on existing literature.

Event	Time slice	References
	(kyr)
CR0	8.2	Barber et al. (1999)
CR1	6.4–6.2	Wanner et al. (2011)
CR2	5.3–5.0	Magny and Haas (2004);
		Roberts et al. (2011)
CR3	4.2–3.9	Walker et al. (2012)
CR4	2.8–3.1	Chambers et al. (2007);
		Swindles et al. (2007)
CR5	1.45–1.65	Wanner et al. (2011)
CR6	0.65–0.15	Wanner et al. (2011)

During the warm middle Holocene, several significant episodes were identified at 6.4, 5.3, and 4.2 ka cal. BP (CR1, CR2, and CR3, respectively; Table 1) (Walker et al., 2012; Wanner et al., 2014). The 4.2 ka event, which may have played a role in the collapse of various civilizations (Magny et al., 2013), is characterized by increased drought in North America, Asia, and the southern Mediterranean region.

During the cooler late Holocene, the 2.8 ka cal. BP cold period (CR4, Table 1) might be at the origin of the decline of the Late Bronze Age civilization (Do Carmo and Sanguinetti, 1999; Weiss, 1982), and CR5 (Table 1) matches the Migration Period that occurred around 1.4 ka cal. BP (Wanner et al., 2014). The late Holocene cooling trend reached the highest development during the Little Ice Age (LIA, CR6 in Table 1) between the 13th and 19th centuries (Frezza and Carboni, 2009).

In the Rhône catchment area, these CRs (or at least CR0 and CR6) are marked by increased river runoff. CR0, or the so-called 8.2 ka event (Alley et al., 1997), is indeed marked by high lake levels generated by more intense rainfall (Magny and Begeot, 2004; Magny et al., 2003). The LIA (CR6), which is well documented in the Rhône watershed, is marked by a period of significant Alpine glacier advance (Goehring et al., 2011; Ivy-Ochs et al., 2009), high lake levels (Magny et al., 2010), and an increase in Rhône River floods (Pichard, 1995; Pichard and Roucaute, 2014).

Material and methods

The present paper is based on a multi-proxy study of a 7.38 m long piston core (RHS-KS55) collected in front of the paleo-“Rhône de Bras de Fer” and “Grand Passon” channels at 67 m water depth (latitude $=$ 43${}^{\circ }$14.35${}^{\prime }$ N; longitude $=$ 4${}^{\circ }$40.96${}^{\prime }$ E) during the RHOSOS oceanographic cruise (R/V Le Suroît, September 2008).

The detailed architecture of the Rhône subaqueous delta was determined from high-resolution seismic data (2000–5200 Hz Chirp system), together with core data constrained by ${}^{14}$C dates (Fanget et al., 2014). Radiometric dates were measured with accelerator mass spectrometer (AMS) ${}^{14}$C on benthic foraminifera or shells at Poznan Radiocarbon Laboratory (PRL, Poland) and at Laboratoire de Mesure ${}^{14}$C (LMC14) at Commissariat à l'Energie Atomique (CEA, France). Nevertheless, because of the low quantity of biogenic carbonates in the proximal part of the Rhône prodelta, we experienced difficulties in dating core RHS-KS55 and observed three age inversions from 336 to 122 cm. Based on seismic and lithofacies correlations at the regional scale (Fanget et al., 2014), we excluded these three ${}^{14}$C dates for core RHS-KS55. We also excluded the three older ${}^{14}$C dates (unit U500, base of the core) since measured ages are considered older than true ages (due to incorporation of reworked benthic foraminifera; for more information, see Sects. 4.1 and 5.1). As a result, eight robust ${}^{14}$C dates were used to create the age model for the last ca. 9 kyr. Ages are ${\mathit{\delta }}^{13}$C-normalized conventional ${}^{14}$C years and are corrected for an assumed air–sea reservoir effect of 400 years. Calendar ages were calculated using the program Clam (version 2.2, Blaauw, 2010) and the Marine13 calibration curve (Reimer et al., 2013). The age model was based on the linear interpolation between the dated levels using basic (non-Bayesian) age–depth modeling software (Blaauw, 2010).

Core RHS-KS55 was split lengthwise, imaged, and visually described for identifying sedimentological facies. Core RHS-KS55 was then subsampled for benthic microfaunal analyses (i.e., ostracods and foraminifera) and clay mineralogy. Three-centimeter-thick slides were collected using a 10 cm sampling step through the core, except from the bottom of the core to ca. 500 cm, where thick slides and sampling step were slightly modified. A total of 79 samples were sieved through a 63 $\mathrm{\mu }$m mesh screen and the residues were dried and dry-sieved using a 125 and 150 $\mathrm{\mu }$m mesh screens. Ostracods and benthic foraminifera were hand-sorted from the > 125 and > 150 $\mathrm{\mu }$m fractions, respectively, and stored in Plummer slides. To illustrate the diversity of benthic meiofauna, total abundance (values normalized for a 100 cm${}^{\mathrm{3}}$ sample volume), species richness ($S)$, Shannon index ($H)$, and evenness index ($E)$ (Hayek and Buzas, 1997; Shannon, 1948) were calculated, as described in Murray (2006), for each level. To highlight vertical patterns in benthic meiofaunal communities, hierarchical clustering was also performed on the 79 samples and on the major species (i.e., occurring with more than 5 % in at least one sample) by mean of the PAST^© software (version 2.09, 2011; Hammer et al., 2001). Cluster analyses were based on the arcsine values of the square of root “pi”, where “pi” is the relative abundance (%) of the species $i$ divided by 100. A tree diagram was constructed according to the Ward's method based on the squared Euclidean distances.

The fraction less than 63 $\mathrm{\mu }$m was used to perform clay mineral analyses. X-ray diffraction (XRD) on oriented mounts of carbonate-free clay-sized (< 2 $\mathrm{\mu }$m) particles was conducted to identify clay minerals with the PANalytical diffractometer, following the routine of the GEOPS Laboratory (Paris Sud University, France) (Z. Liu et al., 2004; Liu et al., 2008). Three XRD runs were carried out, following air drying, ethylene–glycol solvation over 24 h, and heating at 490 ${}^{\circ }$C over 2 h. The position of the (001) series of basal reflections on the three XRD diagrams was used to identify clay minerals. Semi-quantitative estimates of peak areas of the basal reflections for the main clay mineral assemblages of illite (10 Å), smectite (including mixed layers) (15–17 Å), and kaolinite $/$ chlorite (7 Å) were performed on the glycolated curve using the MacDiff software (Petschick, 2000). Relative proportions of kaolinite and chlorite were determined using the ratio 3.57 $/$ 3.54 Å of the peak areas.

Results

Seismic stratigraphic framework, age model, and sedimentation rates

Seismic discontinuities and seismic units described in the following section are based on Fanget et al. (2014). Several deglacial and Holocene seismic units bounded by well-marked discontinuities are identified at the core site (Fig. 2). Seismic data show that core RHS-KS55 goes through surface D500 and represents an expanded record of seismic units U500, U600, and U610. Seismic unit U620a is missing in the studied core (Fig. 2).

“Dip” (NS) chirp seismic profile across the Grand Passon and Bras de Fer subaqueous delta (position in Fig. 1). Surface D500 corresponds to the flooding surface (which coincides here with a wave ravinement surface) that separates the Younger Dryas (seismic unit U400) deposits from the Preboreal deposits (seismic unit U500). The downlap surface D600 is the maximum flooding surface formed during the turnaround between retrogradation and progradation. D610 and D620 are erosional surfaces that mimic flooding surfaces and separate the middle and late Holocene sedimentary wedges (seismic units U600, U610, and U620a) formed in response to the successive shifts of Rhône Channel. Twtt: two-way travel time.

[Figure omitted. See PDF]

The age of seismic unit U500, which corresponds to the early Holocene transgressive parasequence (Fanget et al., 2014), is poorly constrained in core RHS-KS55 (Figs. 3 and 4). Considering the age of the underlying deposits of seismic unit U400 in this area (i.e., paleo-deltaic complex of the Rhône (ERDC), ca. 10.5 ka cal. BP; Berné et al., 2007), we assume that the three ${}^{14}$C dates obtained at 732, 672, and 512 cm (i.e., 12.7, 13.3, and 12.2 ka cal. BP, respectively) are biased (Fig. 4) because of reworking occurring in shallow-water environments. This interpretation is supported by the nature of sediments that composed this unit. Indeed, at the regional scale, seismic unit U500 is made of tempestite deposits (Facies 1, Sect. 4.2.) which mainly contain infra-littoral (i.e., upper shoreface) benthic foraminifera (e.g., Elphidium crispum; see Table 2 and Fanget et al., 2014). These infra-littoral benthic foraminifera are likely reworked by high-energy hydrodynamic processes that regularly winnowed the seafloor. The three older ages obtained within seismic unit U500 are thus probably the result of reworking during the transgression of an underlying Younger Dryas/Preboreal delta front. It is likely that infra-littoral benthic foraminifera picked for these three ${}^{14}$C dates are reworked from older deposits (for more details, see Sect. 5.1).

Summary of ${}^{14}$C dates. Absolute dates were obtained with accelerator mass spectrometer (AMS) ${}^{14}$C dating of well-preserved shells and benthic foraminifera at Laboratoire de Mesure ${}^{14}$C (LMC14) at Commissariat à l'Energie Atomique (CEA, Saclay) and at Poznan Radiocarbon Laboratory (PRL). The ages reported herein are delta ${}^{13}$C-normalized conventional ${}^{14}$C years, corrected for an assumed air-sea reservoir effect of 400 years. Calendar ages were calculated using Clam software (Blaauw, 2010) and the Marine13 calibration curve (Reimer et al., 2013).

Depth	Material	Weight	Sample	${}^{14}$C conventional	Calibrated age	Mean calibrated
(cm)		(mg)	number	age (yr BP)	(yr cal. BP)	age (yr cal. BP)
90–93	Benthic foraminifera	9.5	SacA 27201	1335 $\pm$ 30	790–944	867
120–123	Benthic foraminifera	9.9	SacA 27202	1840 $\pm$ 30	1300–1477	1389${}^{*}$
150–153	Benthic foraminifera	11	SacA 23204	2080 $\pm$ 35	1551–1761	1656${}^{*}$
200–203	Benthic foraminifera	10.6	SacA 23205	1655 $\pm$ 30	1154–1282	1218
300–303	Benthic foraminifera$+$Turritella sp.	10.9	SacA 23206	1900 $\pm$ 30	1362–1527	1445
335–336	Benthic foraminifera $+$Turritella sp.	10	SacA 23207	1705 $\pm$ 30	1185–1318	1252${}^{*}$
	$+$ mixed bivalves
350–353	Benthic foraminifera	11.3	SacA 27203	2760 $\pm$ 35	2351–2619	2485
417–420	Turritella sp.	896	Poz-35061	4335 $\pm$ 35	4375–4574	4475
430–433	Benthic foraminifera	10.5	SacA 27204	6190 $\pm$ 40	6513–6735	6624
440–443	Nucula sp.	11.2	SacA 27205	7830 $\pm$ 40	8192–8374	8283
470–473	Benthic foraminifera	10.2	SacA 23208	8565 $\pm$ 35	9075–9333	9204
510–513	Elphidium crispum	10.1	SacA 23209	10790 $\pm$ 40	12058–12467	12 263${}^{*}$
670–673	Benthic foraminifera	13	SacA 23210	11855 $\pm$ 45	13215–13433	13 324${}^{*}$
730–733	Elphidium crispum	10.6	SacA 23211	11280 $\pm$ 40	12644–12870	12 757${}^{*}$

Mean calibrated ages with an asterisk were excluded from the age model since they are considered biased (see text for more details).

Age model for the upper part (i.e., from 472 to 0 cm; middle and late Holocene) of core RHS-KS55 based on linear interpolation using the Clam 2.2 software (Blaauw, 2010).

[Figure omitted. See PDF]

Sediment features and sedimentation rates of core RHS-KS55. Correlations with seismic units and discontinuities and with the Holocene chronology are shown on the left side of the figure. Black rectangles represent the most significant periods of climate deterioration (known as cold relapses (CRs); Wanner et al., 2014) during the Holocene (e.g., the 8.2 ka event (CR0) and the Little Ice Age (CR6)). AMS ${}^{14}$C dates are mentioned to the right; italic values correspond to age inversions. ${}^{14}$C dates are summarized in Table 2.

[Figure omitted. See PDF]

Only the uppermost part of this unit is confidently dated between ca. 9.2 and 8.3 ka cal. BP (Table 2). The upper boundary of seismic unit U500, called D600 and interpreted as a MFS (Fanget et al., 2014), corresponds to a condensed section formed between ca. 8.3 and 4.5 ka cal. BP in the studied core (Fig. 3). It is characterized by a very low sedimentation rate of ca. 0.01 cm yr${}^{-\mathrm{1}}$ (Fig. 4). Seismic unit U600 progrades on D600, and is related to the marine component of the St Ferréol and Ulmet lobes (Fanget et al., 2014). ${}^{14}$C dates indicate that seismic unit U600 was deposited between ca. 4.5 and 0.9 ka cal. BP (Fig. 3). Sedimentation rates through this interval oscillated between 0.03 and 0.4 cm yr${}^{-\mathrm{1}}$ (Fig. 4). Highest sedimentation rates are recorded along the uppermost part of this unit (i.e., between 300 and 110 cm, Fig. 4). Finally, seismic unit U610, corresponding to the activity of the Grand Passon and Bras de Fer channels, was formed between ca. 900 and 280 yr cal. BP (Fig. 3) (Fanget et al., 2014). As seismic unit U620a is missing within core RHS-KS55, we estimate that the top of the core has an age of ca. 280 yr cal. BP. Sedimentation rate through seismic unit U610 is estimated at ca. 0.11 cm yr${}^{-\mathrm{1}}$ (Fig. 4).

Sedimentary features

Based on lithological description (including lithofacies, sedimentary structures, bioturbation, and color) of core RHS-KS55 (Fanget et al., 2014), three main sedimentary facies (i.e., Facies 1, 2, and 3 (including 3a and b)) are identified and summarized as follows:

• Facies 1: from 738 (core bottom) to 460 cm, core RHS-KS55 consists of numerous silt or very fine sand laminae (in the sense of Campbell, 1967) interlaminated with grayish and beige silty clay with spacing of millimeters to several centimeters (Fig. 4). Within these very thin laminae (millimeters to a few centimeters thick), which are characterized by erosional basal contacts, no sedimentary structures can be identified.

• Facies 2: from 460 to 430 cm, a peculiar interval consisting of heterolithic content in a grayish silty clay matrix is observed (Fig. 4). Turritella sp., as well as several bivalves (e.g., Acanthocardia echinata, Arca tetragona, Nucula sp.) and bryozoans, is identified.

• Facies 3a: from 430 to 320 cm, sediment consists of beige silty clay without visible sedimentary structures. Diffuse veneers of yellowish lighter levels and spots of darker deposits (richer in hydrotroilite), clearly obliterated by intense bioturbation, are observed (Fig. 4). Scattered bryozoans debris, Turritella sp., and bivalves are encountered in this interval.

• Facies 3b: from 320 to 0 cm, sediment consists of grayish and beige silty clay and contains abundant hydrotroilites and bioturbation.

Clay mineralogy

Clay mineral assemblages are dominantly composed of illite, with values ranging from ca. 55 to 85 % (Fig. 5). Illite content exhibits high and quite constant values (ca. 70 %) from the bottom of the core to 465 cm. From 455 to 360 cm, illite content decreases to ca. 60 %. A strong increase is observed at 340 cm, with values reaching ca. 80 %. Illite content remains high between 340 and 110 cm. At 110 cm, illite content drops to ca. 55 % and is generally lower from 100 to 30 cm. Along the uppermost 30 cm of the core, a progressive increase in illite content is identified. Smectite content is inversely correlated to illite content throughout core RHS-KS55 (Fig. 5). It has very low values (ca. 2.5 %) from 738 to 465 cm and exhibits higher (from ca. 5 to 25 %) and erratic values along the uppermost 460 cm of the core. Chlorite content reaches ca. 20 % in the lower part of the core (from the bottom to 350 cm) and drops near 0 % between 350 and 120 cm (Fig. 5). Chlorite content increases to ca. 20 % along the uppermost 120 cm of the core. Kaolinite content is very low throughout core RHS-KS55 with values close to 0 % from the bottom of the core to 360 cm, and ranging from ca. 7 to 10 % between 360 cm and the top of the core (Fig. 5).

Distribution of the main clay minerals in core RHS-KS55. Correlations with seismic units and discontinuities are shown on the left side of the figure. The corresponding positions of seismic units and discontinuities (U500, U600, U610, D600, and D610) are also reported on the plots (gray bands and dashed lines). Black rectangles on the left side of the figure represent the most significant periods of climate deterioration (i.e., cold relapses (CRs)) during the Holocene.

[Figure omitted. See PDF]

Ecological indices (total abundance per 100 cm${}^{\mathrm{3}}$, number of species ($S$), Shannon ($H$), and evenness ($E$) index) describing benthic foraminiferal and ostracod populations. Correlations with seismic units and discontinuities are shown on the left side of the figure. The corresponding positions of seismic units and discontinuities (U500, U600, U610, D600, and D610) are also reported on the plots (gray bands and dashed lines). Black rectangles on the left side of the figure represent the most significant periods of climate deterioration (i.e., cold relapses (CRs)) during the Holocene.

[Figure omitted. See PDF]

Ostracod fauna

Ostracod density and diversity indices

The number of ostracods per sample varies greatly from 68 to 12 821 ind./100 cm${}^{\mathrm{3}}$ in core RHS-KS55 (Fig. 6). From the base of the core to 482 cm, density ranges from 68 to 1266 ind./100 cm${}^{\mathrm{3}}$. The lowest values are observed from 642 to 632 cm, whereas four peaks of 1141, 1243, 1006, and 1266 ind./100 cm${}^{\mathrm{3}}$ are identified at 702, 622, 571, and 541 cm, respectively. From 471 to 352 cm, the number of counted ostracods increases strongly with values reaching up to 12 821 ind./100 cm${}^{\mathrm{3}}$ and 10 539 ind./100 cm${}^{\mathrm{3}}$ at 432 and 412 cm, respectively. The number of ostracods per sample drops to 222 ind./100 cm${}^{\mathrm{3}}$ at 332 cm and ranges from 305 to 1182 ind./100 cm${}^{\mathrm{3}}$ between 322 and 162 cm. Density decreases between 152 and 82 cm and then increases progressively along the uppermost 82 cm of the core.

Species richness ($S)$ oscillates between 7 and 31 species per sample through the core and approximately follows the same trend as density (Fig. 6). Lowest values of $S$ are recorded when ostracods abundances are minimal, i.e., at 642 and 632 cm and between 152 and 82 cm. $S$ is maximal within the interval between 482 and 442 cm. The Shannon index ($H)$ varies between 1.2 and 3.0 through the core (Fig. 6). $H$ increases progressively from the base of the core to 492 cm. $H$ reaches maximal values (ca. 3.0) between 482 and 442 cm and decreases progressively from 433 to 82 cm. Along the uppermost 82 cm of the core, $H$ gradually increases (from ca. 1.4 to 2.2). The evenness index ($E)$, which is between 0.2 and 0.7, follows approximately the same trend as $H$ (Fig. 4). The lowest values of $E$ are observed from 738 to 492 cm, whereas the maximal values are observed from 482 to 422 cm. $E$ is relatively constant along the uppermost 412 cm of the core, with values oscillating around 0.5.

R-mode cluster analyses of the 16 major (more than 5 % in at least one sample) benthic foraminiferal species, and of the 21 major ostracod species according to Ward's method, based on standardized percentages pi of these species.

[Figure omitted. See PDF]

Cluster analysis

R-mode cluster analysis allows us to identify six ostracod clusters plus one single species when a cut-off level of 1.4 is applied (Figs. 7 and 8).

Cluster A is made of Semicytherura incongruens, Pontocythere elongata, and Semicytherura sp. It has a maximal contribution from 738 to 512 cm, with erratic values ranging from ca. 5 to 25%. It decreases strongly from 512 to 442 cm and disappears completely along the uppermost 432 cm of the core.

Cluster B is composed of Propontocypris pirifera, Cytherissa sp., Eucythere sp., Aurila sp., and Cytheridea neapolitana. It shows a low contribution through the core, with values generally < 10 %. It increases only between 482 and 442 cm, where a maximal value of 25 % is reached.

Cumulative percentages of taxa composing the defined ostracod clusters along core RHSKS55. Correlations with seismic units and discontinuities are shown on the left side of the figure. The corresponding positions of seismic units and discontinuities (U500, U600, U610, D600, and D610) are also reported on the plots (gray bands and dashed lines). Holocene cold relapses (CRs) are also indicated on the $y$ axis by the black rectangles (left side of the figure). Periods of increased fluvial discharge are highlighted by peaks in Leptocythere (hatched rectangles). These periods of intensified runoffs correlate to the 2.8 ka event (i.e., CR4) and the LIA (i.e., CR6). Two periods of regional climate deterioration are also observed at 1.3 and 1.1 ka cal. BP. Conversely, the Migration Period Cooling (i.e., CR5) and the Medieval Climate Anomaly (MCA) correspond to periods of increased dryness.

[Figure omitted. See PDF]

Cluster C is constituted by Cytherella sp., Cytheropteron alatum, Cytheropteron monoceros, and Carinocythereis carinata. It exhibits a low contribution (< 10 %) from 738 to 492 cm. It increases strongly between 492 and 372 cm, to reach up to ca. 40 % of the ostracod fauna at 432 cm. Cluster C decreases progressively, and has a minimal contribution along the uppermost 362 cm of the core.

The single species corresponds to Leptocythere spp. This species dominates ostracod fauna from 738 to 492 cm, with values ranging from ca. 41 to 73 %. Along the uppermost 492 cm of the core, Leptocythere spp. oscillates between low (< 5 %) and high (ca. 30–40 %) values. The lowest contributions of this species are recorded from 472 to 432, 332 to 302, 232 to 182, and 122 to 82 cm.

Cluster D is made of Cytheropteron rotundatum and Krithe spp. (juvenile Krithe and K. pernoides). From 738 to 362 cm, it has a low contribution (< 5 %). It increases strongly along the uppermost 352 cm of the core, with values ranging from ca. 30 to 60 %, and reaches a peak of ca. 77 % at 92 and 82 cm.

Cluster E is composed of Argilloecia spp., Loxoconcha laevis, and Paradoxostoma sp. It exhibits a moderate contribution through the core, with values oscillating generally between ca. 5 and 15 %. However, four peaks of ca. 37, 45, 35, and 27 % are recorded at 492, 342–332, 112, and 72 cm, respectively.

Cluster F is constituted by Sagmatocythere sp., Bosquetina dentata, and Pterigocythereis jonesii. It shows a minimal contribution (< 5 %) from 738 to 492 cm. It increases between 482 and 162 cm, with values ranging from ca. 20 to 35 %, except between 352 and 332 cm, where a decrease (< 15 %) is observed. From 162 to 72 cm, Cluster F exhibits erratic values oscillating between ca. 2 and 20 %. Along the uppermost 72 cm of the core, it increases slightly.

Benthic foraminiferal fauna

Benthic foraminiferal density and diversity indices

The number of benthic foraminifera per sample varies greatly from 404 to 74 642 ind./100 cm${}^{\mathrm{3}}$ through the core (Fig. 6). From 738 to 482 cm, density shows erratic values oscillating between 404 and 7130 ind./100 cm${}^{\mathrm{3}}$ and increases progressively. The number of counted specimens increases strongly between 472 and 382 cm, with values ranging from 12 386 to 74 642 ind./100 cm${}^{\mathrm{3}}$. From 382 to 352 cm, density decreases rapidly and drops to 2995 ind./100 cm${}^{\mathrm{3}}$. Along the uppermost 342 cm of the core, benthic foraminiferal densities remain quite constant and below 3000 ind./100 cm${}^{\mathrm{3}}$.

Species richness ($S)$ oscillates between 13 and 49 species per sample through the core (Fig. 6). The lowest values of $S$ (from 13 to 26) are recorded from 738 to 632 cm. $S$ increases from 622 to 492 cm, except between 542 and 532 cm, where a decreased is observed. The highest values of $S$ are encountered from 482 to 212 cm, with values oscillating between 34 and 49 species per sample. $S$ slightly decreases along the uppermost 202 cm of the core and remains quite constant. The Shannon index ($H)$ varies between 1.2 and 3.2 through the core RHS-KS55 (Fig. 6). From 738 to 642 cm, $H$ decreases progressively from 2.3 to 1.2. Erratic values, oscillating between 1.6 and 2.7, are observed from 632 to 492 cm. Even if $H$ slightly decreased between 452 and 402 cm, the highest values are recorded between 482 and 352 cm. From 342 to 132 cm, $H$ decreases progressively, whereas it increases slightly along the uppermost 122 cm of the core. The evenness index ($E)$ exhibits relatively low values through the core, ranging from 0.2 to 0.5, and follows exactly the same trend as $H$ (Fig. 6).

Cluster analysis

R-mode cluster analysis allows us to distinguish four benthic foraminiferal clusters when a cut-off level of 2.4 is applied (Figs. 7 and 9).

Cluster 1 is composed of Elphidium spp. (including E. advenum, E. crispum, E. decipiens, E. granosum, E. incertum, E. macellum, and E. margaritaceum), Nonionella turgida, and Quinqueloculina lata. From 738 to 492 cm, Cluster 1 dominates strongly with abundances oscillating between ca. 58 and 91 %. The contribution of Cluster 1 drops to ca. 15 % at 472 cm. It has a minimal contribution along the uppermost 472 cm of the core, and especially between 342 and 62 cm, where it exhibits values < 10 %.

Cluster 2 is constituted by Cassidulina carinata, Bulimina marginata, and Valvulineria bradyana. It is characterized by a low contribution (less than 10 %) from 738 to 492 cm. It increases progressively from ca. 15 to 72 % between 482 and 202 cm and remains quite constant along the uppermost 192 cm of the core, with values oscillating between ca. 50 and 72 %.

Cluster 3 is made of Haynesina depressula, Ammonia beccarii, and Eggerella scabra. It has a low contribution (< 15 %) through the core RHS-KS55. From 738 to 492 cm, it exhibits relatively erratic values, and three peaks of ca.  12, 5, and 6 % are observed at 699, 632, and 532 cm, respectively. Cluster 3 has a minimal contribution (< 1 %) between 482 and 342 cm and increases slightly from 332 to 162 cm. From 162 to 82 cm, it increases progressively to reach a value of ca. 9 %. Cluster 3 decreases again between 82 and 32 cm and increases slightly along the uppermost 32 cm of the core.

Cluster 4 is constituted by Pseudoeponides falsobeccarii, Textularia agglutinans, Hyalinea balthica, Melonis barleeanus, Bulimina aculeata, Sigmoilopsis schlumbergeri, and Cibicides lobatulus. It shows a low contribution from 738 to 492 cm, with erratic values ranging from ca. 5 to 22 %. It increases strongly between 482 and 402 cm, where values of ca. 55 % are reached. From 402 to 292 cm, Cluster 4 decreases progressively. It remains quite constant along the uppermost 292 cm of the core, with values oscillating between ca. 13 and 25 %.

Discussion

Benthic meiofauna reworking processes in subaqueous deltaic environment

In subaqueous deltaic environments, reworking processes appear to be common within transgressive deposits (Cattaneo and Steel, 2003). In the Rhône subaqueous delta, transgressive deposits consist of tempestite deposits (seismic unit U500), which are the result of regular occurrence of high-energy hydrodynamic processes (including combined storm and flood events; Fanget et al., 2014). These processes regularly winnowed the seafloor and generate erosion, reworking, and transport of sediments. Thus, it is likely that benthic calcareous meiofauna are reworked from older deposits into modern deposits having the same faunal assemblages (Cearreta and Murray, 2000). These reworked benthic meiofauna cannot be considered as in situ, but it appears impossible to distinguish them from the unreworked modern tests and carapaces. This will directly affect AMS dating with measured ages older than true ages, as observed within the transgressive seismic unit U500. Such phenomena have been observed in Denmark (Heier-Nielsen et al., 1995) and Spain (Cearreta and Murray, 2000), and these highlight the difficulty to obtain reliable AMS dates from high-energy transgressive deposits.

Within the most recent prograding units of the highstand systems tract (4.5 to 0.3 ka cal. BP in the present study), we also observe the regular occurrence of reworking and transport of benthic meiofauna. Reworking processes are regularly encountered in shallow-water environments (e.g., Frenzel and Boomer, 2005; Loureiro et al., 2009; Fanget et al., 2013a). Conversely to the transgressive systems tract, reworked benthic meiofauna are easier to identify since they originate from shallow-water environments and deposit into deeper settings. Reworking processes in AMS dating are thus considered as less important and problematic in the highstand systems tract. It is likely than these allochthonous benthic meiofauna are transported and redeposited further offshore within the river plume during periods of increased river discharge (Fanget et al., 2013a). Thus, it can be relevant to use allochthonous meiofauna as biomarkers for a better understanding of transport and reworking processes (Cronin, 1983; van Harten, 1986; Zhou and Zhao, 1999; Fanget et al., 2013a; Angue Minto'o et al., 2015) and to study paleo-hydrology. The distribution pattern of reworked benthic meiofauna through highstand deposits is likely to reflect hydrological fluctuations in the past (see Sect. 5.3).

Cumulative percentages of taxa composing the defined benthic foraminiferal clusters along core RHS-KS55. Correlations with seismic units and discontinuities are shown on the left side of the figure. The corresponding positions of seismic units and discontinuities (U500, U600, U610, D600, and D610) are also reported on the plots (gray bands and dashed lines). Holocene cold relapses (CRs) are also indicated on the $y$ axis by the black rectangles (left side of the figure).

[Figure omitted. See PDF]

Record of Holocene sea-level rise and Rhône delta evolution

Seismic stratigraphy, sedimentological (including clay minerals) and benthic meiofauna data described in the previous section (Sect. 4) allow the subdivision of the studied core into three main intervals. These intervals fairly match the tripartite subdivision of the Holocene (Walker et al., 2012; Wanner et al., 2014) and might be linked to the Holocene sea-level history as well as the Rhône deltaic system evolution.

Interval 1 (ca. 10.5–8.3 ka cal. BP)

This interval encompasses most of the early Holocene. The age of the bottom of the core up to ca. 460 cm cannot be dated precisely because it corresponds to a transgressive parasequence (seismic unit U500, Fig. 2) that formed in a context of shallow-marine environment. Based on the age of the underlying deposits (i.e., the ERDC, seismic unit U400) and on ${}^{14}$C dates, this interval was deposited between ca. 10.5 and 8.3 ka cal. BP (i.e., the early Holocene; Walker et al., 2012). Tempestite (storm-induced) deposits, which are commonly found in lower to middle shoreface environments during periods of storm decelerating flows (Myrow, 1992; Myrow and Southard, 1996; Pérez-López and Pérez-Valera, 2012), characterize this interval (Fig. 4). The intercalation of fine clay and silt layers suggests that these deposits are distal tempestites (i.e., turbidite-like deposited below the storm wave base; Myrow, 1992; Pérez-López and Pérez-Valera, 2012). This facies is interpreted to correspond to an hydrodynamic regime resulting from the combination of E-SE storm waves and flood events (i.e., “wet storms” of Guillén et al., 2006), which regularly winnow the seafloor (Fanget et al., 2014). Tempestite deposits mainly contain foraminifera belonging to Cluster 1 (Elphidium spp., N. turgida, Q. lata) and ostracods belonging to Cluster A (S. incongruens, P. elongata, Semicytherura sp.), as well as to the genus Leptocythere (Figs. 8 and 9). Benthic foraminiferal species, like N. turgida, are typical of shallow prodeltaic environment enriched in organic matter of continental origin (e.g., Barmawidjaja et al., 1992; De Rijk et al., 2000; Diz and Francés, 2008; Van der Zwaan and Jorissen, 1991). Elphidium spp. and Q. lata are commonly reported in sandy silty substrates subject to strong hydrodynamic processes (e.g., Donnici and Serandrei Barbero, 2002; Jorissen, 1988; Rossi and Vaiani, 2008; Sgarrella and Moncharmont Zei, 1993). Similar observations are described in the modern Rhône subaqueous delta (Goineau et al., 2011, 2015; Mojtahid et al., 2009). Similarly, ostracods content of Cluster A are represented by littoral to sublittoral/phytal marine forms (e.g., Bonaduce et al., 1975; Cabral et al., 2006; Carbonel, 1980; Peypouquet and Nachite, 1984; Zaïbi et al., 2012). The genus Leptocythere is commonly found in brackish and shallow-water environments, and many Leptocythere are known to be euryhaline species (e.g., Anadon et al., 2002; Boomer and Eisenhauer, 2002; Carbonel, 1973, 1980; Frenzel and Boomer, 2005; Gliozzi et al., 2005; Van Morkhoven, 1963).

According to paleoenvironmental reconstruction based on benthic meiofauna from core RHS-KS55, the early Holocene is characterized, in the Rhône subaqueous delta, by high-energy hydrodynamic processes and significant organic matter input of continental origin typical of shallow infra-littoral setting. This interpretation is in agreement with the occurrence of tempestite deposits, as well as the global estimates of sea-level rise during the early Holocene (e.g., Bard et al., 1996; Fairbanks, 1989; Smith et al., 2011). Based on the sea-level curve of Stanford et al. (2011), the base of the core (estimated at ca. 10.5 ka cal. BP) corresponds to a sea level of ca. 50 m below its present-day position. Due to the location of core RHS-KS55 at a water depth of 67 m and its length of 7.38 m, a paleo-water depth of ca. 24 m can be estimated at the base of the core (subsidence and compaction being considered as negligible). At the top of the tempestite facies (i.e., at ca. 460 cm), which is dated at ca. 9.2 ka cal. BP, a paleo-water depth of ca. 52 m is estimated. Thus, the resulting rate of sea-level rise within this interval (738–460 cm) is ca. 20 mm yr${}^{-\mathrm{1}}$. This value matches the one found by Stanford et al. (2011) for the early Holocene. Thus, we consider that tempestite deposits, preserved within this transgressive interval (seismic unit U500), are formed at water depth ranging from ca. 20 to 50 m. In the Rhône subaqueous delta, we consider the tempestite facies as a relatively good paleo-bathymetric marker and we have been able to correlate it over a large prodelta area (see cores RHS-KS40, RHS-KS22, and RHS-KS39 in Fanget et al., 2014).

Interval 2 (ca. 8.3–4.5 ka cal. BP)

The interval between 460 and 430 cm corresponds to a period ranging from ca. 8.3 to 4.5 ka cal. BP (i.e., the middle Holocene, Fig. 3), when the Rhône outlet(s) was situated 10 to 30 km landward from the modern shoreline. Considering the resolution of our seismic data (in the order of ca. 0.5 m), it corresponds to the position of the MFS (surface D600, Fig. 2) that marks the transition between retrogradation and progradation. Very low sediment accumulation (ca. 0.01 cm yr${}^{-\mathrm{1}}$, Fig. 4), abundant shell concentration, and very rich microfossil content (up to $\sim$ 13 000 ostracods/100 cm${}^{\mathrm{3}}$ and $\sim$ 75 000 foraminifera/100 cm${}^{\mathrm{3}}$, Fig. 6) indicate sediment starvation within this condensed section which separates transgressive (below) from regressive (above) deposits. This condensed section consists of a silty clay matrix incorporating coarse-grained sediments with reworked shoreface material and shell hash. Clay mineralogy assemblages indicate a clear change in clay mineralogy proportions with a significant decrease in illite content (from 80 to 60 %) and a sharp increase in smectite content (from 0 to 20 %; Fig. 5). Benthic foraminifera belonging to Cluster 4 (P. falsobeccarii, T. agglutinans, H. balthica, M. barleeanus, B. aculeata, S. schlumbergeri, C. lobatulus), and ostracods belonging to Cluster B (P. pirifera, Cytherissa sp., Eucythere sp., Aurila sp., and C. neapolitana), Cluster C (Cytherella sp., C. alatum, C. monoceros, and C. carinata), and Cluster F (Sagmatocythere sp., B. dentata, and P. jonesii) are dominant within this interval (Figs. 8 and 9). Except for C. lobatulus, which is preferentially found in a high-energy shallow-water setting (e.g., Bartels-Jónsdóttir et al., 2006; Javaux and Scott, 2003; Milker et al., 2011; Murray, 2006), foraminifera assemblage is mainly composed of species thriving under a stable environment characterized by marine-derived organic matter input and well-oxygenated sediments (e.g., De Rijk et al., 2000; Debenay and Redois, 1997; Fontanier et al., 2008; Goineau, 2011; Goineau et al., 2011, 2015; Mendes et al., 2004; Mojtahid et al., 2009). Clusters B and F are mainly composed of shallow infra-littoral ostracods (e.g., Bonaduce et al., 1975; Frenzel and Boomer, 2005; Guernet et al., 2003; Ruiz et al., 1997; Zaïbi et al., 2012), whereas Cluster C is primarily made of circa-littoral and epi-bathyal species (e.g., Bonaduce et al., 1975; El Hmaidi et al., 2010; Yamaguchi and Norris, 2012).

The middle Holocene condensed section is very well identified thanks to benthic microfossils indicating mixed assemblages belonging to diverse environments, from infra-littoral to epi-bathyal settings. Shallow-water species highlight incorporation of the previous shoreface and delta mouth sediments that were left in situ during the transgressive submersion. Circa-littoral and epi-bathyal species indicate an abrupt increase in water depth (peak of transgression) and mark the time of maximum landward shift of the shoreline.

Interval 3 (ca. 4.5–0.3 ka cal. BP)

The most recent interval (from 430 cm to the top of core) corresponds to seismic units U600 and U610 (Fig. 2), which formed during the late Holocene a series of regressive deltaic lobes that make up the highstand systems tract in the sequence stratigraphic terminology. They consist of fine-grained prodeltaic deposits and are related to the activity of the St Ferréol and Ulmet distributaries (seismic unit U600), as well as to the synchronous, and thus successive, activity of the Grand Passon and Bras de Fer channels (seismic unit U610) (Fanget et al., 2014). At the core site, clay minerals are dominated by illite, as elsewhere in the Rhône prodelta (Chamley, 1971). Indeed, the Rhône River, receiving its detrital material primarily from the Alps, is particularly rich in illite, associated with some chlorite (Chamley, 1971) that tends to be trapped in sandy sediments during deposition (Chamley, 1971; Giresse et al., 2004). Both minerals represent the relative contribution of physical weathering to sedimentation, since they are resistant to degradation and transport (Chamley, 1971). Relative contents of illite are changing simultaneously with changes in sedimentary facies and activity of different distributaries (Fig. 5). Smectite content is low as a whole but higher than within underlying intervals when sea level was lower. The onset of seaward progradation of the Rhône deltaic lobes corresponds, by definition, to the age of deposits situated immediately above the condensed section. This age is ca. 4.5 ka cal. BP according to the age model. It corresponds to a marked increase in smectite content. It also corresponds to a marked increase in benthic foraminifera belonging to Cluster 2 (C. carinata, B. marginata, and V. bradyana) and ostracods belonging to cluster D (C. rotundatum and Krithe spp.) (Figs. 8 and 9). Foraminifera assemblage is constituted by typical species living in the distal part of the Rhône prodelta, with fine-grained sediments enriched in both terrestrial and marine organic matter (Goineau et al., 2011, 2015; Kruit, 1955; Mojtahid et al., 2009). They are also reported as opportunistic species which respond quickly to fresh phytodetritus input by higher reproduction (De Rijk et al., 2000; Fontanier et al., 2003; Goineau et al., 2011; Jorissen, 1987). Ostracods content of Cluster D is known as common assemblage of circa-littoral to epi-bathyal environments (Bonaduce et al., 1975; Coles et al., 1994; Cronin et al., 1999; Didié et al., 2002; Yamaguchi and Norris, 2012). In the Rhône subaqueous delta, we hypothesize that these species can be tolerant to moderate river influence (Fanget et al., 2013b). At the core site, strong decreases in Cluster 1 and Cluster B (shallow infra-littoral species), as well as increases in Cluster 2 and Cluster D, reveal the establishment of prodeltaic conditions since 4.5 ka cal. BP (Figs. 8 and 9). More precisely, they correspond to the progradation of the St Ferréol and Ulmet lobes. A similar pattern is identified on boreholes in the Rhône delta plain, where the onset of prodeltaic sedimentation is marked by the dominance of V. bradyana around 4 ka cal. BP (Amorosi et al., 2013).

Within Interval 3, we note also the presence of benthic foraminifera belonging to Cluster 3 ($H.$ depressula, A. beccarii, E. scabra), and ostracods belonging to the genus Leptocythere and to Cluster E and F (Figs. 8 and 9). The vertical pattern of these ostracods in this interval will be discussed in further details in the next section (5.3). Foraminifera constituting Cluster 3 are typical of very shallow-water environments, and E. scabra is notably known to be adapted to thrive in organic-matter-enriched and hypoxic sediments (Diz and Francés, 2008; Donnici and Serandrei Barbero, 2002; Jorissen, 1987; Mendes et al., 2004). This assemblage increases in the uppermost 300 cm of the core, in concomitance with increased hydrotroilite content. Authigenic minerals generated by sulfate reduction (hydrotroilite) can be related both to high sedimentation rate (as observed in the core), leading to reducing conditions, and to high organic matter supply. These observations suggest increased river influence that can be linked to the progressive progradation of Rhône delta and to the beginning of the Bras de Fer and Grand Passon channels' activity, located in front of the studied core.

Record of Holocene cold events (CRs)

Ostracods belonging to Clusters E and F, and especially to the genus Leptocythere, show well-marked peaks within highstand prodeltaic deposits (Fig. 8). As previously described, the genus Leptocythere is widely distributed in brackish and shallow marine water environments (Anadon et al., 2002; Boomer and Eisenhauer, 2002; Carbonel, 1973, 1980; Frenzel and Boomer, 2005; Gliozzi et al., 2005; Van Morkhoven, 1963). In the Po delta, the occurrence of Leptocythere sp. is notably related to local increase in fluvial influence (Rossi, 2009), and in the Rhône delta, few valves of Leptocythere are encountered in restricted environmental areas characterized by estuarine conditions (Amorosi et al., 2013). Thus, the distribution pattern of Leptocythere through the highstand deposits would reflect hydrological fluctuations. During the Holocene, high fluctuations in precipitations are recorded, notably during the CRs (Mayewski et al., 2004; Wanner et al., 2014). In Europe, these CRs (or at least CR0, i.e., the 8.2 ka event, and CR6, i.e., the LIA) are characterized by intensified rainfalls (Arnaud et al., 2012; Magny et al., 2010; Magny and Begeot, 2004).

Two intervals of increased occurrence of Leptocythere (and therefore increased rainfall) are identified between ca. 400 and 350, and ca. 70 and 0 cm (Fig. 8). They correspond to ages between ca. 4.0 and 2.2 and between 0.6 and 0.2 ka cal. BP, respectively. These intervals are close to CR4 and CR6 (i.e., the LIA) that are dated between ca. 3.1 and 2.8 and ca. 0.65 and 0.15 ka cal. BP, respectively (Wanner et al., 2014). The hypothesis of increased rainfall and river runoff, in the Rhône watershed, during the cooler late Holocene is supported, at least for the LIA, by the observed advance of the Rhône Glacier (Goehring et al., 2011), the high level in Lake Bourget (France) (Arnaud et al., 2012), higher soil erosion in the French Pre-Alps (Simonneau et al., 2013), increased detritism in the Rhône delta plain (Bruneton et al., 2001; Provansal et al., 2003), and increased Rhône River floods (Pichard, 1995). In contrast, CR5, the so-called Migration Period Cooling, is not characterized by any increase in Leptocythere, suggesting drier conditions in the Rhône watershed compared to CR4 and CR6.

The signature of CR0 (the 8.2 ka event), CR1, CR2, and CR3 is difficult to discriminate since the events are incorporated within a condensed section with very low accumulation rate (Fig. 8).

On the other hand, we also notice that two other periods of increased Leptocythere are recorded between ca.  300 and 230 cm and between ca. 180 and 120 cm, i.e., between ca. 1.4 and 1.3 ka cal. BP and ca. 1.2 and 1.0 ka cal. BP (Fig. 8). These periods are not related to global events (CR-like) but might correlate to periods of regional climate deterioration as attested by a high level of Lake Bourget (Arnaud et al., 2012), as well as by periods of increased detritism in the Rhône delta plain (Provansal et al., 2003).

Conversely, a strong decrease in Leptocythere is observed from 120 to 80 cm (Fig. 8). This suggests drier conditions during this interval, which corresponds to the MCA (ca. AD 950–1250). In the Northern Hemisphere, the MCA is generally described as a warm period characterized by intense dryness. In the Mediterranean, several studies have highlighted drier conditions during this event (e.g., Wassenburg et al., 2013; Martinez-Ruiz et al., 2015; Bassetti et al., 2016). The same signature is also observed in the Alps and the Rhône watershed, with periods of low lake level and low flood frequency, respectively (e.g., Magny, 2004; Wilhelm et al., 2016). Thus, the hypothesis of increased drought at the studied site during the MCA fits well with regional and local observations.

Within the late Holocene interval, the distribution pattern of Cluster E is slightly offset of the single species Leptocythere (Fig. 8). This cluster is constituted by the shallow infra-littoral Paradoxostoma and Loxoconcha species (Bonaduce et al., 1975; El Hmaidi et al., 2010) and by the epi-bathyal Argilloecia species. Argilloecia sp., in the Rhône subaqueous delta, appears to be tolerant to fluvial influence and potentially responds to organic matter supply (Fanget et al., 2013b). Nevertheless, an increase in Cluster E is not recorded within the periods characterized by higher river supply such as the CR4 and CR6 (i.e., the LIA) but slightly after these periods. This possibly indicates that Leptocythere is a better competitor during periods of increased detritism and fluvial discharge.

Conclusion

Our study shows that some environmental and sea-level changes during the Holocene can be clearly depicted from sedimentological and benthic meiofauna proxies.

During the early Holocene (11.7 to 7–8 ka cal. BP), sea-level rise led to the deposition of tempestite sediments that contain shallow infra-littoral benthic meiofauna. These deposits are thought to be formed between ca. 20 and 50 m water depth, and we believe that this feature can be used as a good regional scale paleobathymetric marker.

The middle Holocene (7–8 to 4–5 ka cal. BP) corresponds to a phase of very low sedimentation at the core site, resulting in the formation of a condensed section (i.e., the MFS in sequence stratigraphic terminology) reflecting the further landward position of the shoreline and Rhône outlet(s). This condensed section contains reworked shoreface material within a fine-grained matrix. It displays mixed faunal assemblages, ranging from infra-littoral to epi-bathyal environments, which are the result of erosion processes that occurred during the period of transgressive submersion and, thus, mark the peak of transgression and the subsequent sediment starvation.

Following the transgressive maximum, the late Holocene (4–5 ka cal. BP to 19th century AD) sediment deposits are influenced by a combination of allocyclic and autocyclic factors. The progressive shoreline progradation and prodeltaic lobes' switching is characterized by changes in clay mineralogy content, by the setting-up of benthic meiofauna adapted to thrive in the distal part of the Rhône River influence (i.e., distal St Ferréol and Ulmet lobes), and by the presence of very shallow-water species (i.e., proximal Grand Passon and Bras de Fer lobes).

Within the late Holocene deposits, ostracod assemblages emphasize fluctuations in the Rhône River hydrological activity. In particular, the occurrence of the ostracod genus Leptocythere highlights periods of increased fluvial discharge. These periods of intensified runoffs can be attributed to the 2.8 ka event (CR4) and the Little Ice Age (CR6) that are known to be at the origin of regional climate deterioration in western Europe, as well as periods of regional climate deterioration at ca. 1.3 and 1.1 ka cal. BP. In contrast, the signatures of the early and middle Holocene cold relapses are difficult to explore in the Rhône subaqueous delta since they correspond respectively to (a) a phase of rapid sea-level rise at the origin of shoreline reworking and deposition of tempestite and (b) a period of very low sedimentation at the core site resulting in a condensed interval with low temporal resolution.

Finally, our study demonstrates that prodeltas may provide interesting expanded archives of climate changes at the land/sea interface, with accumulation rates reaching 0.4 m yr${}^{-\mathrm{1}}$. On the other hand, such resolution can be achieved at one single site for only short time intervals, since depot centers migrate rapidly as a consequence of sea-level changes, high sediment fluxes, and lateral shifting of delta lobes. This highlights the need for acquiring series of long cores/boreholes, parallel and orthogonal to deltaic systems.

Data availability

For access to the proxy data, the respective authors should be contacted personally.

Acknowledgements

Core RHS-KS55 was collected during the RHOSOS cruise (2008) on board R/V Le Suroît. We thank the captain and crew of this cruise together with the Genavir technical staff and the scientific parties. Special thanks are due to Bernard Dennielou (IFREMER, Brest) for his commitment at sea and during the processing of the data in the laboratory. We thank the Laboratoire de Mesure du Carbone 14, UMS 2572, ARTEMIS in Saclay for ${}^{14}$C measurements by SMA in the frame of the National Service to CEA, CNRS, IRD, IRSN, and Ministère de la Culture et de la Communication. This work was partly supported by the CNRS-INSU “Mistrals-Paleomex”. We are grateful to Bertil Hebert (CEFREM, University of Perpignan), Gilbert Floch, Angélique Roubi, and Mickaël Rovere (IFREMER, Brest) for their technical support. We would like to thank A. Amorosi and L. Goisan for their comments and suggestions that helped to improve the manuscript. S. Luening is also thanked for commenting on the manuscript during the open discussion. Edited by: N. Combourieu Nebout Reviewed by: A. Amorosi and L. Giosan