Early- to mid-Holocene paleohydrology in northeast Mediterranean: The detrital record of Aliakmon River in Loudias Lake,Greece

Abstract

One of the prominent features of northeast Mediterranean (NEM) Holocene climate are recurrent phases of cold and aridity; their impacts on the hydrological cycle remain at large unknown, as few existing paleohydrological records are either restricted to lake-level fluctuations or focus near the ‘8.2 kyr BP’ event. Here, we present the detrital record of Aliakmon River in Lake Loudias between 9500 and 3000 cal. BP. Magnetic susceptibility (MS) exhibits high correlation with mean grain size (r = 0.7) of silt-sized fractions and is used as a proxy of the distal clastic input of Aliakmon River, whereas organic matter (OM), carbonate content (CaCO₃), water content (WC), and clay concentrations decipher sedimentological and biological processes in Lake Loudias. Periods of high hydrological activity were interrupted by short intervals of low river discharge at c. 9400, 8600, 8200, 7500, 7000, 6200, 5300, and 4500 cal. BP and during a multi-century event centered at 3500 cal. BP, in agreement with marine and terrestrial paleoclimatic reconstructions from NEM. With exception of the wet period between c. 8.6 and 7.9 cal. BP, periods of increased hydrological activity are synchronous to contraction of Pinus forests and increased sea surface temperatures and silt transport in the Aegean Sea. The long-term (~580 years) variability of MS is in-phase with southeast Europe pollen-inferred annual temperature variations and with North Atlantic Ice Rafted Debris (IRD) events.

Keywords

Aliakmon River climate fluvial clastic input Greece Holocene Ice Rafted Debris events northeast Mediterranean paleohydrology

Introduction

A wide array of proxies derived from terrestrial and marine archives has shown that during the Holocene, the climate of eastern Mediterranean experienced prominent cold and arid deteriorations, of multi-centennial to millennial duration, mainly associated with a weakening of the polar vortex and variations in solar irradiance (e.g. Marino et al., 2009; Mayewski et al., 2004; Rohling et al., 2002; Schmiedl et al., 2010). The climate of eastern Mediterranean is also known to be affected by low-latitude (monsoonal) and high-latitude (North Atlantic) atmospheric forcing (Lionello et al., 2006a). To a smaller geographical extent, the subregion of northeast Mediterranean (NEM and also defined as Northern Borderlands of Aegean Sea (NBAS) by Kotthoff et al., 2008a), bounded by the transboundary watersheds of southern Balkan rivers to the north and by the North Aegean Sea to the south (Figure 1), appears to show much faster response to large-scale atmospheric forcing (Tzedakis et al., 2009). Numerous Holocene paleoclimatic and paleoenvironmental reconstructions are available mainly from marine (e.g. Ehrmann et al., 2007; Gogou et al., 2007; Hamann et al., 2008; Kotthoff et al., 2008a, 2008b; Kuhnt et al., 2007; Marino et al., 2009; Schmiedl et al., 2010; Triantaphyllou et al., 2013, 2016) sedimentary sequences (Figure 1), constituting NEM as one of the best studied areas, in terms of Holocene oceanographic dynamics. On the terrestrial realm, past environmental conditions are mainly inferred from lacustrine settings (e.g. Digerfeldt et al., 2007; Francke et al., 2013; Lawson et al., 2005; Panagiotopoulos et al., 2012; Peyron et al., 2011; Pross et al., 2009; Zhang et al., 2014), while sparse studies have focused on other archives, such as fluvial aggradation sequences (e.g. Macklin et al., 2010) and the geoarchaeological record (e.g. Berger et al., 2016; Ghilardi et al., 2012, and references therein). Reconstructed variables from NEM major lacustrine systems (Lake Ohrid, Lake Prespa, Lake Dorjan, Tenaghi Philippon peatland, and Lake Xinias; Figure 1), also considered in this study, include lake-level fluctuations (e.g. Digerfeldt et al., 2007; Francke et al., 2013), palynological records, and pollen-derived precipitation and temperature reconstructions (e.g. Lawson et al., 2005; Panagiotopoulos et al., 2012; Peyron et al., 2011; Pross et al., 2009), as well as organic, inorganic, and isotopic sediment proxies (e.g. Francke et al., 2013; Zhang et al., 2014). These lacustrine systems involve karstic lakes that lack sedimentary input from major fluvial systems, and even though they provide valuable records of environmental change, the signals of past changes can be biased from local processes, such as temperature inversions and interannual variability (Wu et al., 2007). However, sedimentary archives of freshwater lakes along river deltas, situated on the fringe between terrestrial and marine environments, are representative of processes that take place in the alpine foreland where fluvial systems originate and exhibit continuity with the ocean shelf and margin stratigraphic sequences. The fluvial systems of NEM (Evros River, Nestos River, Strymonas River, Axios River, Aliakmon River, Pinios River, and the rivers of the northwest coast of Turkey; Figure 1), involve relatively small watersheds (<40.000 km²) and drain areas of contrasting lithology and climatic regimes. Nowadays, all northeastern Mediterranean rivers discharge directly to the north Aegean Sea, but during the late glacial and Holocene, Aliakmon and Axios rivers were discharging into a shallow marine embayment that now comprises the Plain of Thessaloniki. Massive deposition of deltaic sediments during the Holocene led to the formation of bird-foot-type deltas of the major rivers and of alluvial-fan-type deltas of minor rivers. Delta progradation caused a gradual isolation of the bay from the sea, resulting in temporarily estuarine conditions and finally to the formation of a backshore freshwater lake, Lake Loudias.

Figure 1.

Map of the main fluvial systems, mountain ranges, Aliakmon River watershed (patched in gray color), and ocean circulation patterns of northeastern Mediterranean (NEM) region and of the sites cited in the text. Marine cores SL-148, SL-152, MNB-3, SL-21, M-4G, and MD01-2430 that have been used for the paleoenvironmental reconstruction of northern Aegean Sea and the Sea of Marmara are also shown.

In this paper, we present the early- to mid-Holocene detrital record of Aliakmon River in Lake Loudias. Three sediment cores (NN3, NN4, and S1) recovered from the southern part of the lake (Figure 2) are used for stratigraphic correlations and for the establishment of the age–depth relationship, while selected physical sediment properties of core S1 are used to reconstruct the depositional history of Aliakmon River. The cores were extracted during a multi-year geoarchaeological project that aimed to evaluate the geomorphological changes in the western part of Thessaloniki Plain (Ghilardi et al., 2008a, 2008b) and the environmental evolution of the lake-dwelling Neolithic settlement of Nea Nikomedeia, which represents the oldest and most important Neolithic settlement in northern Greece (Ghilardi et al., 2012). The detrital input of Aliakmon River is expressed by the magnetic susceptibility (MS) record, whereas the environmental variability of Lake Loudias is described by the records of organic matter (OM), calcite, water, and clay contents. The overall objective of this work is to reconstruct the early- to mid-Holocene hydrological history of Aliakmon River through the study of sediment properties and to correlate the novel results with paleoclimatic and paleoenvironmental records from the broader NEM region.

Figure 2.

General setting of the study area illustrated on combined LANDSAT imagery and SRTM 90 topography (elevation scale: red color: 30–20 m, yellow color: 20–10 m, and blue color: 10–0 m). Holocene extents of Lake Loudias (black dashed line), Aliakmon River delta (white dashed line), and main lobes (yellow color) are visible, together with the respective location of cores NN4, S1, and NN3. Core S1 has been presented in Ghilardi et al. (2008a), while cores NN4 and NN3 and additional information of core S1 in Psomiadis et al. (2014). Cutoff meanders and late-Holocene abandoned channels, as well as the locations of the Neolithic settlement of Nea Nikomedeia, the ancient city of Vergina and the former Macedonian capital of Pella, are also shown.

Study site

Aliakmon River originates on the northwest part of continental Greece; the river catchment constitutes a west-east trending sigmoidal feature with surface area of natural basin (in the 1920s, 2200 km² were artificially added to Aliakmon River drainage; Styllas, 2016) of 6100 km² (Figure 1). The highest altitude of Aliakmon River watershed is 2520 m, and its northern tributary originates just south of Prespa Lake (Figure 1). Upper basin (>600 m elevation) lithology is composed of felsic and mafic rocks (24%), but soft sediments (flysch–molasse, neogene, and quaternary deposits) are the most abundant formation (60%), resulting in high Holocene sediment yields that is according to the estimates of delta accommodation space reduction, average 6.6 × 10⁶ t/yr (Styllas, 2017). Carbonate formations (16%) occupy the lower part of Aliakmon drainage and the entire catchment of Vermion Mountain (Figure 2).

The climate of Aliakmon watershed and delta undergoes continental and marine influence, respectively, with a strong seasonal cycle. Average annual precipitation at the limit between upper and lower catchment is 765 mm, average annual temperature is 12°C, while annual evaporation is 340 mm, representing 44% of annual precipitation (period of observations AD 1961–2000, data source: Greek Public Operation of Electricity). Annual precipitation in the delta plain is 490 mm (data source: Hellenic National Meteorological Society, Station Trikala Imathias), while average annual discharge at the river delta is 70.6 m³/s (period observations from 1963 to 1999, data source: Greek Public Operation of Electricity) and peaks during December to March, in contrast to a double precipitation peak occurring during October–December and April–June (Styllas, 2016). Winter precipitation patterns of Aliakmon River drainage are related to western and central Mediterranean cyclonic activity, while spring and summer precipitation is mainly controlled by local convective processes. During winter, moisture-laden air masses from Ionian Sea are pushed onto Pindos Range under the influence of westerlies, resulting in considerable amounts of precipitation along Aliakmon River headwaters (Xoplaki et al., 2004). This pattern is observed during negative phases of the North Atlantic Oscillation (NAO), expressed by the negative correlation coefficient (r = −0.43) between winter (December through March (DJFM)) NAO index (Hurrell and NCAR Staff, 2016) and Aliakmon River average annual discharge. Monthly correlations between NAO index and monthly discharge become more significant by the end of the winter season, especially during March (r = −0.6), as rain-on-snow events along higher watershed elevations, combined with increasing spring temperatures, contribute to high water and sediment discharge in the delta.

The Holocene delta of Aliakmon River has an area of 570 km² and is part of the broader Thessaloniki Plain (Figure 2). Lake Loudias existed since the Pleistocene and was more related to the coalescence of Pleistocene fans of Almopia and Aliakmon rivers (Figure 2). It was also located close to Mount Vermion marginal travertines with additional runoff from subterraneous springs, while its more recent form resulted by the confluence of Axios and Aliakmon main lobes since the late glacial. The initially shallow freshwater lake became gradually deeper during the early stages of Holocene between 9000 and 8400 cal. BP (Ghilardi et al., 2008a, 2012). From 8400 to 7800 cal. BP, the vicinity of Lake Loudias was intruded by the rising sea level, a process that resulted in locally brackish conditions, which together with the globally recognized cold and arid ‘8.2 ka RCC event’ (Alley et al., 1997) likely forced the abandonment of the Neolithic settlement of Nea Nikomedeia (Ghilardi et al., 2012). Lake Loudias persisted until mid-Holocene (~4000 cal. BP), when massive fluvial sedimentation resulted in the gradual infilling of the lake and the formation of small lacustrine bodies and marshes. During late-Holocene, Lake Loudias was transformed into an extensive marsh (Ghilardi et al. 2012; Psomiadis et al., 2014), to be artificially drained in the 1920s. The southern shore of Lake Loudias receives sediments mainly from Aliakmon River and to a lesser extent from small rivers draining Vermion Mountain. The location of cores S1, NN4, and NN3 (Figure 2) is representative of distal fluvial sedimentation in a low-energy environment (Ghilardi et al., 2008b, 2012; Psomiadis et al., 2014).

Analytical methods

The 675-cm-long core S1 was recovered from the margin of Aliakmon River delta (40°38′40″N, 22°20′51″E). The one half-split of the core was analyzed for a number of physical sediment parameters in the sedimentology laboratory of Bjerknes Centre for Climate Research (BCCR), University of Bergen, Norway. The other half-split of the core was analyzed in the Centre Européen de Recherche et d’Enseignement en Géosciences de l’Environnement (CEREGE); the results are presented in Ghilardi et al. (2008b, 2012) and Psomiadis et al. (2014).

Methods applied to the core in BCCR involved MS measurements (contiguous 1 cm throughout the core) carried out by a Bartington MS2E sensor. The weight loss-on-ignition (LOI), dry and wet bulk density, and water content (WC) were measured at 1-cm intervals between 100 and 675 cm of core length. The LOI was determined by continuously retrieving volumetric samples of 1 cm³ and drying them overnight at 105°C, before the dry weight was measured. The weight LOI (%) was calculated following the methodology of Dean (1974).

Grain-size distributions of the sediment were realized at a 5-cm interval, as there was not adequate material left following the retrieval of sediment for LOI and bulk density determinations, realized by a Micromeretics Sedigraph 5100 Particle Size Analysis System (Webb and Orr, 1997). Approximately 2 cm³ of sediment from each 5-cm sampling interval was wet sieved through a 125-µm mesh with the use of hexametaphosphate dispersing agent. The analyses were run on a grain-size range between 1 and 63 µm (coarse silt to clay) as sand fractions in the remaining material had minor (1–2%) contributions.

The chronology is based on a synthetic age–depth model of five previously published accelerator mass spectrometry (AMS) ¹⁴C dates, from cores S1 and NN4 (Table 1). Sediment samples with peat, OM, plant remnants, and charcoal were considered, as deposition of organic rich layers within the short distance (5.5 km) between the two cores in low-energy conditions is expected to form under similar environmental forcing. The correlation between the cores was also based on their MS stratigraphy.

Table 1.

Radiocarbon dating results used for the construction of core S1 composite age–depth model (from Ghilardi et al., 2012, and Psomiadis et al., 2014).

Core ID	Sample type	Altitude relative to mean sea level (m)	Laboratory reference	Age ¹⁴C BP	cal. BP
NN4	Plant	+1.68	SacA-21947	2930	3140 ± 30
NN4	Organic matter	+0.13	SacA-21951	5095	5885 ± 30
NN4	Charcoal	−2.4	Poz-34965	6750	7650 ± 40
NN4	Organic matter	−2.7	SacA-21952	8095	9020 ± 35
S1	Peat	−2.35	Poz-16762	8400	9462 ± 50

Results

Core lithology

The middle and upper parts of core S1 are homogeneously composed of silty clay sediments, intercalated with coarse silt and fine sand thin layers. Observations exclude the existence of hiatuses and/or major flood events, suggesting low-energy sedimentary deposition. In total, four distinct lithostratigraphic units (U1–U4) were defined from bottom to top, on the basis of sediment color variations, grain size, MS, OM, calcium carbonate (CaCO₃), WC, and clay concentration downcore trends (Figure 3).

Figure 3.

Synthetic lithological description of core S1 with the stratigraphic units; their associated sedimentary environments and downcore variations of magnetic susceptibility (MS); clay weight (%); and the contents (%) of organic matter (OM), carbonate (CaCO₃), and water (WC). The early-Holocene period of the record is characterized by fluvial and/or alluvial conditions as illustrated from the high MS values (Psomiadis et al., 2014).

Lithostratigraphic unit U1 (675–600 cm) is composed of fluvial sediments of brown color, with fine sand/coarse silt (56–68 µm) horizons alternating with layers of silt and clay. Unit 1 shows high MS (90.5 × 10⁻⁵ SI) and low OM (3.9%) values, negatively correlated with each other (Figure 3 and Table 2). U1 exhibits high correlation between MS and mean grain size (Psomiadis et al., 2014; Figure 7) that implies that the increase in the magnetic signal results from the dominance of Aliakmon River detrital fraction in Lake Loudias sedimentation. High values of MS and calcite (CaCO₃ > 10%) are positively correlated with each other (r = 0.40; Table 2) and coincide with OM and clay minimum values. The up-core decrease in MS is indicative of a changing environment from fluvial to lacustrine and to lower energy depositional regime.

Table 2.

Correlation coefficients between selected physical parameters for stratigraphical units A, B, C, and D. Bold numbers correspond to statistical significant correlations.

Unit A	MS	OM	CaCO₃	WC	Unit B	MS	OM	CaCO₃	WC
MS	–	−0.67	0.40	−0.57	MS	–	−0.51	0.02	−0.32
OM		–	0.13	0.85	OM		–	−0.22	0.64
CaCO₃			–	0.00	CaCO₃			–	0.23
WC				–	WC				–
Unit C	MS	OM	CaCO₃	WC	Unit D	MS	OM	CaCO₃	WC
MS	–	−0.53	−0.48	−0.18	MS	–	−0.28	0.12	−0.03
OM		–	0.70	0.57	OM		–	−0.38	0.37
CaCO₃			–	0.21	CaCO₃			–	0.02
WC				–	WC				–

MS: magnetic susceptibility; OM: organic matter; WC: water content.

Unit 2 (600–320 cm) has a light gray color and constitutes the lacustrine part of the record consisted of mainly silty sediments, with variable clay content and interbedded centimeter-scale layers of fine sand (Figure 4). Coarsening upward (CU) trends denote increases in sediment transporting energy and correspond to increases in grain size. The highest concentrations of clay (28%), WC (50%), and OM (5.0%) of the entire record are indicative of a calm sedimentary environment. OM is negatively correlated (r = −0.22) with CaCO₃ and positively correlated with WC and clay contents (Table 2 and Figure 3).

Figure 4.

(a) Illustration of a half-split meter of core S1 (5–5.5 m depth) with coarse silt layers indicated (black circles), together with the MS profile for the same depth, and the graphic relation between (b) MS and (c) mean grain size for stratigraphic units 2 and 3.

A marked change in color from light gray to a light brown sediment with dark clay horizons, between U2 and the overlying U3, attests shallower conditions and the transition from lacustrine to floodplain depositional environment. Unit U3 (320–150 cm) demonstrates the highest carbonate content (10.5%), which is positively correlated with OM (Table 2), likely a result of enhanced biological productivity from water column stratification that led to precipitation of authigenic calcium carbonate and deposition of OM. High CaCO₃ values can also be attributed to carbonate clastic sediment deposition from the streams of Vermion Mountain.

The upper part of the floodplain deposits is characterized by the formation of marshes (U4, 150–100 cm), with dark gray fine-grained layers that display high OM (5%), the lowest MS (10.2 × 10⁻⁵ SI) values, and pronounced peaks of calcite. MS is negatively correlated with CaCO₃ (r = −0.48), suggesting a biological origin of calcite and/or increasing influence of the streams from Vermion Mountain (Figure 2). MS peaks are also negatively correlated with sediment WC (%), as clastic inputs during periods of increased river flow reduce the available space for water to be stored.

MS as a proxy for detrital input

The MS signal of fluvial sediments depends on the relative amounts of diamagnetic (e.g. calcite and quartz, MS negative), paramagnetic (clay minerals, MS weak and positive), and ferromagnetic minerals (magnetite and titanomagnetite, MS high and positive) minerals (Arnaud et al., 2005). Figure 4 presents a half-meter section of core S1 between 500 and 550 cm (Figure 4a), together with the MS profile (Figure 4b). The coherence between silt concentrations and MS is also evident both from the graphic (Figure 4c) and statistical correlation (r = 0.7), between mean grain size and MS.

Given the high correlation between MS and grain size and the structure of Aliakmon River catchment, which consists of a ‘carbonate buffer’ below 600 m, we assume that the magnetic fractions transported to the delta to be representative of the upper watershed hydrological conditions. Increased MS values in Lake Loudias sediments can be the result of high precipitation and flooding activity, vegetation changes, poor soil development and enhanced erosion, or from a combination of these factors. Low MS values are related to decreasing hydrological activity of the upper watershed and to either increasing clastic input from the lower Aliakmon watershed and the small rivers that drain Vermion Mountain, as detrital carbonate fractions exhibit low MS values or from deposition of authigenic calcite during warm periods with stagnant waters that also dampen MS concentrations.

Lake Loudias age–depth relationship

The age–depth relationship was established using a composite chronology contemplated by stratigraphic cross-correlations of core S1 with the lithostratigraphic columns of cores NN4 and NN3, as presented in Ghilardi et al. (2012) and Psomiadis et al. (2014), but were also refined here from their MS profiles (Figure 5a). The age–depth relationship for Lake Loudias was then constructed by fitting a linear regression line between the radiocarbon dates (Figure 5b). The ages of core NN3 (Figure 5a) were ignored because of the absence of characteristic MS peaks along the middle part of the core (3.8–4.9 m of sediment depth) where available datings of cores NN4 and S1 exist.

Figure 5.

(a) MS correlation of cores NN4, S1, and NN3, showing the tie points of pronounced magnetic susceptibility (MS) peaks, used to correlate the Lake Loudias sedimentary sequences, and (b) the radiocarbon dates and composite age–depth linear relationship for Lake Loudias.

Discussion

A 6500-year record of Aliakmon River detritism

For the purpose of this work, we select the time interval between 9500 and 3000 cal. BP because it is constrained chronologically by the age model, it represents the lacustrine part of the Lake Loudias record and is characterized by low-energy conditions and fine-grained deposits. Given the considerations of the previous sections, we assume MS to be a reliable proxy of Aliakmon River clastic load in Lake Loudias. MS has been successfully used as a proxy for flooding activity in the case of Rhone River sediments in Lake Bourget (Arnaud et al., 2005; Chapron et al., 2005; Debret et al., 2010). MS shows a long-term increasing trend from 9500 to 4500 cal. BP (Figure 6), related with increasingly early- to mid-Holocene wet conditions in agreement with high lake levels, of selected lakes across Greece and the Balkans between 8000 and 5500 cal. BP (Digerfeldt et al., 2007). In contrast to MS, OM exhibits a long-term decreasing trend from 8800 to 3000 cal. BP. Periods of increased hydrological activity and high clastic input cause mixing of the water column and reduce lake productivity because of increased turbidity (Wolfe et al., 2005), and/or simply result in the dilution of OM, explaining the negative correlations between MS and OM. The relative amount of OM increases under low (MS minima) hydrological periods, except for two occasions at c. 6200 and 4400 cal. BP, when both proxies simultaneously decrease (Figure 6). WC (%) of sediments is fluctuating near a constant value of 50% between 8800 and 5000 cal. BP, decreasing from then onward, reflecting the transition from lacustrine to floodplain sedimentary environment. Prominent lowerings of WC correspond to increased MS values as a response to the increase in coarser fractions in lake sediments. CaCO₃ is nearly constant throughout the record (~10%), with characteristic peaks at 9400 and 3500 cal. BP, respectively.

Figure 6.

Core S1 magnetic susceptibility (MS), together with the percentages (%) of organic matter (OM), calcite (CaCO₃), and water content (WC). Gray vertical bars correspond to MS minima and are associated with periods of minimal fluvial activity. Significant environmental and climatic incidents from NEM terrestrial and marine archives are also noted.

Despite their long-term trends, all physical parameters demonstrate considerable variability. MS minima represent Aliakmon River phases of low sediment transporting capacity and occur at c. 9400, 8600, 8200, 7500, 7000, 6200, 5300, and 4500 cal. BP and during a multi-century event around c. 3500 cal. BP (Figure 6, gray bars). The periods of reduced discharge recorded in Aliakmon River sedimentary record, match terrestrial and marine pollen records of aridity from NEM, within dating uncertainties. For example, the North Aegean pollen sequence of core SL-152 records centennial-scale decreases in Quercus and nonsaccate arboreal pollen (AP) at c. 9300, 8600, 8100, 7400, 6500, 5600, 4700, and 4100 cal. BP (Kotthoff et al., 2008a), while south Aegean core LC-21 shows the dominance of cold species benthic foraminiferal from 8600 to 8200 cal. BP and at c. 6200 and 3500 cal. BP, respectively (Rohling et al., 2002), that coincide with Aliakmon River periods of low discharge.

Characteristic peaks in OM (6.5%), calcite (CaCO₃ > 20%), and MS minima occur approximately at c. 9400 cal. BP (Figure 6), in timing with north Atlantic ‘9.3ka – Bond 6’ event (Bond et al., 1997) that in Tenaghi Philippon was expressed by reduced annual precipitation (Peyron et al., 2011). This short-term interval of Aliakmon River low flows but increased CaCO₃ and OM was likely associated with warm and humid conditions that caused intense thermal stratification of Lake Loudias water column and led to enhanced precipitation of authigenic calcite in the epilimnion and to deposition of OM in the hypolimnion (Figure 6). Such conditions are better explained by increased summer precipitation in the lower watershed with high lake levels and increased temperatures since wet winter conditions in Aliakmon River drainage are associated with cool conditions, extensive snow cover, and high annual discharge. Accordingly, in Lake Prespa, the sedimentary record shows a remarkable increase in δ¹⁸O_carb, after 9400 cal. BP, which has been attributed to increased summer precipitation (Leng et al., 2010), while in Sidari archaeological site located on the island of Corfu (western Greece), increased gully erosion and vertical aggradation rates during the same period are also explained by intense summer precipitation events (Berger et al., 2016). The early-Holocene arid episode in Aliakmon River catchment was followed by a phase of high hydrological activity between 9400 and 8700 cal. BP that fits well with phases of high water levels in Lake Dorjan between 9500 and 8700 cal. BP (Francke et al., 2013). In north Aegean Sea, this wet phase coincides with anoxic bottom conditions that led to the formation of sapropel S1a (Kotthoff et al., 2008a, 2008b; Schmiedl et al., 2010; Triantaphyllou et al., 2016), which was largely fostered by enhanced fluvial transport of NEM rivers, as manifested by the distribution of clay minerals along many North Aegean sites (Ehrmann et al., 2007; Roussakis et al., 2004).

Following the early-Holocene wet phase in Aliakmon catchment, a rapid decrease in hydrological activity between 8800 and 8600 cal. BP is combined with a prominent decrease in CaCO₃ that reached record minima (1.5%) at approximately 9000 cal. BP, attesting increasingly colder conditions and fluvial transport incompetence. The continental scheme is in-phase with the beginning S1a interruption in the northern Aegean Sea (Kotthoff et al., 2008a, 2008b; Rohling and Pälike, 2005; Schmiedl et al., 2010; Triantaphyllou et al., 2016). Overall, the period from 8600 to 7900 cal. BP is characterized by increasing fluvial deposition in Lake Loudias (Figure 6) that was briefly interrupted at c. 8200 cal. BP but reached maxima at 7900 cal. BP. This short-term phase of reduced fluvial activity is associated with the ‘8.2 ka’ event (Alley et al., 1997), but for Aliakmon River, watershed appears to be of moderate intensity. On the contrary, the impacts of the ‘8.2 ka’ event appear to be of higher magnitude in the terrestrial archives of Tenaghi Philippon, expressed by a reduction in winter and summer temperatures by 4°C and 2°C, respectively (Pross et al., 2009), and to decreases in annual precipitation by 100–150 mm (Peyron et al., 2011), although marine core SL-152 pollen data suggest a weaker cooling (Kotthoff et al., 2008a).

The gradual increase in Aliakmon River fluvial activity between 8600 and 7900 cal. BP is coupled by high sediment transporting energy in Lake Prespa (Panagiotopoulos et al., 2013), while the recorded aridity along the eastern part of NEM (i.e. Tenaghi Philippon) is further supported from vegetation changes in southern Bulgaria, with the dominance of coniferous forest between 8500 and 7800 cal. BP in higher elevations of Rila Mountains (Tonkov et al., 2013). Increasingly wet conditions of Aliakmon River between c. 8600 and 7900 cal. BP are consistent with high lake levels in European mid-latitudes (period 12 of Magny, 2004). Situated just north of 40°N, the Aliakmon River hydrological regime shows similarities with the middle zone of the hydrological tri-partition model of Europe proposed by Magny et al. (2003). This contrasting hydrological pattern raises the question of a climatic/hydrological zonation across NEM with a western, NAO-influenced part including the Aliakmon River watershed and the Prespa Lake and with an eastern part mainly influenced by climatic patterns from northeast Europe and Asia (i.e. Siberian High (SH)) that includes Lake Dorjan, Nisi Fen, and Tenaghi Philippon; a situation not very different from present conditions.

The 8600–7900 cal. BP wet phase terminated within a multi-centennial time interval, as Aliakmon River activity decreased from 7700 to 7300 cal. BP (Figure 6). This phase of decreasing fluvial activity likely corresponds to the ‘7.4 ka’ event, which is linked to a weakening in North Hemisphere summer insolation (Bond et al., 2001). This distinct cool and arid episode that fits well with Aliakmon River reduced activity has been reported in Lake Ohrid (Vogel et al., 2010a) and the northern Aegean Sea (Kotthoff et al., 2008b; Triantaphyllou et al., 2016).

Overall, wet mid-Holocene conditions prevailed between 7000 and 4500 cal. BP across NEM, as inferred from Lake Dorjan high water and trophic levels that resulted from local enhanced catchment runoff and nutrient erosion (Francke et al., 2013; Zhang et al., 2014). Our MS record suggests that the wet period between 7400 and 6100 cal. BP is characterized by short periods of low fluvial activity at c. 7000 and 6200 cal. BP, alternating with multi-centennial wet periods. The period of minimal clastic input at c. 6200 cal. BP coincides with very low OM values and with minor CaCO₃ and WC changes (Figure 6), pointing to arid conditions and probably to the absence of littoral vegetation. The highly variable hydrological regime of Aliakmon River between 7400 and 6100 cal. BP is manifested by other regional records of environmental change. The occurrence of extensive alluviation episodes in the vicinity of Mount Olympus between 7800 and 6200 cal. BP, with the most voluminous alluviation episode taking place between 6200 and 5900 cal. BP (Krahtopoulou and Veropoulidou, 2014), is indicative of intense hydrological changes and overall wetness. Pollen-based winter precipitation and temperature reconstructions from North Aegean core SL-152 exhibit considerable variability as well, characterized by a cooling trend, with lowest temperature and highest precipitation values occurring between 6100 and 5400 cal. BP (Kotthoff et al., 2008b), in agreement with Aliakmon River periods of increased clastic deposition.

The period characterized by the most intense hydrological variability of Aliakmon River watershed is between 5800 and 4500 cal. BP and coincides with pronounced climatic and fluvial instability in Tenaghi Philippon region between 5800–5700, 5400, and 5000–4900 cal. BP, respectively, as recorded in Aggitis River paleoecological data (Lespez et al., 2016), and with warm and humid conditions in North Aegean between 5500 and 4000 cal. BP (Triantaphyllou et al., 2013). This period of high-frequency detrital variability was followed by a hydrologically active period between 4200 and 3750 cal. BP (Figure 6). In southern Greece, higher frequency of erosive rainfall events confined between 4800 and 4200 BP was traced both in Peloponnese (Pope and Van Andel, 1984) and in the island of Crete (Macklin et al., 2010). Finally, the most recent period of low hydrological activity in our record is the multi-century dry event at c. 3500 cal. BP, expressed by MS minima and peaks in OM and CaCO₃ concentrations, which is similar to the ‘9.3 ka’ event, denoting warm temperatures, intense thermal stratification of the lake’s water column, and/or increased clastic influence of Aliakmon lower watershed and Vermion Mountain small rivers. In local scale, this phase of aridity is synchronous with reduced water levels, increased influx of soil material, and high MS values in Nisi Fen (Lawson et al., 2005) and matches the beginning of the 3500–2500 BP globally recognized Rapid Climate Changes (Mayewski et al., 2004).

Timing with climatic events in NEM

In contrast to the well-established marine dynamics, there is a complete lack of continuous reconstructions of the terrestrial system hydrological changes during the Holocene in NEM. The detrital record of Aliakmon River in Loudias Lake is one of the first continuous high-resolution records of hydrological variability, consistent with several key events registered in terrestrial and marine records across NEM. More specific, periods of increased clastic input and high MS values show an opposite sign with Pinus pollen concentrations from marine core SL-152 (Figure 7a and b). Increased regional dryness expressed by the expansion of Pinus forests in Rhodope and Rila mountains generally coincides with periods of low discharge. One of the most pronounced droughts of Aliakmon River watershed at c. 5300 cal. BP is concomitant with the north Aegean Pinus pollen highest values (Figure 7a and b). Minor time discrepancies between core S1 arid intervals and the marine response lie within age uncertainties.

Figure 7.

Comparison between the hydrological record of Aliakmon River from core S1 (this study) with selected proxy records from northeastern Mediterranean region and the southern Aegean Sea: (a) Aliakmon River magnetic susceptibility (MS), (b) Pinus pollen counts from marine core SL-152 (Kotthoff et al., 2008a), (c) 12-point running mean of organic matter (OM, continuous line, this study) and total organic carbon (TOC, line with dots) of Lake Loudias and Prespa Lake (Panagiotopoulos et al., 2013), (d) sand/silt concentrations of marine core SL-148 (Hamann et al., 2008), (e) sand/silt concentrations of core S1 (this study), and (f) warm versus cold species benthic foraminiferal abundances (100% dominance of warm species, lower percentages correspond to cooler sea surface temperatures) of south Aegean marine core LC-21 (Rohling et al., 2002). Gray vertical bars correspond to MS minima and are associated with periods of minimal fluvial activity, similar to Figure 6.

The temporal variability of Lake Loudias OM and Prespa Lake total organic carbon (TOC) exhibits similar behavior for the first part of the record until 6500 cal. BP, despite their different physical characteristics (i.e. surface area, input sources, and water depth) and the temporal resolution of the considered records (Figure 7c). In both cases, production and/or preservation of organic carbon is reduced during periods of increased hydrological activity of Aliakmon River, likely due to intense water column mixing. Even though Prespa Lake lacks a major fluvial input and is highly karstic, high amounts of precipitation, catchment runoff and seepage flow, changes in vegetation cover, and cold temperatures can be potential causes of water mixing, bottom ventilation, and reduced amounts of TOC. After 6500 cal. BP Lake Loudias is more affected by local fluvial processes and continuous reduction of the water depth (transition from U3 to U4), while Lake Prespa was still affected by regional climate forcing, but local processes are not discarded. The concentrations of Aliakmon River fine sand and silt fractions match closely the sand and silt concentrations in marine core SL-148 sediments. In both records, sand is a minor contributor to the signals (1–2%). The close resemblance between fluvial and oceanographic processes (Figure 7d and e) reveals a link between upland sediment dynamics and deltaic deposition with the north Aegean continental margin. Pulses of silt on the ocean bottom (core SL-148; Hamann et al., 2008) may be the response of the marine system to detrital inputs and/or to a thermohaline reactivation with higher bottom current velocities.

A lagged response between high fluvial energy and marine sediment transport dynamics cannot be ruled out given the extensive low-gradient shelf area between the river mouth and core SL-148 position. The overall coherence of the sand/silt pulses enhances the possibility of a common forcing mechanism for Aliakmon River and north Aegean shelf dynamics since core SL-148 location is influenced from Thermaikos Gulf hydrodynamic circulation, especially when surface circulation is controlled by north winds (Karageorgis and Anagnostou, 2003).

The timing of Aliakmon River hydrological activity with the mid-latitudes and NAO-like patterns appears to be partially decoupled from the southern Aegean Sea thermal regime, as recorded by the ratio of warm versus cold benthic foraminiferal in marine core LC-21 (Rohling et al., 2002). Lower assemblages of warm species (Figure 7f) have been related to an increase in northerly outbreaks of cold and dry air in NEM, through an intensification of the SH (Marino et al., 2009; Rohling et al., 2002). However, during these cold stages, the Aliakmon River does not always show decreased detrital activity, as, for example, occurs around the ‘8.2 kyr event’. This in turn implies that periods of increased hydrological activity of Aliakmon River associated with negative phases of NAO may have coincided with periods of enhanced SH activity. Such large-scale complex atmospheric circulation is not a contradictory condition, as recent geochemical studies from high-latitude lakes in Finland suggest NAO and SH to be operating individually and to not entirely block each other during mid-Holocene (Saarni et al., 2016).

MS cyclicity and external forcing

A final thing to consider is the cyclic, quasi-periodic variation of the Aliakmon River MS signal that appears to have a period of approximately 500 years (Figure 8). If the present teleconnection between the North Atlantic and Mediterranean wet season precipitation (Xoplaki et al., 2004) has been active throughout early- to mid-Holocene, the observed hydrological cyclicity of Aliakmon River may be linked processes of similar period (~500 years) in the North Atlantic. As shown in Figure 8, the periods of Aliakmon River increased activity are in-phase with pollen-derived annual temperature anomaly variations in southeast Europe (Davis et al., 2003) and with the Ice Rafted Debris (IRD) episodes in North Atlantic (Bond et al., 1997).

Figure 8.

Correlation between the (a) 12-point running mean (~580 years) MS record of Aliakmon River with (b) SE Europe pollen-reconstructed annual temperature anomalies compared with the AD 1890–1990 reconstructed temperatures (Davis et al., 2003) and (c) the record of North Atlantic stacked IRD events from Bond et al. (2001). Time lags between the signals range between 250 and 500 years, a temporal discrepancy that lies within the age uncertainties of the records.

The origin of the IRD episodes during the Holocene is not straightforward (Mayewski et al., 2004), and the atmospheric responses to the variation of solar irradiance can be multiple, affecting the North Atlantic drift ice regime, the North Atlantic thermohaline circulation, and the Arctic oscillation (AO)/NAO system (Bond et al., 2001; Chapman and Shackleton, 2000). Despite the series of globally short-term cold and arid events during the Holocene, the long-term increases of Aliakmon River hydrological activity appear to be coupled with intervals of high IRD activity and increasing temperatures in SE Europe, while periods of low river discharge mostly correspond to intervals of air temperature decline and low IRD activity. The striking similarity of the curves illustrated in Figure 8 confirms the existence of the teleconnection between NEM and the northern latitudes in sub-millennial timescales during early- to mid-Holocene, but further conclusions will require additional high-resolution hydrological proxies with better chronological constraints and a larger spatial extent.

Conclusion and perspectives

A study based on sediment properties was used to provide new insights on the paleohydrology of northwestern Greece. The continuity of core S1 and its ability in recording the deposition of Aliakmon River sediments in Lake Loudias, at a distal location of major delta depocenters, proves to be a key element for deciphering the river’s early- to mid-Holocene hydrological history. Of equal importance is the fact that the upper (>600 m elevation) catchment lithology is dominated by igneous rocks and sedimentary formations (84%), which contain highly magnetic minerals, whereas the lower basin is composed by carbonate formations (16%) that exhibit low magnetism. The Aliakmon River MS record represents sediment fractions transported to the delta from the upper watershed during stages of increased fluvial discharge. Consequently, the evolution of Aliakmon River detritism between 9500 and 3000 cal. BP constitutes a reliable proxy for Holocene paleohydrology of NEM region. Periods of enhanced hydrological activity during wet climatic phases were interrupted by short periods of low discharge at c. 9400, 8600, 8250, 7500, 7000, 6200, 5300, 4500, and 3500 cal. BP, in agreement with reconstructions from marine and lacustrine records from across NEM. These periods coincide with cold and arid intervals, triggered from variations in solar irradiance and large-scale atmospheric circulation, as is also indicated here by the in-phase variation of Aliakmon River hydrological activity with southeastern Europe annual temperature variations and the North Atlantic IRD events.

The main advantage of utilizing the Aliakmon River detrital record is that according to our composite, age–depth model provides high sedimentation rates (1 cm–50 years), which in the absence of high-resolution terrestrial paleoenvironmental records opens new perspectives on the timing and climate forcing mechanisms of NEM hydrology. In order to examine the interaction of contrasting large-scale climatic patterns during these phases of cold and aridity, future studies in Lake Loudias should involve (a) better chronological constraints, (b) sampling locations closer to Aliakmon River lobes that will ensue even higher sedimentation rates, (c) use of geochemical proxies and trace-elemental analyses for direct comparisons with the MS signal, and (d) use of biomarker isotopic (δD and δ¹⁸O) compositions, which will bring new information on basin-scale vegetational changes, permit direct comparisons with isotopic signals of marine archives, and provide more accurate information about relative precipitation sources that drive the observed hydrological variability.

Footnotes

Acknowledgements

Mike Styllas is particularly thankful to Professor Atle Nesje of University of Bergen/Bjerknes Centre for Climate Research who ‘opened’ his eyes toward the use of sediment physical parameters in lacustrine environments as paleoclimate proxies.

Funding

The majority of the laboratory work presented in this paper was realized as part of an EU Marie Curie early-stage researcher fellowship. Radiocarbon dating was funded within the framework of Matthieu Ghilardi’s PhD (2004–2007) and by the ARTEMISCNRS French National Programme for 14C dating. This article is also part of the PALEOMEX-INEE program from CNRS and was funded by the ARCHEOMED workshop (Directors Laurent Carozza and Laurent Lespez).

References

Alley

Mayewski

Sowers

et al . (1997) Holocene climate instability: A prominent widespread event 8200 years ago. Geology 25: 483–486.

Arnaud

Revel

Chapron

et al . (2005) 7200 years of Rhone River flooding activity in Lake Le Bourget: A high-resolution sediment record of the NW Alps hydrology. The Holocene 15(3): 420–428.

Berger

Lespez

Kuzucuoğlu

et al . (2016) Interactions between climate change and human activities during the early to mid-Holocene in the eastern Mediterranean Basins. Climate of the Past 12: 1847–1877.

Bond

Kromer

Beer

et al . (1997) A pervasive millennial-scale cycle in North Atlantic Holocene and glacial climates. Science 278: 1257–1266.

Bond

Kromer

Beer

et al . (2001) Persistent solar influence on North Atlantic climate during the Holocene. Science 294: 2130–2136.

Chapman

Shackleton

(2000) Evidence of 550-year and 1000-year cyclicities in North Atlantic circulation patterns during the Holocene. The Holocene 10(3): 287–291.

Chapron

Arnaud

Noël

et al . (2005) Rhone River flood deposits in Lake La Bourget: A proxy for Holocene environmental changes in the NW Alps. Boreas 34: 404–416.

Davis

BAS

Brewer

Stevenson

et al . (2003) The temperature in Europe during the Holocene reconstructed from pollen data. Quaternary Science Reviews 22: 1701–1716.

Dean

Jr (1974) Determination of carbonate and organic matter in calcareous sediments and sedimentary rocks by a loss-on-ignition: Comparison with other methods. Journal of Sedimentary Petrology 44: 242–248.

10.

Debret

Chapron

Desmet

et al . (2010) North western Alps Holocene paleohydrology recorded by flooding activity in Lake Le Bourget, France. Quaternary Science Reviews 29: 2185–2200.

11.

Digerfeldt

Sandgren

Olsson

(2007) Reconstruction of Holocene lake-level changes in Lake Xinias, central Greece. The Holocene 17(3): 361–367.

12.

Ehrmann

Schmiedl

Hamman

et al . (2007) Clay minerals in Late Glacial and Holocene sediments in northern and southern Aegean Sea. Palaeogeography, Palaeoclimatology, Palaeoecology 249: 36–57.

13.

Francke

Wagner

Leng

et al . (2013) A Late-Glacial to Holocene record of environmental change from Lake Dorjan (Macedonia, Greece). Climate of the Past 9: 481–498.

14.

Ghilardi

Fouache

Queyrel

et al . (2008a) Human occupation and geomorphological evolution of the Thessaloniki plain (Greece) since mid-Holocene. Journal of Archaeological Science 35(1): 111–125.

15.

Ghilardi

Kunesch

Styllas

et al . (2008b) Reconstruction of Mid-Holocene sedimentary environments in the central part of the Thessaloniki Plain (Greece), based on microfaunal identification, magnetic susceptibility and grain-size analyses. Geomorphology 97(3–4): 617–630.

16.

Ghilardi

Psomiadis

Cordier

et al . (2012) The impact of rapid early- to mid-Holocene palaeoenvironmental changes on Neolithic settlement at Nea Nikomedia, Thessaloniki Plain, Greece. Quaternary International 266: 1–15.

17.

Gogou

Bouloubassi

Lykousis

et al . (2007) Organic geochemical evidence of abrupt late glacial-Holocene climate changes in the North Aegean Sea. Palaeogeography, Palaeoclimatology, Palaeoecology 256: 1–20.

18.

Hamann

Ehrmann

Schmiedl

et al . (2008) Sedimentation processes in the eastern Mediterranean Sea during the Late Glacial and Holocene revealed by end-member modelling. Marine Geology 248: 97–114.

19.

Hurrell

and NCAR Staff (eds) (2016) The climate data guide: Hurrell North Atlantic Oscillation (NAO) index (station-based). Available at: https://climatedataguide.ucar.edu/climate-data/hurrell-north-atlantic-oscillation-nao-index-station-based.

20.

Karageorgis

Anagnostou

(2003) Seasonal variation in the distribution of suspended particulate matter in the northwest Aegean Sea. Journal of Geophysical Research 18: 1–8.

21.

Kotthoff

Müller

Pross

et al . (2008a) Lateglacial and Holocene vegetation dynamics in the Aegean region: An integrated view based on pollen data from marine and terrestrial archives. The Holocene 18: 1019–1032.

22.

Kotthoff

Pross

Müller

et al . (2008b) Climate dynamics in the borderlands of the Aegean Sea during formation of Sapropel S1 deduced from a marine pollen record. Quaternary Science Reviews 27: 832–845.

23.

Krahtopoulou

Veropoulidou

(2014) Linking coastal and inland records: Landscape and human histories in Pieria, Macedonia, Greece. In: Touchais

Laffineur

Rougemont

(eds) PHYSIS: L’ environment naturel et la relation homme-milieu dans le monde égéen protohistorique: AEGAEUM, 14e Rencontre Égéne Internationale. Leuven: Peeters Publishers, pp. 153–158.

24.

Kuhnt

Schmiedl

Ehrmann

et al . (2007) Deep-sea ecosystem variability of the Aegean Sea during the past 22 kyr as revealed by benthic foraminifera. Marine Micropaleontology 64: 141–162.

25.

Lawson

Al-Omari

Tzedakis

et al . (2005) Lateglacial and Holocene vegetation history at Nisi Fen and the Boras mountains, northern Greece. The Holocene 15: 873–887.

26.

Leng

Baneschi

Zanchetta

et al . (2010) Late Quaternary palaeoenvironmental reconstruction from Lakes Ohrid and Prespa (Macedonia/Albania border) using stable isotopes. Biogeosciences 7: 3109–3122.

27.

Lespez

Glais

Lopez-Saez

et al . (2016) Middle Holocene rapid environmental changes and human adaptation in Greece. Quaternary Research 85: 227–244.

28.

Lionello

Malanotte-Rizzoli

Boscolo

et al . (2006a) The Mediterranean climate: An overview of the main characteristics and issues, Mediterranean climate variability. Elsevier B.V. 20: 1–26.

29.

Macklin

Tooth

Brewer

et al . (2010) Holocene flooding and river development in a Mediterranean steepland catchment: The Anapodaris Gorge, south central Crete, Greece. Global and Planetary Change 70: 35–52.

30.

Magny

(2004) Holocene climatic variability as reflected by mid-European lake-level fluctuations, and its probable impact on prehistoric human settlements. Quaternary International 113: 65–79.

31.

Magny

Bégeot

Guiot

et al . (2003) Contrasting patterns of hydrological changes in Europe in response to Holocene climate cooling phases. Quaternary Science Reviews 22: 1589–1596.

32.

Marino

Rohling

Sangiorgi

et al . (2009) Early and middle Holocene in the Aegean Sea: Interplay between high and low latitude climate variability. Quaternary Science Reviews 28: 3246–3262.

33.

Mayewski

Rohling

Stager

et al . (2004) Holocene climate variability. Quaternary Research 62: 243–255.

34.

Panagiotopoulos

Aufgebauer

Schäbitz

et al . (2013) Vegetation and climate history of the Lake Prespa region since the Lateglacial. Quaternary International 293: 157–169.

35.

Peyron

Goring

Dormoy

et al . (2011) Holocene seasonality changes in the central Mediterranean region reconstructed from the pollen sequences of Lake Accesa (Italy) and Tenaghi Philippon (Greece). The Holocene 21(1): 1–16.

36.

Pope

Van Andel

(1984) Late Quaternary alluviation and soil formation in southern Argolid: Its history, causes and archaeological implications. Journal of Archaeological Science 11: 281–306.

37.

Pross

Kothoff

Müller

et al . (2009) Massive perturbation in the terrestrial ecosystems of the Eastern Mediterranean region associated with the 8.2 kyr B.P. climatic event. Geology 37: 887–890.

38.

Psomiadis

Ghilardi

Demory

et al . (2014) Late Pleistocene to Mid-Holocene landscape reconstruction in the western part of Thessaloniki Plain, Greece: Evidence for environmental changes, and implications for human occupation. Zeitschrift für Geomorphologie 58: 067–087.

39.

Rohling

Pälike

(2005) Centennial-scale climate cooling with a sudden cold event around 8,200 years ago. Nature 434: 975–979.

40.

Rohling

Mayewski

Abu-Zied

et al . (2002) Holocene atmosphere-ocean interactions: Records from Greenland and the Aegean Sea. Climate Dynamics 18: 587–593.

41.

Roussakis

Karageorgis

Conispoliatis

et al . (2004) Last glacial-Holocene sediment sequences in N. Aegean basins: Structure, accumulation rates and clay mineral distribution. Geo-Marine Letters 24: 97–111.

42.

Saarni

Saarinen

Dulski

(2016) Between the North Atlantic Oscillation and the Siberian High: A 4000-year snow accumulation history inferred from varved lake sediments in Finland. The Holocene 26(3): 423–431.

43.

Schmiedl

Kuhnt

Ehrmann

et al . (2010) Climatic forcing and regional differentiation of eastern Mediterranean deep-water formation and benthic ecosystems during the past 22000 years. Quaternary Science Reviews 29: 3006–3020.

44.

Styllas

(2017) Natural processes versus human impacts during the last century: A case study of the Aliakmon River Delta. In: Barceló

Kostianoy

(eds) Handbook of Environmental Chemistry Special Edition. Berlin: Springer, pp. 1–19.

45.

Tonkov

Bozilova

Possnert

(2013) Postglacial vegetation history as recorded from the subalpine Lake Ribno (NW Rila Mts) Bulgaria. Open Life Sciences 8: 64–77.

46.

Triantaphyllou

Gogou

Bouloubassi

et al . (2013) Evidence for a warm and humid Mid-Holocene episode in the Aegean and the Levantine Seas (Greece, NE Mediterranean). Regional Environmental Change 14: 1697–1712.

47.

Triantaphyllou

Gogou

Dimiza

et al . (2016) Holocene climate optimum centennial-scale paleoceanography in the NE Aegean (Mediterranean Sea). Geo-Marine Letters 36: 51–66.

48.

Tzedakis

Palike

Roucoux

et al . (2009) Atmospheric methane, southern Europe vegetation and low mid-latitude links on orbital and millennial time scales. Earth and Planetary Science Letters 277(3–4): 307–371.

49.

Vogel

Wagner

Zanchetta

et al . (2010a) A paleoclimate record with tephrochronological age control for the last glacial-interglacial cycle from Lake Ohrid, Albania and Macedonia. Journal of Paleolimnology 44: 295–310.

50.

Webb

Orr

(1997) Analytical Methods in Fine Particle Technology. Norcross, GA: Micromeritics Instrument Corporation.

51.

Wolfe

Karst-Riddoch

Vardy

et al . (2005) Impacts of climate and river flooding on the hydro-ecology of a floodplain basin, Peace-Athabasca Delta, Canada since A.D. 1700. Quaternary Research 64: 147–162.

52.

Guiot

Brewer

et al . (2007) Climatic changes in Eurasia and Africa at the last glacial maximum and mid-Holocene: Reconstruction of polled data using inverse vegetation modelling. Climate Dynamics 29(2): 211–229.

53.

Xoplaki

González-Rouco

Luterbacher

et al . (2004) Wet season Mediterranean precipitation variability: Influence of large-scale dynamics and trends. Climate Dynamics 23: 63–78.

54.

Zhang

Reed

Wagner

et al . (2014) Lateglacial and Holocene climate and environmental change in the northeastern Mediterranean region: Diatom evidence from Lake Dojran (Republic of Macedonia/Greece). Quaternary Science Reviews 103:51–66.