Multiscale monsoon variability during the last two climatic cycles revealed by spectral signals in Chinese loess and speleothem records

The East Asian Monsoon (EAM) exhibits a significant variability on timescales ranging from tectonic to centennial as inferred from loess, speleothem and marine records. However, the relative contributions and plausible driving forces of the monsoon variability at different timescales remain controversial. Here, we spectrally explore time series of loess grain size and speleothem δO records and decompose the two proxies into intrinsic components using the empirical mode decomposition method. Spectral results of these two proxies display clear glacial and orbital periodicities corresponding to ice volume and solar cycles, and evident millennial signals which are in pace with Heinrich rhythm and Dansgaard–Oeschger (DO) cycles. Five intrinsic components are parsed out from loess grain size and six intrinsic components from speleothem δO records. Combined signals are correlated further with possible driving factors including the ice volume, insolation and North Atlantic cooling from a linear point of view. The relative contributions of components differ significantly between loess grain size and speleothem δO records. Coexistence of glacial and orbital components in the loess grain size implies that both ice volume and insolation have distinctive impacts on the winter monsoon variability, in contrast to the predominant precessional impact on the speleothem δO variability. Moreover, the millennial components are evident in loess grain size and speleothem δO records with variances of 13 and 17 %, respectively. A comparison of the millennial-scale signals of these two proxies reveals that abrupt changes in the winter and summer monsoons over the last 260 kyr share common features and similar driving forces linked to highlatitude Northern Hemisphere climate.


Introduction
The East Asian Monsoon (EAM), as a significant part of Asian monsoon circulation, plays an important role in driving the palaeoenvironmental changes in East Asia (An, 2000).The EAM fluctuations can be quantified at different time intervals ranging from thousands of years to intraseasonal periodicities, and the primary driving force of the monsoon variability on each timescale is not unique (An et al., 2015).Multiscale monsoon variability has been inferred from numerous proxies generated from deep-sea sediments (e.g.Wang et al., 1999Wang et al., , 2005)), eolian deposits (e.g.An, 2000;Sun et al., 2012), and speleothem records (e.g.Wang et al., 2001Wang et al., , 2008)), which provide valuable insights into the changing processes and potential driving forces of the EAM variability.In particular, Chinese loess has been investigated intensively as a direct and complete preserver of the EAM changes, with great effort spent on deciphering the EAM variability at both orbital and millennial scales (e.g.An et al., 1990;Ding et al., 1994Ding et al., , 2002;;Porter and An, 1995;Guo et al., 1996;Chen et al., 1997;Liu and Ding, 1998;Liu et al., 1999;An, 2000;Chen et al., 2006).
On the orbital timescale, the EAM variation recorded by Chinese loess-paleosol sequences was characterized by an alternation between the dry-cold winter monsoon and the wet-warm summer monsoon (Liu and Ding, 1998;An, 2000).A strong 100 kyr periodicity was detected in the Chinese loess particle size record, implying an important impact of glacial boundary conditions on the EAM evolution (Ding et al., 1995).Obliquity and precession signals were also clear in loess-based proxies (Liu et al., 1999;Ding et al., 2002;Sun et al., 2006).In addition to these dominant periodicities, some harmonic periodicities related to orbital parameters were also found in the EAM records, such as the ∼ 75, ∼ 55, and ∼ 30 ka spectral peaks (Lu et al., 2003;Sun et al., 2006;Yang et al., 2011).In contrast, absolute-dated speleothem δ 18 O records revealed an evident 23 kyr cycle, implying a dominant role of summer insolation in driving the summer monsoon variability (Wang et al., 2008;Cheng et al., 2009).Different variances of obliquity and precession signals in monsoonal proxies suggest that the responses of the winter and summer monsoons to the orbital forcing were dissimilar (Shi et al., 2011).The various patterns of orbitalscale monsoon fluctuations between the loess proxies and speleothem δ 18 O records most likely reflected the sensitivity of various archives and proxies to the EAM variability (Clemens et al., 2010;Cheng et al., 2012;Sun et al., 2015;Cai et al., 2015).
On the millennial timescale, the rapid monsoon oscillations inferred from Chinese loess were not only persistent during the last two glacial cycles (Porter and An, 1995;Guo et al., 1996;An and Porter, 1997;Chen et al., 1997;Ding et al., 1999;Sun et al., 2010;Yang and Ding, 2014), they were also evident during early glacial extreme climatic conditions (Lu et al., 1999).The millennial-scale monsoon variability during the last glacial period was strongly coupled to climate changes recorded in the Greenland ice core and North Atlantic sediments, indicating a dynamic connection between the EAM variability and the high-latitude Northern Hemisphere climate (Porter and An, 1995;Guo et al., 1996;Chen et al., 1997;Fang et al., 1999).Recently, a combination of proxies from Chinese loess, speleothem, and Greenland ice core with modelling results indicated that the Atlantic meridional overturning circulation might have played an important role in driving the rapid monsoon changes in East Asia during the last glaciation (Sun et al., 2012).
Although previous studies have revealed that past EAM variability principally comprise a mixture of forcing signals from ice volume, solar radiation, and North Atlantic climate, the relative contributions of glacial, orbital and millennial forcing to the EAM variability remain unclear.In this study, we conducted a comprehensive investigation of multiscale EAM variability over the last 260 kyr by analysing the mean grain size (MGS) record from a Gulang loess sequence (a proxy indicator of the East Asian winter monsoon intensity) and speleothem δ 18 O record of Hulu and Sanbao caves (a debatable indicator of the summer monsoon intensity).Our objectives are to evaluate the relative contributions of glacial-interglacial-to millennial-scale signals registered in these two widely used monsoon proxies, and to emphasise the glacial-interglacial discrepancy and millennial similarity between loess and speleothem records during the last two glacial cycles.

Data and methods
The data for the loess sequence was collected at a section in Gulang, Gansu Province, China (37.49• N, 102.88 • E, 2400 m a.s.l.), which is situated in the northwestern part of the Chinese Loess Plateau.It is about 10 km to the southwest, on the margin of the Tengger Desert (Fig. 1).In this region, the average annual precipitation and temperature over the last 20 years are 350 mm and 5.7 • C.About 70 m of loess accumulated at Gulang during the last two climate cycles.A high sedimentation rate and weak pedogenesis in this region made the Gulang loess sequence very sensitive to orbital and millennial monsoon changes (Sun et al., 2012(Sun et al., , 2015)).
The samples used in this study were collected at 2 cm intervals, corresponding to 50-100-year resolution for the loess-paleosol sequence.The grain size data of the upper 20 m were from a 20 m pit near Gulang (Sun et al., 2012), and the lower part spanning the last two glacial cycles were from another 50 m section (for the full grain size data, see Supplement Table S1).The mean grain size data of the composite 70 m section have been employed for a chronological reconstruction (for a detailed description, see Sun et al., 2015).The Gulang chronology was evaluated by comparison with a 249 kyr Chinese loess millennial-scale oscillation stack (CHILOMOS) record in the northern Loess Plateau (Yang and Ding, 2014) (Fig. 2); the close matches between these two records indicate a high reliability of our Gulang age construction.Unlike previous studies (Sun et al., 2012(Sun et al., , 2015)), here we performed spectral and decomposing analy-  Sun et al., 2015) and CHILOMOS stack median grain size (Md, green; Yang and Ding, 2014) with the benthic δ 18 O (black; Lisiecki and Raymo, 2005) and Sanbao/Hulu speleothem δ 18 O (magenta; Wang et al., 2008;Cheng et al., 2009) records.The red and black dashed lines denote tie points derived from optically stimulated luminescence (OSL) dating and benthic δ 18 O correlation, respectively.sis on the mean grain size time series in order to decipher multiscale variability and dynamics of the winter monsoon.
The absolute-dated speleothem δ 18 O records from Sanbao/Hulu caves (0-224 ka, Wang et al., 2008) and the Sanbao cave (224-260 ka, Cheng et al., 2009) (Fig. 1; Supplement Table S2) were selected to infer the summer monsoon variability spanning the last two glacial-interglacial cycles.Compatible with the analysis by Wang et al. (2008), we plotted the Hulu δ 18 O data 1.6 ‰ more negative than that from the Sanbao cave (Fig. 2).Interpretation of the Chinese speleothem δ 18 O records remains debatable as a direct indicator of the summer monsoon intensity since various factors like seasonal changes in precipitation amount, moisture sources, and circulation patterns would influence the speleothem δ 18 O composition (e.g.Yuan et al., 2004;Wang et al., 2001Wang et al., , 2008;;Cheng et al., 2009;Clemens et al., 2010;Dayem et al., 2010;Pausata et al., 2011;Maher and Thompson, 2012;Caley et al., 2014).Nevertheless, the high similarity between millennial events in Chinese speleothem and Greenland ice core revealed that speleothem δ 18 O is a reliable indicator of seasonal monsoon change (Wang et al., 2001;Clemens et al., 2010).More recently, a modeldata comparison suggested that Chinese speleothem δ 18 O can be regarded as a monsoon proxy to reflect the southerly wind intensity rather than the precipitation change (Liu et al., 2014).Thus, spectral and decomposed results of the composite speleothem δ 18 O record time series were used in this study to address multiscale variability and dynamics of the summer monsoon.
To detect the presence of glacial-millennial periodicities, we performed spectral analysis on the 260 ka records of Gulang MGS and speleothem δ 18 O using both the Multi-taper (MTM, implemented in the SSA toolkit, Vautard et al., 1992) (http://www.atmos.ucla.edu/tcd/ssa/)and REDFIT (Schulz and Mudelsee, 2002) methods, which are related to the empirical orthogonal function and Lomb-Scargle Fourier transform, respectively.The MTM method has the advantages of quantified and optimized trade-off between the spectral leakage reduction and variance reduction, and is suitable for time series affected by high-noise levels (Lu et al., 1999), but MTM requires equally spaced data and therefore an interpolation is needed.The REDFIT program estimates the firstorder autoregressive (AR1) parameter from unevenly sampled time series without interpolation, which avoids a too "red" spectrum (Schulz and Stattegger, 1997), but uses the weighted overlapped segment averaging (WOSA) method for the spectral leakage reduction and variance reduction, which makes the trade-off not quantifiable.The similar spectral periodicities derived from both REDFIT and MTM methods were regarded as dominant frequencies at glacialmillennial bands.
The decomposed components of loess MGS and speleothem δ 18 O records were parsed out using the technique of empirical mode decomposition (EMD; Huang et al., 1998) (for the MATLAB code of EMD method, see the code in the Supplement).EMD directly extracts energy which is associated with intrinsic timescales in nonlinear fluctuations, and iteratively decomposes the raw complex signal with several characteristic timescales coexisting into a series of elementary intrinsic model function (IMF) components, avoiding any arbitrariness in the choices of frequency bands in this multiscale study.The EMD method has been widely applied to various palaeoclimate databases such as ice cover (Gloersen and Huang, 2003), North Atlantic oscillation (Hu and Wu, 2004), solar insolation (Lin and Wang, 2006), and temperature under global warming (Molla et al., 2006).This approach has also been used to decipher the multiscale variations of Indian monsoon (Cai et al., 2015).However, the application of the EMD method on the loess record remains poorly investigated with limited understanding of decomposed components at glacial and orbital timescales due to the low-resolution proxy variations (Yang et al., 2011(Yang et al., , 2008)).In this study, we applied EMD on linearly interpolated loess and speleothem data with a 100-year interval to quantify the relative contributions of both orbital and millennial components.

Multiscale monsoon variability
The highly comparable spectral results between REDFIT and MTM methods showed that apparent periods identified in the MGS spectrum are at ∼ 100, ∼ 41, ∼ 23, ∼ 15, ∼ 7, ∼ 5, ∼ 4, and ∼ 3-1 kyr, for REDFIT and MTM methods, over the 80 and 90 % confidence levels, respectively (Fig. 3).It is shown that the potential forcing of the glacial-interglacial and orbital EAM variability is part of external (e.g. the orbitalinduced summer insolation; An, 1991;Wang et al., 2008)  and the internal factors (e.g. the changes in the ice volume and CO 2 concentrations; Ding et al., 1995;Lu et al., 2013;Sun et al., 2015).The coexistence of the ∼ 100, ∼ 41, and ∼ 23 kyr periods in the Gulang MGS record confirms the dynamic linkage of the winter monsoon variability to glacial and orbital forcing.Based on the spectral results, many millennial frequencies are detected, which can be mainly divided into two groups, namely, ∼ 7-4 and ∼ 3-1 kyr, which possibly correspond to the Heinrich (∼ 6 kyr) rhythm and the DO (∼ 1.5 kyr) cycles recorded in the North Atlantic sediments and Greenland ice core (Bond et al., 1993;Dansgaard et al., 1993;Heinrich, 1988).Taking into account the sampling resolution and surface mixing effect at Gulang, the residual component (< 1 kyr) might contain both centennial and noisy signals, which was excluded for further discussion in this study.
The different oscillation patterns composing loess MGS and speleothem δ 18 O time series are separated out using the EMD method as presented in Figs. 4 and 5. REDFIT spectral analysis is further carried out on each IMF with dominant periods as shown.Five IMFs are generated for the Gulang MGS data on the glacial-to-millennial timescale.The variability of Gulang MGS is dominated by the lowest fre-  (Berger, 1978) to ascertain plausible impacts of glacial and orbital factors on the EAM variability (Fig. 6).The low-frequency component of the Gulang MGS record is well correlated to the global ice volume change inferred from the benthic δ 18 O record (Lisiecki and Raymo, 2005) with a correlation coefficient (R 2 ) of 0.56, reinforcing the strong coupling between the winter monsoon variation and ice volume changes, particularly in terms of glacialinterglacial contrast (Ding et al., 1995).However, fine MGS signals at the precessional scale seem more distinctive than those in the benthic δ 18 O stack.For example, the remarkable peaks in the MGS around 85, 110, and 170 ka have no counterpoints in the benthic δ 18 O record.By comparing MGS data with the summer insolation record, the overall ∼ 20 kyr periodicity is damped but still visible during both glacial and interglacial periods except for the insolation maxima around 150 and 220 ka (Fig. 6).The coexistence of the glacial and orbital cycles in loess MGS indicates that both the ice volume and solar insolation have affected the winter monsoon variability, and their relative contributions are 32 and 55 %, respectively, as estimated from variances of the glacial (IMF5) and orbital (IMF4 and IMF3) components.Wang et al., 2008;Cheng et al., 2009) records with summer insolation at 65 • N (red; Berger, 1978) and benthic δ 18 O record (black; Lisiecki and Raymo, 2005).The vertical gray bars represent the interglacial periods.
The speleothem δ 18 O record varies quite synchronously with the July insolation, characterized by a dominant precession frequency (Fig. 6).This in-phase change is thought to support a dominant role of summer insolation in the Northern Hemisphere in driving the summer monsoon variability at the precession period (Wang et al., 2008), given that the palaeoclimatic interpretation of the speleothem δ 18 O is quite controversial (Wang et al., 2001(Wang et al., , 2008;;Yuan et al., 2004;Hu et al., 2008;Cheng et al., 2009;Peterse et al., 2011).
The different contributions of glacial and orbital variability in the loess MGS and speleothem δ 18 O records indicate that the driving forces associated with these two proxies are different.The loess grain size is directly related to the northwesterly wind intensity, reflecting that atmospheric process is linked to the Siberian-Mongolian High (Porter and An, 1995).The speleothem δ 18 O might be influenced by multiple factors such as the isotopic depletion along the vapour transport path (Pausata et al., 2011), changes in δ 18 O values of meteoric precipitation or the amount of summer monsoon precipitation (Wang et al., 2001(Wang et al., , 2008;;Cheng et al., 2009), and seasonality in the amount and isotopic composition of rainfall (Clemens et al., 2010;Dayem et al., 2010;Maher and Thompson, 2012).Even at the orbital timescale, proxymodel comparison suggests that the response of the winter and summer monsoon to obliquity and precession forcing are dissimilar (Shi et al., 2011) It is quite clear that the EAM is formed by the thermal gradient between the Asian continent and the Pacific Ocean to the east and southeast (Halley, 1986;Xiao et al., 1995;Lestari and Iwasaki, 2006).In winter, due to a much larger heat capacity of water in the ocean than that on the land surface, a higher barometric pressure forms over the colder Asian continent with a lower pressure over the warmer ocean.This gradient is the driving force for the flow of cold and dry air out of Asia, forming the winter monsoon (Gao, 1962).On the glacial-interglacial timescale, the buildup of the northern high-latitude ice sheets during the glacial periods strengthens the barometric gradient which results in intense winter monsoons (Ding et al., 1995;Clark et al., 1999).The contemporaneous falling sea level and land-ocean pressure gradient further enhances winter monsoon circulation during glacial times (Xiao et al., 1995).The other factor that influences the land-ocean differential thermal motion is the orbitally induced changes in solar radiation.The precession-induced insolation changes can lead to regional land-ocean thermal gradients whilst obliquity-related insolation changes can result in meridional thermal gradients, both of which can substantially alter the evolution of the Siberian and Subtropical Highs and the EAM variations (Shi et al., 2011).

Impacts of high-latitude cooling on millennial EAM oscillations
The EAM variations are persistently punctuated by apparent millennial-scale monsoon events (Garidel-Thoron et al., 2001;Wang et al., 2001;Kelly et al., 2006).The millennialscale events of the last glacial cycle were firstly identified in Greenland ice cores (Dansgaard et al., 1993;Meese et al., 1997).Subsequently, well-dated loess grain size and speleothem δ 18 O records in China have been found to apparently correspond with rapid climate oscillations in the North Atlantic (Porter and An, 1995;Guo et al., 1996;Chen et al., 1997;Ding et al., 1998;Wang et al., 2001).The most striking evidence is the strong correlation between the loess grain size, speleothem δ 18 O and Greenland ice core δ 18 O records during the last glaciation (Ding et al., 1998;Wang et al., 2001;Sun et al., 2012).The findings of these abrupt changes have been extended to investigate glacial-interglacial cycles using loess and speleothem records (Ding et al., 1999;Cheng et al., 2006Cheng et al., , 2009;;Wang et al., 2008;Yang and Ding, 2014) and North Atlantic sediments (McManus et al., 1999;Channell et al., 2012).Unlike previous comparisons based on original proxy variability, here we combine the IMF1 and IMF2 components of the loess MGS and IMF1, IMF2, and IMF3 components of speleothem δ 18 O records as robust reflections of millennialscale signals of the winter and summer monsoons, with variances of 13 and 17 %, respectively.The combination of the two millennial signals of the loess MGS and speleothem δ 18 O records are compared further with the North Atlantic cooling events over the last two glacial cycles to reveal the dynamic links of abrupt climate changes in East Asia and the North Atlantic (Fig. 7).The Younger Dryas (YD) and Hein- rich events (H 1 -H 6 ) are well detected in loess and speleothem records around 12, 16, 24, 31, 39, 48, 55, and 60 ka.Most of the millennial-scale events in the loess MGS and speleothem δ 18 O records are well aligned with comparable timing and duration during the last two glacial cycles.However, some MGS valleys such as A17, A23, B17, B18, and B22 are not well matched with the speleothem δ 18 O minima, possibly due to uncertainties in the loess chronology.The comparable millennial-scale events between grain size of Gulang and the CHILOMOS stack (Yang and Ding, 2014) show the nature of replication of Gulang MGS record within the dating uncertainty, confirming the persistent millennial-scale winter monsoon variability spanning the last two glacial cycles (Fig. 7).
The millennial-scale monsoon signals over the last two glacial cycles have been well compared with the cooling events recorded in the North Atlantic sediments, demonstrating a dynamic link between abrupt climate changes in East Asia and the North Atlantic.As identified in Chinese speleothem records, the magnitudes of abrupt climate events are identical between the last and the penultimate glacial cycles (Wang et al., 2008).However, the duration and amplitude of these millennial events are quite different between the glacials and interglacials.The duration of the millennial monsoon events is relatively shorter and the amplitude larger during glacial periods, suggesting a plausible glacial modulation of rapid climate changes (McManus et al., 1999;Wang et al., 2008).The potential driving mechanism for rapid EAM changes has been attributed to changing climate in the highlatitude Northern Hemisphere, e.g. the reduction of the North Atlantic deep water circulation triggered by fresh water inputs from melting icebergs (Broecker, 1994).The North Atlantic cooling can affect the zonal high pressure systems, including the Azores-Ural-Siberian-Mongolian high (Palmer and Sun, 1985;Rodwell et al., 1999;Yuan et al., 2004), which can further transmit the abrupt cooling effect into East Asia and result in significant EAM changes (Porter and An, 1995;Wang et al., 2001).Apart from the geological evidence, numerical modelling also suggests that the Atlantic meridional overturning circulation might affect abrupt oscillations of the EAM, while the westerly jet is the important conveyor introducing the North Atlantic signal into the EAM region (Miao et al., 2004;Zhang and Delworth, 2005;Jin et al., 2007;Sun et al., 2012).

Conclusions
The EAM displays variabilities on timescales ranging from thousands of years to intraseasonal periodicities as inferred from numerous proxy indicators.Multiscale signals were spectrally detected and naturally decomposed from Chinese loess and speleothem records over the last two climatic cycles in this study, permitting an evaluation of the relative contributions of glacial, orbital and millennial components in the EAM record from a linear point of view.The spectra of Gulang MGS and speleothem δ 18 O data show similar periodicities at glacial-to-orbital and millennial timescales, corresponding to the rhythms of changing ice volume, orbitally induced insolation, and North Atlantic cooling (i.e.Heinrich rhythm and Dansgaard-Oeschger cycles).Amplitude variances of the decomposed components reveal significant glacial and orbital impacts on the variation in loess grain size and the dominant precession forcing in the speleothem δ 18 O variability.The millennial components are evident in the loess and speleothem proxies with variances of 13 and 17 %, respectively.Millennial IMFs were combined to decode the synchronous nature of rapid changes of these two proxies.High similarity of millennial-scale monsoon events both in terms of the rhythms and duration between the loess and speleothem proxies implies that the winter and summer monsoons share common millennial features and similar driving forces.
The Supplement related to this article is available online at doi:10.5194/cp-11-1067-2015-supplement.

Figure 1 .
Figure 1.Map showing the loess distribution and locations of the Gulang loess section, Sanbao, and Hulu caves.Dotted lines indicate the precipitation isohyets (PI).

Figure 4 .Figure 5 .
Figure 4. IMFs of Gulang MGS series (a) and corresponding spectra (b).Numbers in black are dominant periods and dotted lines represent the 90 % confidence level.

Figure 7 .
Figure 7.Comparison of millennial-scale variations among Gulang MGS (blue), CHILOMOS stack Md (green; Yang and Ding, 2014) and Sanbao/Hulu speleothem δ 18 O (magenta; Wang et al., 2008; Cheng et al., 2009) records over the last two glacial-interglacial cycles.Cyan dotted lines are the YD and the Heinrich events identified among the three records and gray bars indicate interglacial periods.The numbers represent well-correlated Chinese interstadials identified among the three records.