Cold-Water Coral Mound Archive Provides Unique Insights Into Intermediate Water Mass Dynamics in the Alboran Sea During the Last Deglaciation

The Alboran Sea is widely recognized to host numerous cold-water coral ecosystems, including the East Melilla Coral Province. Yet, their development through time and response to climatic variability has still to be fully understood. Based on a combined investigation of benthic foraminiferal assemblages, foraminiferal stable isotope compositions, grain size analysis, sediment geochemistry, and macrofaunal quantification, this study identifies key events and processes having governed cold-water coral development at the East Melilla Coral Province between Greenland Stadial 2.1 and the Early Holocene. The transition from Greenland Stadial 2.1 to Greenland Interstadial 1 is associated to a decline of bryozoan communities and their replacement by cold-water corals, together with changes in benthic foraminiferal assemblages and a decrease in the sediment mean grain size. These results suggest that a rapid decrease in bottom currents and the establishment of dysoxic and mesotrophic conditions at the seafloor, possibly associated to enhanced fluvial input, resulted in the decline of bryozoans as the dominant suspension feeding organisms and their replacement by a thriving cold-water coral community. This transition from a bryozoan to a coral dominated environment is concomitant with the beginning of the African Humid Period, confirming that increasing fluvial input could have been a main factor triggering the establishment of cold-water corals in the East Melilla Coral Province during Greenland Interstadial 1. A change in benthic foraminiferal communities and an increase in the sediment mean grain size mark the passage from the Early to Late Greenland Interstadial 1. The current velocity of intermediate water masses is suggested to have increased during the Early to Late Greenland Interstadial 1, whilst simultaneously fluvial input would have reduced. Such changes suggest that the climate became more arid during the second phase of Greenland Interstadial 1.

The Alboran Sea is widely recognized to host numerous cold-water coral ecosystems, including the East Melilla Coral Province. Yet, their development through time and response to climatic variability has still to be fully understood. Based on a combined investigation of benthic foraminiferal assemblages, foraminiferal stable isotope compositions, grain size analysis, sediment geochemistry, and macrofaunal quantification, this study identifies key events and processes having governed coldwater coral development at the East Melilla Coral Province between Greenland Stadial 2.1 and the Early Holocene. The transition from Greenland Stadial 2.1 to Greenland Interstadial 1 is associated to a decline of bryozoan communities and their replacement by cold-water corals, together with changes in benthic foraminiferal assemblages and a decrease in the sediment mean grain size. These results suggest that a rapid decrease in bottom currents and the establishment of dysoxic and mesotrophic conditions at the seafloor, possibly associated to enhanced fluvial input, resulted in the decline of bryozoans as the dominant suspension feeding organisms and their replacement by a thriving cold-water coral community. This transition from a bryozoan to a coral dominated environment is concomitant with the beginning of the African Humid Period, confirming that increasing fluvial input could have been a main factor triggering the establishment of cold-water corals in the East Melilla Coral Province during Greenland Interstadial 1. A change in benthic foraminiferal communities and an increase in the sediment mean grain size mark the passage from the Early to Late Greenland Interstadial 1. The current velocity of intermediate water masses is suggested to have increased during the Early to Late Greenland Interstadial 1, whilst simultaneously fluvial input would have reduced. Such changes suggest that the climate became more arid during the second phase of Greenland Interstadial 1.

INTRODUCTION
Scleractinian cold-water corals (CWCs) shape worldwide, unique and highly diverse environments (Freiwald et al., 2004;Roberts et al., 2009). They settle in areas favorable for their suspension feeding life style, such as steep cliffs, rock overhangs or pre-existing seamounts, forming with time large CWC reefs (Freiwald, 2002;Roberts et al., 2009). The most prevalent scleractinian CWC species in the North Atlantic Ocean and Mediterranean Sea are Desmophyllum pertusum (formerly known as Lophelia pertusa (Linnaeus, 1758), see Addamo et al., 2016) and Madrepora oculata (Freiwald et al., 2004;Wheeler et al., 2007). Their development is promoted by high food availability coupled to enhanced hydrodynamic regimes, whilst the physico-chemical properties of the ambient seawater (e.g., density, salinity, temperature) also play a role (White et al., 2005;Mienis et al., 2007;Dullo et al., 2008;Davies et al., 2009;Davies and Guinotte, 2011;Hanz et al., 2019). A key characteristic of CWCs are their capacity to baffle suspended sediments entering their framework (Wilson, 1979;Mienis et al., 2009). The interplay between this sediment baffling and coral growth leads, over time, to the formation of CWC mounds . These mound structures are essentially confined to continental margins between 200 and 1000 m depth, where they can form build-ups of over 300 m in height (Mienis et al., 2007;Roberts et al., 2009;Wienberg and Titschack, 2015;Hebbeln et al., 2016Hebbeln et al., , 2019. The development of coral mounds depends primarily on the presence of dynamic bottom current regimes. Indeed, vigorous bottom currents increase lateral food supply, enabling corals to thrive, and support mound aggradation through sediment supply (White et al., 2005;Mienis et al., 2007;Hebbeln et al., 2016). The equilibrium between sediment and nutrient supply is fundamental for a coral mound to maintain sustained development through time. Although coral mounds demonstrate a discontinuous record, often punctuated by a number of hiatuses, they present much higher sedimentation rates than surrounding off-mound settings where sedimentation is limited (Rüggeberg et al., 2007;Eisele et al., 2008;Margreth et al., 2011;Raddatz et al., 2011;Stalder et al., 2018;Wienberg et al., 2018). As such, they possess a strong potential to recognize short term paleoenvironmental changes and may serve as unique high-resolution sedimentary records, preserving sediments that are otherwise eroded or non-deposited in neighboring open environments (Thierens et al., 2013).
Living scleractinian CWCs, mostly represented by M. oculata, are presently found in the Mediterranean Sea. They have, for example, been observed in the Northern Ionian Sea (Savini et al., 2014), in the Sicily channel , in the Tyrrhenian Sea (Taviani et al., 2017(Taviani et al., , 2019, in the Adriatic and Aegean Seas (Vafidis et al., 1997;Freiwald et al., 2009), in the Gulf of Lions (Zibrowius, 2003;Orejas et al., 2007;Fabri et al., 2014), and at the summit of ridges in the Alboran Sea Fink et al., 2013;Lo Iacono et al., 2014;Stalder et al., 2015Stalder et al., , 2018Wang et al., 2019). The oldest occurrences of fossil CWC mounds in the Mediterranean Sea date back to the Late Pliocene/Early Pleistocene and are related to the recolonization of the Mediterranean Sea by CWCs after the Messinian Salinity crisis (Taviani et al., 2005;Vertino et al., 2014). Mediterranean coral mounds are concentrated in the Alboran Sea, along the Moroccan margin van Rooij et al., 2013;Lo Iacono et al., 2014). The largest CWC mound field, the Melilla mound field, covers an area in excess of 500 km 2 parallel to the continental margin (van Rooij et al., 2013; Figure 1). It can be subdivided into two coral provinces, the West and East Melilla Coral Provinces, respectively 35 km northwest and 15 km northeast of the Cape Tres Forcas (Lo Iacono et al., 2014; Figure 1). Coral mounds situated in the West Melilla Coral Province dwell at depths ranging from 300 to 430 m and reach heights between 8 and 21 m above the seafloor (Lo Iacono et al., 2014). Coral mounds that are part of the East Melilla Coral Province demonstrate higher morphological variability. The mounds that are located to the North of the East Melilla Coral Province possess the particularity of being situated along three steep volcanic ridges (Brittlestar Ridges I, II and III). The summits of these Brittlestar Ridges are at ca. 250 m depth, whilst their base is found at ca. 450 m depth Lo Iacono et al., 2014;Wang et al., 2019;Figure 1). They are between 100 and 250 m wide and extend up to 3 km in length Fink et al., 2013).
The recent (last 20 ka) history of coral proliferation of the East Melilla Coral Province has been documented by a number of studies in the past years (Fink et al., 2013(Fink et al., , 2015Stalder et al., 2015Stalder et al., , 2018Wang et al., 2019;Feenstra, 2020). The conclusions drawn from these different studies indicate that corals experienced their most intense recent proliferation phase during Greenland Interstadial 1 (GI-1) and during the Early Holocene, whilst coral growth was reduced or halted during the Mid to Late Holocene, Younger Dryas and Last Glacial Maximum. Fink et al. (2013Fink et al. ( , 2015 and Stalder et al. (2015Stalder et al. ( , 2018 suggested that increased productivity was the main factor driving coral growth during GI-1 and the Early Holocene. Wang et al. (2019), based on high benthic foraminiferal accumulation rates, further suggested that high productivity is a decisive parameter sustaining mound development during GI-1 and the Early Holocene. Moreover, these authors suggest that a rearrangement of surface and intermediate water masses would also have triggered sustained coral growth at these times, whilst Feenstra (2020) further stress the importance of fluvial input and changes in water mass configuration for the initiation of coral growth during GI-1 and the Early Holocene. Based on these combined observations, it appears clear that GI-1 and its upper and lower boundaries are key periods for recent coral growth history in the Southeastern Alboran Sea. However, previous studies lack the sufficient resolution to enable a detailed investigation of the interaction between environmental changes and coral growth during GI-1. The aims of this study are thus: (1) to provide a high-resolution paleoenvironmental reconstruction of mound formation during the last deglaciation at the East Melilla Coral Province, with a particular focus on GI-1; (2) to assess in detail the variability of intermediate water masses during GI-1, and at the lower and upper GI-1 boundaries, in the Southeastern Alboran Sea. . The bold brown arrows illustrate the general circulation of the inflowing surface Modified Atlantic Water (MAW), which describes two anticyclonic gyres: the Western Alboran Gyre (WAG) and the Eastern Alboran Gyre (EAG). The Atlantic Anticyclonic Gyre (AAG) is illustrated by the thin brown arrows. The surface water circulation is based on a compilation of results from Lanoix (1974), Tintore et al. (1988), Viúdez andTintoré (1995) andL'Helguen et al. (2002), whist the bathymetry uses the GEBCO_2019 gridded bathymetric data. (C) Precise location of the Banc des Provençaux and Brittlestar Ridge I and II (BR I and BR II). The location of the studied core MD13-3455G is illustrated by the red dot, whereas the two yellow dots indicate the location of the neighboring cores MD13-3459G and MD13-3462G sampled during the same cruise. This map uses MultiBeam Echosounder data retrieved during the "GATEWAY" cruise (N • 194) onboard the Marion Dufresne II Research Vessel (van Rooij et al., 2013).

Geographical and Geological Setting
The Alboran Sea is the westernmost basin of the Mediterranean Sea and is approximately 400 km long and 200 km wide, with an average depth of 1300 m and a maximum depth of ca. 1800 m (Olivet et al., 1973;Comas et al., 1999). It is limited to the North by the Iberian Peninsula and to the South by North Africa. The Alboran Sea demonstrates an intricate seafloor morphology, with the presence of ridges, seamounts, plateaus, troughs and a number of sub-basins (Figure 1; Comas et al., 1999). This complex morphology is a result of the extension and subsidence that occurred in the region during the Early to Middle Miocene (Comas et al., 1999;Faccenna et al., 2004;Do Couto et al., 2016). The extension processes and associated crustal thinning led to the intrusion of the metamorphic basement by volcanic material (López Casado et al., 2001;Duggen et al., 2008). The shallow Banc des Provençaux, situated at ca. 200 m depth, is one of these volcanic edifices (Figure 1). Its limits extend out to form the three Brittlestar ridges (I, II and II) that are partly built up by CWCs and belong to the East Melilla Coral Province (Figure 1). Brittlestar Ridge I is situated at ca. 50 km from the mouth of the Moulouya river, which is one of North Africa's largest rivers (Emelyanov and Shimkus, 2012). The bases of the Brittlestar Ridges are surrounded by erosional moats, testimony to the strong currents affecting the area (Comas et al., 2009). Video surveys along the ridges revealed a diverse fauna consisting of numerous organisms, such as ophiuroids, soft corals and sponges. However, only scarce unregularly distributed living scleractinian coral colonies (D. pertusum and M. oculata) were observed, concentrated on the crest of the ridges ).

Oceanographic Circulation
The circulation patterns in the semi-enclosed Mediterranean Sea are governed by the exchange between inflowing low salinity Atlantic Water and outflowing highly saline Mediterranean Outflow Water at the Strait of Gibraltar. The Alboran Sea is dominated by three water masses. The surface water consists of low salinity (< 36.5 psu) Modified Atlantic Water (MAW) that extends down to approximately 200 m depth and becomes increasingly saltier as it flows eastwards (La Violette, 1983;Millot, 2009). High saline (ca. 38.5 psu) and warm (ca. 13.3 • C) Levantine Intermediate Water flows beneath the MAW between 200 and 600 m depth. It enters the Western Mediterranean through the Straits of Sicily to then exit through the Strait of Gibraltar, contributing to approximately 70% of the total Mediterranean Outflow Water (Millot, 2013). The Levantine Intermediate Water circulates along the Spanish margin, whilst the Southern Alboran Sea would be rather dominated by Shelf Water, i.e., a mixture of MAW and Western Mediterranean Deep Water (WMDW; Ercilla et al., 2016). Western Mediterranean Deep Water, colder (ca. 12.8 • C) and slightly fresher (ca. 38.4 psu) than the above lying Levantine Intermediate Water and Shelf Water, further supplies the Mediterranean Outflow Water and makes up the deepest water mass in the region (Millot et al., 2006;Millot, 2013). The WMDW circulates essentially along the Moroccan margin in the Alboran Sea (Hernández-Molina et al., 2002;Ercilla et al., 2016). The surface circulation in the Alboran Sea varies seasonally. The MAW enters the Northeastern Alboran Sea as a jet (1.6 Sv; 1 Sv = 10 6 m 3 .s −1 ) (Katz, 1972;Lanoix, 1974;La Violette, 1983). As a result of the entry of this jet, the surface circulation in the Alboran Sea describes two anticyclonic gyres, namely the Western Alboran Gyre and the Eastern Alboran Gyre (Heburn and La Violette, 1990;Lafuente et al., 1998; Figure 1). The strength of these two Alboran gyres fluctuates through time as a function of atmospheric and oceanographic parameters, such as the strength of westerly winds and the exchange rates between Modified Atlantic Water inflow and Mediterranean Outflow Water outflow (Heburn and La Violette, 1990). During winter, the Eastern Alboran Gyre tends to weaken and is replaced by a jet circulating eastward along the African coast (Tintore et al., 1988;Heburn and La Violette, 1990). A number of eddies develop along the gyres, such as the Almeria Cyclonic Gyre or the Atlantic Anticyclonic Gyre (L'Helguen et al., 2002; Figure 1). The East Melilla Coral Province is directly influenced by the Eastern Alboran Gyre and the Atlantic Anticyclonic Gyre (Lanoix, 1974;Viúdez and Tintoré, 1995;L'Helguen et al., 2002; Figure 1). Both gyres contribute to mixing between surface and intermediate water masses down to ca. 300 m depth (Heburn and La Violette, 1990;Lafuente et al., 1998).

Sample Collection
This study is based on the detailed investigation of the 491 cm long MD13-3455G gravity core ( Figure 1). A 1 cm slice of each half core section was taken every 5 cm for the investigation of benthic foraminiferal assemblages. The same 5 cm sampling interval was chosen to pick planktonic and benthic foraminifera for the analysis of their stable oxygen and carbon isotopic composition and for the particle-size analysis of the sediment matrix. The core was sampled every 10 cm to measure the Total Organic Carbon of the sediment matrix. All results were plotted with the ggplot2 package for R (Wickham, 2016;R Core Team, 2018).

Benthic Foraminiferal Assemblages
After drying and weighing, samples were wet sieved through a series of 200, 125, and 63 µm mesh sieves and dried at room temperature. Samples did not possess the same volume. The resulting fractions were then weighed separately. A minimum target number of 300 individuals was picked and identified per sample on the fraction > 125 µm. The sample was split with a microsplitter when necessary. Benthic foraminiferal densities were calculated by dividing the total number of foraminifera (number of individuals per aliquot multiplied by the split) by the dry fraction's weight. The benthic foraminiferal relative abundances were calculated for each sample and converted to percentages. The Shannon diversity index was calculated for each sample using the PRIMER6 software (Clarke and Gorley, 2006). Multivariate statistical analysis was performed using the software PRIMER6 (Clarke and Gorley, 2006). Foraminiferal composition data was root transformed in order to down-weigh the contribution of the most abundant species, thus increasing the contribution of rarer species (Field et al., 1982). Bray-Curtis similarity coefficients were then calculated (Bray and Curtis, 1957) and the resulting similarity matrix was used to perform a hierarchical cluster analysis. The hierarchical cluster analysis groups samples that share high similarity levels and is used to define benthic foraminiferal assemblages that are specific to a given depth range/epoch. A Similarity Percentage (SIMPER) analysis was computed at two different similarity levels in order to precisely identify the species responsible for the similarity among and the dissimilarity between clusters identified in the hierarchical cluster analysis (Clarke and Gorley, 2006).

Oxygen and Carbon Stable Isotope Analysis
A portion of the material prepared for the assessment of benthic foraminiferal assemblages was dry sieved through 250 and 212 µm mesh sieves. Between 5 and 12 specimens of the planktonic foraminifera Globigerina bulloides and the benthic foraminifera Cibicides lobatulus were picked from the 212-250 µm size fraction in order to measure their stable oxygen and carbon isotope compositions. This fraction was chosen to avoid any ontogenic effect on the measurements (Schiebel and Hemleben, 2017). Specimens were first cleaned three times with distilled water in an ultrasonic bath for 2 s to remove any attached particles or infilling. The measurements were then made using a Thermo Fisher Scientific GasBench II connected to a Thermo Finnigan Delta Plus XL isotope ratio mass spectrometer at the Stable Isotope Laboratory of the University of Lausanne according to a method adapted from Spötl and Vennemann (2003). Results are reported in the conventional δ-values in permil ( ) relative to the Vienna Pee Dee Belemnite (VPDB) standard. Analytical standard deviations (1σ) average 0.04 for δ 13 C and 0.06 for δ 18 O values based on 8 replicate analyses of standards in each sequence of 40 samples.

Grain Size Analyses
Grain size was determined on the bulk material and on the siliciclastic fraction using the Malvern Mastersizer 3000 at the department of Geology of the University of Ghent. Large clasts (e.g., coral and bryozoan fragments) were first hand picked out of each sample. Samples were then placed in 35% H 2 O 2 and boiled until the reaction stopped to remove organic matter. One half of each sample was boiled for 2 min in 10% HCl. This allowed dissolution of the carbonate material, hence only the siliciclastic fraction remained. Before measurement, all samples were placed in 2% sodium polymetaphosphate and boiled to assure total particle separation. A 2 mm mesh sieve placed above the receiver of the Malvern Mastersizer 3000 retained any remaining large particles. Each sample was measured three times and an average value was taken as the final result. The mean particle-size of the bulk and siliciclastic fractions, together with the mean sortable silt size (i.e., the mean of the 10-63 µm siliciclastic fraction; McCave et al., 1995), were calculated using the Rysgran Package for R (Gilbert et al., 2015;R Core Team, 2018). The mean sortable silt size served as a proxy for bottom water current velocity (McCave and Hall, 2006), in combination with changes in benthic foraminiferal assemblages. Indeed, using solely mean sortable silt size as a bottom current velocity indicator is possibly biased in CWC environments, since coral framework can locally reduce bottom current speed and lead to the deposition of fine material (Huvenne et al., 2009;Eisele et al., 2011;Fentimen et al., 2020), thus leading to an underestimation of bottom current velocities.

Organic Geochemistry
Total Organic Carbon (TOC) and Mineral Carbon (MinC) content were determined at the Laboratory of Sediment Geochemistry at the University of Lausanne using the Rock-Eval6 pyrolysis technique. The RockEval6 pyrolysis technique determines a Hydrogen Index (HI) and an Oxygen Index (OI), respectively corresponding to the quantity of pyrolyzable organic compounds relative to TOC and the quantity of CO 2 relative to TOC (Baudin et al., 2015). The graphical representation of HI vs. OI (i.e., the modified Van Krevelen diagram) has regularly been used to determine the type of kerogen and the origin of the organic matter present in samples, although the exactness of this method has been subject to discussion (Katz, 1983;Delvaux et al., 1989;van Krevelen, 1993;Tyson, 1995). The measured OI values in the samples presented in this study are extremely high and are probably the result of cracking at low temperature of poorly crystallized CaCO 3 (Baudin et al., 2015). Thus, the OI was not used to determine the type of kerogen. Instead, HI vs. Tmax was plotted in order to determine the kerogen type and the origin of the organic matter (Delvaux et al., 1989). Tmax being the temperature at which the maximum hydrocarbon yield resulting from kerogen cracking occurs (Baudin et al., 2015). Moreover, the high OI implies that TOC values are possibly overestimated (Behar et al., 2001).

Radiometric Dating
Uranium-series isotope measurements ( 230 Th/U) were carried out on 6 cold-water coral fragments (D. pertusum and M. oculata) using a multicollector inductively coupled plasma source mass spectrometer MC-ICPMS (Thermo Fisher Scientific Neptuneplus) coupled with a dissolver (Aridus I) at the Institute for Environmental Physics at Heidelberg University (Table 1). Prior to measurements, pristine-looking coral fragments were selected, physically cleaned by sand blasting and further chemically cleaned using a weak acid leaching. More details about the protocol can be found in Wefing et al. (2019) and in Frank et al. (2004).
Radiocarbon dating was performed at the ETH in Zürich on 5 benthic foraminiferal samples ( Table 2). The epibenthic foraminiferal species Cibicides lobatulus, Cibicides refulgens and Discanomalina coronata were picked in order to obtain between 4 and 10 mg of pure carbonate. The samples were first dissolved in phosphoric acid. The extracted CO 2 was then converted to graphite and measured by Accelerator Mass Spectrometry (AMS) technique using the MICADAS dedicated instrument  Ages are corrected for a reservoir age of 390 ± 80 years (Siani et al., 2000). (Synal et al., 2007). Results were corrected for 13 C and calibrated using the Marine13 calibration curve (Reimer et al., 2013) and the software OxCal v4.2.4 (Ramsey, 2017). A reservoir age of 390 ± 80 years was applied to all ages (Siani et al., 2000).

Chronology
The chronostratigraphic framework of core MD13-3455G was constructed based on absolute radiocarbon and Uranium-series isotope (U/Th) ages (Tables 1, 2) combined with the planktonic (G. bulloides) and benthic (C. lobatulus) oxygen isotope records. High planktonic and benthic δ 18 O values, respectively between ca. 2 to 3 and 3 to 4 , characterize the base of the core from 490 to 440 cm. Considering the U/Th and radiocarbon ages, this interval is interpreted as Greenland Stadial 2.1 as defined by Rasmussen et al., 2014 (GS-2.1; Figure 2). Greenland Interstadial 1 is defined as the period with lower planktonic and benthic δ 18 O values between 440 and 15 cm (Rasmussen et al., 2014). The following decrease in planktonic and benthic δ 18 O values (respectively from about 2 to 0.5 and from 3 to 1.5 ) from 15 cm to the top of the core represents the Holocene (Figure 2). The Younger Dryas would have possibly been eroded or not deposited, although a short excursion toward higher benthic δ 18 O values, respectively from about 2 to 3 , is situated between 25 and 15 cm. Low planktonic δ 18 O values (about 2 ) would however suggest that this interval does not correspond to the Younger Dryas and that these deposits are indeed absent in this core. Such discontinuity and the presence of hiatuses is a regularly observed feature in cores recovered from CWC mounds (Dorschel et al., 2005;Rüggeberg et al., 2007;Frank et al., 2009;Raddatz et al., 2011;Stalder et al., 2018). The resulting chronology hence indicates that a great part of the core, from ca. 440 to 15 cm, corresponds to GI-1 (Figure 2). Greenland Stadial 2.1 and Greenland Interstadial 1 can be related respectively to the time periods formerly known as the Oldest Dryas and the Bølling-Allerød. The latest chronology defined by Rasmussen et al. (2014) will be used throughout this study.

Core Description and Macrofauna
Scleractinian CWCs are absent during GS-2.1 (Figure 2). This interval is characterized by a dense accumulation of the bryozoan Buskea dichotoma associated to the brachiopod species Gryphus vitreus (Figure 2). The upper part of the core, corresponding to GI-1 and the Holocene, consists of framework building scleractinian CWCs (D. pertusum and M. oculata) embedded in a matrix of mixed, clay to silt sized carbonate/siliciclastic material (Figure 2; Feenstra, 2020). Sedimentary and macrofaunal variations observed in core MD13-3455G are reported in more detail by Feenstra (2020).

Relative Abundances and Diversity
A total number of 167 benthic foraminifera species were recognized (Supplementary Material). Benthic foraminiferal densities vary between 653 ind.g −1 during the Late GI-1 (32 cm) and 10560 ind.g −1 during the Early GI-1 (302 cm; Supplementary Material). The Shannon diversity index ( Figure 3A) reaches minimum values during the GS-2.1 (1.831) and maximum values during the Holocene (3.589). Higher benthic foraminiferal diversity also characterizes the Early GI-1, whilst the Late GI-1 (from 230 to 25 cm) is marked by lower values (Figures 3A,B). The distribution of the main benthic foraminifera species is illustrated in Figures 3A,B. As all Bolivinid species (Bolivina alata, B. difformis, B. dilatata, B. pseudoplicata, B. pseudopunctata, B. spathulata, B. striatula and B. variabilis) follow the same distribution pattern throughout FIGURE 2 | Stratigraphy, core description and macrofaunal content of core MD13-3455G. The core descriptions and macrofaunal content were previously reported by Feenstra (2020). The stratigraphy is constructed by combining the G. bulloides planktonic (light blue curve) and C. lobatulus benthic (dark blue curve) δ 18 O records ( VPDB), U/Th coral and calibrated median 14 C foraminifera ages (Tables 1, 2). The resulting stratigraphic column is given to the far right (according to Rasmussen et al., 2014). The benthic foraminiferal assemblages (BFA) identified at 47% similarity and at 62% similarity thanks to the hierarchical cluster analysis ( Figure 4) are indicated.

Benthic Foraminiferal Assemblages
The hierarchical cluster analysis allowed to distinguish 2 clusters at 47% similarity and 5 clusters and an outlier sample at 62% similarity (Figure 4). Each cluster relates to a specific benthic foraminiferal assemblage (BFA). At 47% similarity, clusters separate between a GS-2.1 assemblage (BFA.A) and a combined GI-1 and Holocene assemblage (BFA.B). The most contributing species to BFA.A are C. lobatulus (ca. 30%), A. stelligerum (ca. 12%) and T. angulosa (ca. 9%), while C. lobatulus (ca. 14%), C. laevigata (ca. 10%) and M. subrotunda (ca. 6%) contribute the most to BFA.B. More detailed results of the SIMPER analysis are given in Tables 3, 4. At 62% similarity, clusters separate between a GS-2.1 assemblage (BFA.b), a transitional GS-2.1/GI-1 assemblage from 437 to 447 cm (BFA.d), an Early GI-1 assemblage from 432 to 232 cm (BFA.f), a Late GI-1 assemblage from 232 to 17 cm (BFA.e) and a Holocene assemblage from 12 cm to the core top (BFA.c) (Figure 4). The 5 species contributing the most to each cluster are reported in Table 5, whereas the species contributing the most to dissimilarity between clusters are given in Table 6. Overall, the clustering confirms and provides more details about benthic foraminiferal distribution patterns described in section "Relative Abundances and Diversity."

Stable Carbon Isotopic Records
The δ 13 C values of the benthic C. lobatulus have a range between 0.7 during the Late GI-1 (157 cm) and 1.4 during the Late GI-1(337 cm), while the δ 13 C values of the planktonic G. bulloides range between −2.1 during the Late GI-1 (367 cm) and −0.8 during the Holocene (7 cm; Figure 5). The planktonic δ 13 C values have a higher variability compared to the benthic δ 13 C values (Figure 5). Average planktonic and benthic δ 13 C values are comparatively high during GS-2.1 (ca. −1.2 and 1.3 respectively) (Figure 5). The planktonic δ 13 C values have more variability during GI-1. Greenland Interstadial 1 can be divided into two phases; a first phase between ca. 440 and 210 cm and a second phase between 210 and 25 cm, respectively corresponding to average lower (ca. −1.6 ) and higher (−1.2 ) δ 13 C values. The benthic δ 13 C record can also be divided into the same two phases during GI-1, although in contrast to the planktonic δ 13 C values, the early phase corresponds to average higher benthic δ 13 C values and the later phase to average lower values. These two phases correspond to two distinct benthic foraminiferal assemblages: Frontiers in Marine Science | www.frontiersin.org FIGURE 3 | (A) Benthic foraminiferal Shannon diversity and distribution of the main benthic foraminifera species (expressed as percentages). The benthic foraminiferal assemblages (BFA) identified at 47% similarity and at 62% similarity thanks to the hierarchical cluster analysis (Figure 4) are illustrated. The stratigraphic column defined in Figure 2, the G. bulloides planktonic and C. lobatulus benthic δ 18 O records and the U/Th coral and calibrated median 14 C benthic foraminifera ages are further indicated. (B) Distribution of the main benthic foraminifera species (expressed as percentages). The benthic foraminiferal assemblages (BFA) identified at 47% similarity and at 62% similarity thanks to the hierarchical cluster analysis (Figure 4) are illustrated. The stratigraphic column defined in Figure 2, the G. bulloides planktonic and C. lobatulus benthic δ 18 O records and the U/Th coral and calibrated median 14 C benthic foraminifera ages are further indicated.
f and e (Figure 5). Planktonic and benthic δ 13 C values show an overall decoupled trend, noticeably during the Late GI-1 (between ca. 160 and 180 cm; Figure 5). An increase in planktonic δ 13 C values (from ca. −1.6 to −0.8 ) marks the transition toward the Holocene, also corresponding to a change in benthic foraminiferal assemblage (BFA.e to BFA.c). Benthic δ 13 C values remain relatively high and invariable during the Holocene (ca. 1.2 ; Figure 5).

Grain Size
The mean grain size of the bulk material varies between ca. 5 and 60 µm (Figure 5). The coarsest intervals are limited to GS-2.1 (ca. 60 µm), Late GI-1 (ca. 36 µm) and Holocene (ca. 40 µm). A number of coarser intervals are observed during GI-1 (Figure 5). The mean grain size of the siliciclastic fraction ranges from ca. 5.3 to 9.4 µm, i.e., in the fine silt size range (Figure 5). The transition from GS-2.1 to GI-1 is characterized by a decrease in the mean grain size of the siliciclastic fraction from ca. 7.0 to 5.3 µm. Greenland Interstadial 1 is divided into different phases based on the mean grain size of the siliciclastic fraction. Between 440 and 280 cm, mean grain size is low (ca. 6 µm). This first phase is followed by a slight increase from ca. 6 to 7 µm between 280 and 180 cm (Figure 5). An important coarsening takes place at 180 cm (from ca. 7 up to 9.4 µm). This coarser interval (average mean grain size of 9 µm) spans from 180 to 120 cm (Figure 5). At 120 cm, the mean grain size decreases to ca. 8 µm and remains high until the end of GI-1, although the latest part of GI-1 (from 50 to 15 cm) shows a slight increase in mean grain size (ca. 8.5 µm). The passage from GI-1 to the Holocone is marked by a decrease in mean grain size of the siliciclastic fraction, from ca. 8.6 to 6.2 µm (Figure 5). FIGURE 4 | Hierarchical cluster analysis obtained from the benthic foraminifera Bray-Curtis similarity matrix. Benthic foraminifera abundances were root transformed. Two clusters are distinguished at 47% similarity (A,B). At a higher similarity level (62%), the cluster analysis allows to isolate 5 clusters and one outlier sample. Each cluster corresponds to a specific benthic foraminiferal assemblage. The species contributing the most to the similarity among and the dissimilarity between clusters are given in Tables 3-6. TABLE 3 | List of the 5 most contributing species to the similarity within clusters A and B, which separate at 47% similarity on the hierarchical cluster analysis ( Figure 5). The average abundance, similarity, contribution (%) and cumulative contribution (%) are indicated for each species. The average similarity of the cluster is also given.
The sortable silt mean grain size ranges from ca. 16.9 to 21.6 µm (Figure 5). The trend in mean sortable silt matches the one observed for the siliciclastic fraction ( Figure 5). Highest values are recorded during GS-2.1 and the second part of GI-1 (between 180 to 120 cm and 30 to 15 cm), whilst lowest values relate to the first part of GI-1, corresponding to BFA.f. Decreases in the mean sortable silt characterize the transitions between GS-2.1/GI-1 and GI-1/HC, respectively corresponding to the transitions BFA.A/BFA.B and BFA.e/BFA.c (Figure 5). An increase that takes place at 180 cm, mirroring the trend noticed in the mean siliciclastic fraction, can be related to the change in The average abundance in each cluster, average dissimilarity, average dissimilarity/standard deviation ratio (Diss/SD), contribution (%) and cumulative contribution to the total dissimilarity are indicated for each species. A high Diss/SD ratio indicates that the considered species is a good discriminating species. The average dissimilarity between the two clusters is also given.
benthic foraminiferal assemblage (BFA e to f) that separates the Early and Late GI-1.

Organic Geochemistry
Total organic carbon (TOC) varies between 0.18 wt% during GS-2.1 (460 cm) and 0.84 wt% at the transition between GI-1 and the Holocene (20 cm; Figure 5). Lowest values are recorded during GS-2.1 (0.18 wt%) and during the Holocene (0.34 and 0.37 wt%). The transition from GS-2.1 to GI-1 is characterized by an increase followed by a gradual decrease in TOC from 0.77 to 0.44 wt% ( Figure 5). Low TOC values (ca. 0.6 wt%) are recorded for the Early phase of GI-1. It is followed by a gradual increase to higher TOC values, from ca. 0.6 to 0.8 wt%, toward the end of GI-1 (between 280 and 210 cm; Figure 5). This slight increase in TOC matches the first coarsening phase noticed in the sortable silt and siliciclastic fraction mean grain size, together with the change in benthic foraminiferal assemblage (Figure 5). Total organic carbon content remains unchanged and is relatively high (ca. 0.8%) until the end of GI-1. The HI vs. Tmax diagram provides information about the origin of the organic matter and its relation to changes in benthic foraminiferal assemblages (Figure 6). All samples, except for two outliers (15 and 460 cm), are restricted to HI and Tmax values situated respectively between 123 and 286 mgHC/g TOC), and 409 and 428 • C ( Figure 6B). All samples, except sample 15, plot in the mixed marine-terrestrial/terrestrial deltaic domain ( Figure 6A). The Greenland Stadial 2.1 samples, corresponding to BFA.A, are associated to the highest HI values (248 to 286 mgHC/g TOC) and are associated to the most marine related organic matter. The Late GI-1 (corresponding to BFA.e) and the Early GI-1 (corresponding to BFA.f) samples relate to respectively intermediate (188 to 270 mgHC/g TOC) and low (123 to 229 mgHC/g TOC) HI values ( Figure 6A). The Early GI-1 is associated with more terrestrial deltaic organic matter, whereas the Late GI-1 is rather affiliated to more mixed terrestrial/marine organic matter ( Figure 6A). Higher Tmax (426 and 428 • C) and intermediate (173 and 251 mgHC/g TOC) HI values characterize the Holocene (corresponding to BFA.c). Thus the organic matter during the Holocene appears to have a rather mixed terrestrial/marine origin, although the number of samples analyzed does not allow any precise interpretation.

Greenland Stadial 2.1
The Greenland Stadial 2.1 benthic foraminiferal assemblage (BFA.A) resembles the glacial assemblage described by Stalder et al. (2018) at BRI, with high abundances of the epibenthic C. lobatulus and D. coronata and infaunal C. laevigata (Figures 3A,B and Table 3). The filter feeding D. coronata and C. lobatulus preferentially cling to elevated substrata (e.g., coral branches, dropstones, hydroids) in well-oxygenated waters dominated by vigorous bottom currents (Linke and Lutze, 1993;Schönfeld, 1997;Margreth et al., 2009;Schönfeld et al., 2011). Together with the infaunal C. laevigata and H. balthica, they are typically associated with the availability of fresh, high quality marine organic matter (De Rijk et al., 2000;Schmiedl et al., 2000; TABLE 5 | List of the 5 most contributing species to the similarity within clusters b, c, d, e and f, which separate at 62% similarity on the hierarchical cluster analysis ( Figure 5). The average abundance, similarity, contribution (%) and cumulative contribution (%) are indicated for each species. The average similarity of the cluster is also given.  The average abundance in each cluster, average dissimilarity, average dissimilarity/standard deviation ratio (Diss/SD), contribution (%) and cumulative contribution to the total dissimilarity are indicated for each species. A high Diss/SD ratio indicates that the considered species is a good discriminating species. The average dissimilarity between each pair of clusters is also given. Milker et al., 2009;Stalder et al., 2018). This assemblage, together with the HI vs. Tmax diagram, confirms that the seafloor at BRI received the input of high quality marine derived organic matter during GS-2.1 (Figure 6). Surface productivity was important during GS-2.1, as suggested by the high planktonic δ 13 C values ( Figure 5). Based on Al-normalized elemental ratios, Feenstra (2020) suggest that arid continental conditions prevailed at the end of the last glacial period and that eolian dust was an important fertilizer for marine communities during GS-2.1. The presence of vigorous bottom currents is directly recorded by the high sortable silt mean grain size (Figure 5). These increased bottom currents would have been favorable for the large epibenthic foraminifera C. lobatulus and D. coronata, and the infaunal T. angulosa (Figures 3A,B and Table 3). Trifarina angulosa is reported from current-swept areas and is resistant to permanent winnowing (Mackensen et al., 1995;Schönfeld, 2002). The association of C. lobatulus, D. coronata, T. angulosa, Cassidulina spp. and CWCs is reported from the Porcupine Seabight (Margreth et al., 2009;Schönfeld et al., 2011;Fentimen et al., 2018), Norwegian margin (Spezzaferri et al., 2013;Stalder et al., 2014) and the Alboran Sea Stalder et al., 2015Stalder et al., , 2018. However, this assemblage is here observed in relation to a thick accumulation of the bryozoan Buskea dichotoma associated to the brachiopod G. vitreus, whilst CWCs are absent (Figure 2). Thus, despite seemingly favorable conditions such as high surface productivity and vigorous bottom currents, CWCs did not thrive during GS-2.1 (Figure 2). Bryozoans, CWCs and G. vitreus share the same need for sustained bottom currents, although bryozoans and G. vitreus prefer moderate velocities and do not cope with excessively strong bottom currents (Emig and Arnaud, 1988;Scholz and Hillmer, 1995;Bjerager and Surlyk, 2007;Roberts et al., 2009). Likewise, bryozoans and CWCs thrive when fresh organic matter is available, although CWCs are more tolerant to varying food sources than bryozoans (Best and Thorpe, 1994;Holbourn et al., 2002;Duineveld et al., 2004;van Oevelen et al., 2016). Lower temperatures coupled to lower sea level stand and increased upwelling promoted the development of bryozoan reef mounds at the Great Australian Bight (James et al., 2000;Holbourn et al., 2002). Similar cold conditions and low sea level stand during GS-2.1 at BRI are suggested by the high planktonic δ 18 O values and the presence of the brachiopod G. vitreus (Figure 2). Indeed, G. vitreus is found today along the Mediterranean continental margin, at depths between ca. 100 and 250 m (Emig, 1988). It is associated to "white coral" debris between 235 and 255 m depth near the Hyères Islands (SE France), where it reaches densities between 5 and 10 ind.m −2 (Emig and Arnaud, 1988). Below this depth range, G. vitreus becomes rarer and quickly disappears from the faunal community (Emig and Arnaud, 1988). Thus, these observations would suggest that the sea level during GS-2.1 was 75 ± 10 m lower than today. Additionally, A. stelligerum and H. balthica, two foraminifera species affiliated to the GS-2.1 assemblage (Figures 3A,B), are reported to prefer cold conditions (Husum and Hald, 2004;Murray, 2006). Hence, the presence of bryozoans and the absence of corals during GS-2.1 may be partly induced by lower temperatures at BRI. Overall, the environment during GS-2.1 is characterized by cold temperatures, high surface productivity and relatively vigorous but not excessively strong bottom currents. Increased bottom current velocity during glacial periods in the Western Mediterranean has been linked FIGURE 5 | Planktonic (G. bulloides) and benthic (C. lobatulus) δ 13 C values, Total Organic Carbon (%), mean grain size of bulk material (µm), siliciclastic (µm) and sortable silt fraction (corresponding to the 10 to 63 µm siliciclastic grain size range, expressed in µm). The benthic foraminiferal assemblages (BFA) identified at 47% similarity and at 62% similarity thanks to the hierarchical cluster analysis (Figure 4) are illustrated. The stratigraphic column defined in Figure 2, the G. bulloides planktonic and C. lobatulus benthic δ 18 O values and the U/Th coral and calibrated 14 C benthic foraminifera ages are provided as supporting information.

Less than 2 samples
to a stronger intensity of the Levantine Intermediate Water Toucanne et al., 2012). We suggest that such increased intermediate water flow occurred during the cold GS-2.1 at BRI (Figure 7). Moreover, Cacho et al. (2006) reported that deep water density and ventilation reached a maximum during MIS 2. Recently formed dense waters in the Mediterranean may lie at intermediate depths, above 1500 m, hence never reaching the deepest parts of the western Mediterranean (Sparnocchia et al., 1995;Millot, 1999;Ercilla et al., 2016). In addition, WMDW has been traced at depths shallower than 500 m along the Moroccan margin (Ercilla et al., 2016) and contributes significantly to the intermediate waters in the region. Thus, it is conceivable that shoaling of WMDW along the Moroccan margin during the last glacial period, especially during GS-2.1, led to a higher admixture of this water mass to overlying intermediate waters. The contribution of these nutrient-rich and well-ventilated waters at intermediate depths, together with an important input of eolian dust, would have favored the proliferation of bryozoans and associated epibenthic foraminifera.

The Early Greenland Interstadial 1
Transition Between Greenland Stadial 2.1 and Greenland Interstadial 1 The transition from GS-2.1 to GI-1 is marked by a significant change in the benthic macro and microfauna (Figures 2, 3).
The larger epibenthic C. lobatulus and D. coronata retreat, together with the infaunal T. angulosa and the cold water related H. balthica and A. stelligerum. Simultaneously, CWCs become the dominant macro-organism, to the detriment of bryozoans (Figure 2). The benthic foraminiferal assemblage at the transition between GS-2.1 and GI-1 (BFA.d) is characterized by high abundances of the epibenthic agglutinated S. wrightii, together with the small epibenthic hyaline G. praegeri and R. bradyi (Tables 3, 4). Spiroplectammina wrightii has been associated to modern CWC habitats dominated by dead corals at the   (Table 4). The selected foraminifera do not encompass the complete assemblage of each time period and are a simplified interpretation. The illustrated macrofauna and microfauna are not to scale.
Porcupine Seabight, NE Atlantic (Margreth et al., 2009). Its ecological requirements are not well constrained, however a related species, S. biformis, is known to tolerate severe dysoxia and to respond positively to increased nutrient availability in the Norwegian Drammensfjord (Alve, 1991;Gross, 2000). Szymańska et al. (2017) found it to be the most abundant genus in the Adventfjorden and Hornsund fjords of Svalbard, whilst it is also a common occupant of seep environments (Bernhard et al., 2001;McGann and Conrad, 2018). The epiphytic R. bradyi and G. praegeri have previously been observed together with CWCs at BRI (Stalder et al., 2015(Stalder et al., , 2018. Rosalina bradyi is frequently found along shallow shelf sites of the Mediterranean and dwells attached to elevated substrates, such as vegetation (Jorissen, 1987;Cimerman and Langer, 1991;Barras et al., 2014). It also shows a preference for warmer waters (Murray, 2006). Laboratory experiments have demonstrated that it responds positively to sporadic pulses of high amounts of organic matter and can tolerate extreme salinity changes and harsh oxygen depletion (Heinz et al., 2001;Fontanier et al., 2008). Thus, BFA.d indicates that the transition from GS-2.1 to GI-1 was marked by decreasing seafloor oxygenation and eutrophication. Rapidly decreasing bottom current velocities, inferred from the decrease in sortable silt mean grain size, would have further characterized the passage from GS-2.1 to GI-1. We link this decrease in bottom current velocity to a slowdown of intermediate water mass circulation (Figure 7). A substantial reduction in WMDW flow is observed in the deep western Mediterranean at 14.5 ka . We propose that the change in faunal composition and the abrupt decrease in sortable silt mean grain size at the transition between GS-2.1 and GI-1 is linked to this reduction of WMDW flow. A reduced WMDW flow would have weakened the ventilation of intermediate water masses. This reduced ventilation and current velocity at BRI may have in turn resulted in an important faunal turnover. These observations suggest that the decrease in the intensity of WMDW flow in the Southeastern Alboran Sea started as early as ca. 14.8 ka ago. Furthermore, the planktonic and benthic δ 18 O values support that surface and bottom waters also became warmer (Figure 2). This observation is supported by the disappearance of cold water related benthic foraminifera species (A. stelligerum and H. balthica) and the appearance of the warmer water dwelling R. bradyi (Figure 3).

Increased Terrestrial Input During the Early Greenland Interstadial 1
The environmental changes that happened at the GS-2.1/GI-1 transition intensified during the Early GI-1. High abundances of the infaunal Bulimina spp. and C. oolina, together with increased R. bradyi abundances, characterize BFA.f (Figure 3 and Tables 5, 6). Bulimina spp. is a typical genus for dysoxic and eutrophic environments and is dominant in the Mediterranean Sea near the Po and Rhône river deltas (Phleger and Soutar, 1973;Lutze and Coulbourn, 1984;Jorissen, 1987;Schmiedl et al., 2000;Mojtahid et al., 2009). Bulimina spp. is adapted to feed on degraded organic matter (De Rijk et al., 2000;Schmiedl et al., 2000). The intermediate or deep infaunal C. oolina has been described as a species adapted to suboxic conditions and feeding preferentially on degraded organic matter (Corliss, 1985;Sen Gupta and Machain-Castillo, 1993;Bernhard and Sen Gupta, 1999;Schmiedl et al., 2000;Fontanier et al., 2002). Elevated organic matter concentrations in the deeper parts of the sediment offer favorable conditions for the proliferation of this species (Rathburn and Corliss, 1994). Thus, the association of Bulimina spp., C. oolina and R. bradyi suggests that the eutrophication trend that began at the transition between GS-2.1 and GI-1 was reinforced throughout the Early GI-1. Moreover, BFA.f relates to a more terrestrial organic matter origin (Figure 6). Organicrich layers (ORLs), presenting a similar assemblage composition as BFA.f, have previously been described in the deep Alboran basins (Comas et al., 1996;Cacho et al., 2002;Rogerson et al., 2008). Organic-rich layers boast high TOC levels (up to 0.9 wt%), although noticeably lower than sapropel deposits in the eastern Mediterranean (> 2 wt% TOC) (Emeis et al., 2000;Cacho et al., 2002;Martínez-Ruiz et al., 2003;Rogerson et al., 2008). The sediment during the Early GI-1 at BRI does not have higher TOC content compared to GS-2.1. However, based on the benthic foraminiferal assemblage, we propose that it can still be related to ORL1, although not representing an ORL deposit per se.
The warming that started at the transition between GS-2.1 and GI-1 persisted throughout the Early GI-1, as indicated by the planktonic and benthic δ 18 O values (Figure 2). Higher abundances of R. bradyi during the Early GI-1 confirm this trend (Figure 3). Permanent low bottom current velocity during the Early GI-1, implied by the stable low sortable silt mean grain size, suggests that the weakening of intermediate water flow persisted throughout the Early GI-1 (Figure 7). Furthermore, the lower abundance of the epibenthic C. lobatulus, S. wrightii, G. praegeri and the infaunal T. angulosa in BFA.f in comparison to BFA.d may even indicate that the intermediate water mass flow decreased during the Early GI-1 ( Table 4). This would in turn suggest a persistent low or even dwindling WMDW flow.

Environmental Parameters Triggering Coral Proliferation
The African continent experienced important climatic changes during the Late Pleistocene and Holocene (Jolly et al., 1998;Gasse, 2000). A warm and humid period, known as the African Humid Period, generally reported to have lasted from 11 to 5.5 ka, had a major impact on marine and terrestrial systems, such as faunal and floral communities, as well as human migration pathways (Barth, 1857;Kuper and Kröpelin, 2006;deMenocal and Tierney, 2012 and references therein). East African lakes were tens or hundreds of meters higher than today, whilst large river systems flowed through the Western Sahara (Talbot and Livingstone, 1989;Skonieczny et al., 2015). The abrupt change from cold hyper-arid to warm humid conditions began as soon as 14.8 ka BP off Cap Blanc, Mauritania (deMenocal et al., 2000). The composition of BFA.f (Tables 3, 4) and the type of organic matter present in the sediment (Figure 6) could suggest that the East Melilla Coral Province experienced an increase in fluvial input during the GS-2.1/GI-1 transition and the Early GI-1. Moreover, the coral fragment dated at the transition GS-2.1/GI-1 has an age of 14.814 ± 0.052 ka (Table 1 and Figure 2). This age can be considered near to concomitant with the beginning of the African Humid Period as defined by deMenocal et al. (2000) at IODP Site 658C and the beginning of ORL1 (ca. 14.800 ka; Cacho et al., 2002;Jimenez-Espejo et al., 2008;Rodrigo-Gámiz et al., 2011;Martinez-Ruiz et al., 2015). Thus, this would confirm that the switch from cold arid to warm humid conditions over the North African continent, related to the beginning of the African Humid Period, led to a sharp increase in fluvial input at BRI. At the same time hydrodynamic changes were taking place in the Mediterranean Sea. The contribution of MAW in the Alboran Sea increases with rising sea level and deepening of the Gibraltar sill (Sierro et al., 2005). Simultaneously, changes in deep-water mass dynamics occurred. Rogerson et al. (2019) proposed that during the last deglacial, dense waters formed in the Gulf of Lions were injected at intermediate depths rather than at deeper depths. This would have isolated water masses above 2000 m from water bodies situated beneath this limit (Rogerson et al., 2019). Following the observations made by Ercilla et al. (2016), we suggest that this boundary may have been situated at shallower depths along the Moroccan margin, separating BRI from the influence of deeper well-ventilated waters. The important freshwater inflow from the North African continent, combined to the injection of more MAW in the Alboran Sea, further resulted in weak intermediate and deep-water ventilation and increased stratification (Sierro et al., 2005;Rogerson et al., 2008;Toucanne et al., 2012). Thus, the equilibrium at BRI appears to shift between a system primarily influenced by exchanges between intermediate and deep water masses during GS-2.1 to a system dominated by surface water dynamics and separated from the influence of deep water masses during the Early GI-1. The combination of rapidly increasing fluvial input, temperature and sea level rise, and the slowdown of intermediate and deep water mass flow, triggered a change in benthic habitats and the advent of a thriving CWC ecosystem (Figure 7). Weaker bottom water currents, combined to the high fluvial input and the presence of a dense coral framework, led to dysoxic conditions at the sediment/water interface at BRI (Figure 7). These observations confirm that CWCs can withstand and even thrive in weakly ventilated bottom waters if sufficient food is available.

The Late Glacial Interstadial 1
The Late GI-1 assemblage (BFA.e), in comparison to the Early GI-1 assemblage (BFA.f), is characterized by higher C. bradyi, C. laevigata, C. lobatulus, G. subglobosa and Bolivinids. abundances and lower Bulimina spp., G. praegeri, R. bradyi and S. wrightii abundances (Figures 3, 4 and Tables 5, 6). The opportunistic Bolivina spp. has an infaunal living strategy and dwells in highly productive areas associated to low seafloor oxygenation (Lutze and Coulbourn, 1984;Mackensen et al., 1995;Schmiedl et al., 2000;Murray, 2006). The association of Bolivina spp. and C. oolina suggests that the seafloor during the Late GI-1 experienced periods of oxygen depletion. However, in contrast with the Early GI-1, higher abundances of the epibenthic C. lobatulus and the infaunal C. laevigata, C. brady and G. subglobosa indicate that the seafloor was overall better oxygenated during the Late GI-1. The species G. subglobosa is regularly found along the Alboran platform and prefers high energy and well-oxygenated environments (Mackensen et al., 1995;Milker et al., 2009). It responds favorably and quickly to the input of fresh phytodetritus (Gooday, 1993;Suhr et al., 2003;Sun et al., 2006). It can be inferred from the higher C. lobatulus and G. subglobosa abundances and increased sortable silt mean grain size that bottom currents picked up during the Late GI-1 (Figures 3, 6). We propose that these observations can be related to an increase in velocity of intermediate water masses during the later phase of the GI-1, together with a re-establishment of deep water circulation in the Southeast Alboran Sea (Figure 7). In parallel to these changes taking place at intermediate and possibly deep water depths, surface waters would also have been affected by varying conditions. The origin of the organic matter changes from more terrestrial to more marine derived (Figure 6). The higher contribution of phytodetritus feeding benthic foraminifera, such as C. lobatulus and G. subglobosa, confirms that BRI received less terrestrial material and was experiencing an increased input of fresher marine derived organic matter. This change in the origin of the organic matter reaching the seafloor and in the benthic foraminiferal community indicates that fluvial input was overall reduced during the Late GI-1 (Figure 7). Such a reduction of fluvial input may be linked to a transition toward less humid and more arid conditions above the continent. We propose that the Late GI-1, as observed in core MD13-3455G, could represent a precursor phase to the arid conditions that prevailed during the Younger-Dryas in the region (Fink et al., 2015;Stalder et al., 2015). Furthermore, the slightly higher planktonic and benthic δ 18 O values during the Late GI-1 may suggest a minor cooling of surface and intermediate water masses. Coral proliferation phases are less sustained during the Late GI-1 in comparison to the Early GI-1 (Figure 2). Thus, the weaker fluvial input at BRI during the Late GI-1 had a negative effect on the coral community, although it did not lead to a full retreat (Figure 7). Moreover, reduced fluvial input in the Southeast Alboran Sea may have favored the re-establishment of intermediate and deep water mass circulation. It can be inferred from these observations that intermediate and deep water mass flows reached a minimum from ca. 14.8 to 14.3 ka, whilst it picked up after 14.3 ka (Figure 7).
Environmental conditions during the Late GI-1 were not stable. A short interval during the middle of the Late GI-1 (ca. 120 to 60 cm) corresponds to a decrease in sortable silt mean grain size and higher abundance of C. oolina (Figures 3, 6), thus implying a slowdown in bottom current velocity and an increased input of degraded organic matter. We propose that this time interval is associated to a short period of increased humidity over the continent and to reduced intermediate and deep water mass flow, similar to the Early GI-1, although of a shorter duration.

The Holocene
The Holocene is known to show significant short-term climatic variability in the Alboran Sea (Bond et al., 2001;Frigola et al., 2008;Ausín et al., 2015;Fink et al., 2015). However, the reduced size of the Holocene horizon in core MD13-3455G does not allow a precise interpretation of changes having occurred during this period. Thus, we can only track the general environmental trends that characterized the Holocene at BRI. These would correspond to major rearrangements in water masses, which led to important changes in benthic communities. Rapidly increasing surface and bottom water temperatures during the Holocene are inferred from the decrease in planktonic and benthic δ 18 O values and lower H. balthica abundance compared to the Late GI-1 (Figures 2, 3). This matches previous observations made in the Alboran Sea (Cacho et al., 2001;Ausín et al., 2015). This general warming at shallow and intermediate water depths is associated to a shift in the coral community, M. oculata becoming the dominant species instead of D. pertusum, and in the benthic foraminiferal assemblage (Figure 2). This change from D. pertusum to M. oculata dominated environments is documented at a number of different areas in the East and West Melilla Coral Provinces (Fink et al., 2013(Fink et al., , 2015Stalder et al., 2015Stalder et al., , 2018Wang et al., 2019;Feenstra, 2020), and is hence a widespread phenomenon. The decrease in sortable silt mean grain size together with the retreat of the epibenthic foraminifera species C. lobatulus indicates that the seafloor experienced a decrease in bottom current strength during the Holocene. We relate this decrease in bottom water currents, comparable to the passage from GS-2.1 to the Early GI-1, to a reduced flow of intermediate water masses. Reduced rates of Levantine Intermediate Water flow into the Alboran Sea during interglacials have previously been documented in the Corsica Through and in the Sicily Channel (Incarbona et al., 2011;Toucanne et al., 2012). In contrast, the influence of MAW during interglacials increases as a result of sea level rise (Sierro et al., 2005). The increase in MAW flow could explain higher abundances of the phytodetritus feeding infaunal G. subglobosa and U. mediterranea. The species U. mediterranea is considered opportunistic in the Mediterranean Sea, responding positively to the important supply of fresh organic matter (De Rijk et al., 2000;Schmiedl et al., 2000;Fontanier et al., 2006). Thus, increasing MAW inflow would favor surface productivity during the Holocene, in turn promoting the development of phytodetritus feeding foraminifera at the seafloor.

CONCLUSION
A range of measurements made on the MD13-3455G gravity core demonstrate that a combination of short-term fluctuations in continental and marine systems control coral development and demise at the East Melilla Coral Province. Furthermore, it is likely that fluvial input and nutrient supply are the primary environmental factor driving cold-water coral proliferation in the Southeast Alboran Sea, whilst they may influence the dynamics of intermediate and deep water masses in the region. An important change in benthic communities took place at the transition from Greenland Stadial 2.1 to Greenland Interstadial 1. Cold conditions, lower sea-level and stronger bottom currents during Greenland Stadial 2.1 would have favored the development of an ecosystem dominated by the bryozoan Buskea dichotoma associated to the brachiopod Gryphus vitreus. These environmental conditions were tightly linked to enhanced intermediate and deep water mass flow. Rapidly increasing fluvial input and surface and bottom water temperatures, coupled to a weakening of bottom currents, led to the decline of bryozoans and their replacement by cold-water corals during the Early Greenland Interstadial 1. Changes in the benthic foraminiferal fauna, such as the lower abundances of large epibenthic (e.g., C. lobatulus and D. coronata) and increased numbers of infaunal species (e.g., Buliminids and C. oolina), also reflect the environmental changes from Greenland Stadial 2.1 to Greenland Interstadial 1. We propose that this switch is related to the beginning of the African Humid Period, establishing more humid and warmer conditions over the continent and initiating continental runoff. Simultaneously, the flow of intermediate and deep water masses would have rapidly decreased in the Southeast Alboran Sea. Warmer bottom waters and fluvial related eutrophication were detrimental for bryozoans but proved optimal for cold-water coral growth. Benthic foraminiferal assemblages suggest that the sediment/water interface was poorly oxygenated during the Early Greenland Interstadial 1.
During the mid-Greenland Interstadial 1, bottom currents strengthened, clearly separating this period between an early and a later phase. At the same time, the origin of the organic matter becomes more marine, whilst the benthic foraminiferal fauna is subject to compositional changes. We suggest that the impact of fluvial input on benthic communities was reduced during the later phase of Greenland Interstadial 1. More arid continental conditions together with a re-establishment of dynamic intermediate and deep water mass circulation would have characterized this period. However, our results stress that conditions were unstable at the seafloor during the Late Greenland Interstadial 1, with periods of increased fluvial input resulting in eutrophication and low oxygen availability. In contrast, bottom currents were weaker during the Holocene whilst water temperatures rapidly increased. This led to the development of renewed microfaunal and coral communities, M. oculata becoming the dominant species instead of D. pertusum.

DATA AVAILABILITY STATEMENT
All the datasets used in this study are available at the open-access repository PANGEA using the following link: doi.pangaea.de/10.1594/PANGAEA.914702.

AUTHOR CONTRIBUTIONS
RF acquired the data, interpreted the data, wrote and coordinated the manuscript. EF and TA participated in the data acquisition. AR participated in data interpretation and in writing the manuscript. TV and IH participated in data acquisition, interpretation and in writing the manuscript. DV head organizing member of the research cruise, participated in writing the manuscript. AF participated in data interpretation and in writing the manuscript, provided financial support.
the Eastern Mediterranean: cold-water coral ecosystemsand seeps" and FN-200021_149247 "4D-Diagenesis@Mound: Understanding thetemporal and spatial variability of early diagenesis in carbonate mounds" for funding this research. Mr. Efraim Hall is thanked for his help with the CT quantification of the bioclasts and sedimentary core descriptions. We also thank Prof. Norbert Frank (Institute of Environmental Physics, University of Heidelberg, Germany) for performing and providing the Uraniumseries isotopemeasurements. We are grateful for the ship time provided by IPEV on the R/V Marion Dufresne II within the framework of the EuroFLEETS GATEWAYS project (grant agreement 228344). We further thank Dr. Tim Collart for the help provided with the Rysgran Package for R and Prof. Sébastien Bertrand (Department of Geology, Ghent University, Belgium) for the assistance with the grainsize measurements.