Quaternary rainfall variability is governed by insolation in northern China and ice-sheet forcing in the South

Quaternary Asian low-latitude hydroclimate cyclicity has long been attributed to insolation forcing, in contrast to the dominant ice-sheet and CO2 controls identified in mid-high-latitude regions. However, debates exist regarding these rainfall variations and forcings due to inconsistent reconstructions and simulations. Here, by combining rainfall proxy records with multi-model simulations, dominant 23 ka rainfall cycle in northern China and 100 ka rainfall cycle in southern China and Southeast Asia were found. We propose that rainfall mainly occurs in summer in the north, primarily driven by insolation. Rainfall in the south is largely forced by high-latitude ice sheets, with enhanced spring and autumn rainfall in southern China and weakened rainfall in western Maritime Continent during glacial periods. This study highlights the seasonal contributions to orbital-scale rainfall changes, and sheds light on the Asian hydroclimate conditions associated with high-low-latitude climate interactions. The latitude-dependent periodicity of East Asian hydroclimate during the Quaternary was probably a response to summer insolation in the north and ice sheet forcing in the south, suggested by analyses of proxy-based rainfall records and climate model simulations.

M onsoons are seasonal variations in atmospheric circulations accompanied by the reversal of dry and wet surface climate conditions, and are mainly forced by seasonal pressure system shifts driven by the annual/seasonal solar radiative forcing cycle 1 . Theoretically, precession forcing contributes more to the solar radiation amplitude received by Earth than obliquity or eccentricity 2 , thus, the 23 ka cycle should be considered as the dominant mode affecting monsoon climate variations. This hypothesis has been tested in some early rainfall reconstruction and modeling studies in mid-low-latitude regions [3][4][5][6][7] . In contrast, highlatitude hydroclimate conditions tend to respond to ice-sheet and CO 2 forcings as a result of multiple atmospheric processes, thus suggesting strong glacial-interglacial cycles of 100 ka in these regions [8][9][10] .
Such latitude dependent forcing of monsoon hydroclimates has also been proposed in the East Asian monsoon region. Monsoon rainfall variabilities in the Loess Plateau and Yangtze River Valley in northern-central China are more sensitive to ice-sheet and CO 2 forcings and less sensitive to insolation forcing than that in southern China 9,[11][12][13] . This results in dominant 100 ka cycle of summer monsoon climate in northern China but 23 ka cycle in southern China. However, all recent precipitation isotope record reanalysis and model simulations have revealed that glacialinterglacial cycles are dominant in southern China 14,15 , even in tropical Southeast Asia [15][16][17] , highlighting the dominant ice sheet and CO 2 forcing on the lower latitude hydroclimate. Another paradox is that most previous simulation studies have suggested the existence of strong precession forcing on the East Asian monsoon hydroclimate [18][19][20][21][22][23] . Even when ice-sheet and greenhouse gas forcings were included in these simulations, the dominant precession signal can still be observed in northern China 14,24 . This can be supported by recently rainfall reconstructions in the Tengger Desert in northern China, which revealed prominent precession cycles in local rainfall variation, challenging the view that highlatitude ice sheets are the primary driver of northern China monsoon hydroclimate 25 . Such controversial complicates our understanding of hydroclimate in East Asian monsoon domain and surrounding tropical climate zone associated with high-lowlatitude climate interactions.
To address these conundrums, we present a silicate weathering record from Site U1429 (31°37.04'N, 128°59.85'E) of the International Ocean Discovery Program (IODP) in the northern East China Sea ( Fig. 1a and Supplementary Fig. 1a) to quantitatively reconstruct the rainfall changes in northern China over the last 400 ka. We also analyzed a transient climate model simulation (comprising orbital (ORB) + greenhouse (GHG) + ice sheet (ICE) forcings) through the last 300 ka using the fully coupled National Center for Atmospheric Research-Community Climate System Model version 3 (NCAR-CCSM3) model 26,27 , and studied the periodicity of simulated rainfall variations and the drivers of this periodicity in East-Southeast Asia. In addition, we compared the Transient Climate Evolution (TraCE)-21ka simulation (and the individual forcing simulations) performed using the same model 28,29 and iTraCE based on an isotope-enabled Community Earth System Model (iCESM) [30][31][32] to assess the hydroclimatic changes across East-Southeast Asia, focusing on the individual forcing mechanisms and isotope-climate relationship, respectively.

Results and discussion
Dominant precession cyclicity of rainfall variability in northern China. It is feasible to reconstruct the rainfall series in northern China based on silicate weathering record derived from IODP Site U1429, due to the robust sediment provenance tracing results based on Nd isotopic and rare earth element compositions; these results suggest that the Yellow River in northern China is the major supplier of clay-sized sediments in the study region (Supplementary Figs. 1a and 2) (Supplementary Discussion). The Yellow River has long been regarded as a typical river influenced by intense catchment erosion. It is globally known for its high sediment flux (~11 billion ton/yr, before 1980s) sourced dominantly from the Chinese Loess Plateau (~90% of sediments) in the middle reach of this river 33,34 . Such huge loess materials supply makes the effect from eolian dust on Yellow River sediment composition can be ignored. Thus chemical alteration of sediments during transport from Loess Plateau to our study site has the potential to record local hydroclimate changes.
The chemical weathering intensity can be calculated by comparing the concentration of mobilized elements (K) to a nonmobile element (Al) during incongruent weathering of silicates 35 . Lower K-to-Al ratios (K/Al) indicate greater mobile element depletion in sediments and thus higher chemical alteration. The Yellow River basin has an arid to semi-arid continental climate, being more arid and cold in the upper and middle reaches, and more humid and temperate in the lower reaches. Here we selected surface sediment samples distribute from upper to lower reaches of Yellow River ( Supplementary Fig. 1a). High linear correlations have been observed between the K/Al ratios of clay-sized sediments and mean annual temperature and rainfall at each site, suggesting the controlling effects of temperature and rainfall on the silicate weathering intensity ( Supplementary Fig. 1b, c). Sediment samples of Yellow River and Site U1429 were all analyzed on the clay-sized fraction, which avoids the grain-size effect on their chemical composition. Effect of physical erosion on chemical weathering has also been excluded (Supplementary Discussion). Here, we established a model covering K/Al ratio, temperature and rainfall based on the multiple statistical analyses, and further quantitatively calculated the rainfall changes during the last 400 ka (Supplementary Discussion). The calculated mean annual rainfall ranges from 470-717 mm (Fig. 2c). This is close to the modern annual rainfall (~562-648 mm) in the Yellow River middle to lower reaches 36 and our simulated annual rainfall amount (~407-854 mm) in northern China (Fig. 2c, d).
Our rainfall reconstruction result suggests a distinct story compared to the findings of most previous studies. It shows strong precession cycle and is consistent with the 300 ka simulation results of annual rainfall in northern China (Fig. 2c, d). This differs from the loess magnetic susceptibility, a classic summer monsoon rainfall proxy, and from existing records of loess carbon and beryllium isotopes ( Supplementary Fig. 3e, f, g), which are all suggest dominant 100 ka cycle 9,37-39 . The loess magnetic susceptibility has been suggested to be strongly affected by the dilution effect of the dust sedimentation rate 23,40 . Loess carbon isotopes in inorganic carbonates (δ 13 C IC ) have been adopted as a terrestrial vegetation proxy, and further the monsoon hydroclimate indicators in northern China, are controlled by the effects of both temperature and rainfall conditions 9,39 . The beryllium isotope ( 10 Be) concerntration in loess can probably also be influenced by sediment grain size 41 despite the dust flux, recycled 10 Be and geomagnetic effects have been removed 37 . Thus, such complex responses in proxies of loess magnetic susceptibilities, carbon and beryllium isotopes likely diminish the 23 ka cycle recorded in rainfall changes in northern China. The spectral analysis of our rainfall record also suggests stronger power at the precession band than at the eccentricity band (Fig. 2c, d), corresponding to the simulated power ratio distributions of the 23 ka to 100 ka annual rainfall cycles (Fig. 1a) and indicating the dominance of precession cyclicity in affecting the hydroclimate of northern China. Such dominant 23 ka cyclical changes in northern China rainfall have also been observed in previous rainfall records in the Yellow River mouth 42 and Tengger Desert 25 ( Supplementary Fig. 3h, i).
Dominant glacial-interglacial cyclicity of rainfall variability in southern China. Explicit rainfall records cover glacial-interglacial cycle in southern China are currently lacking. Precipitation isotopes (e.g., stalagmite δ 18 O and leaf wax δ 2 H) are informative hydrologic tracers that are commonly used to reconstruct hydroclimate conditions in Asian monsoon regions 12,16,[43][44][45] , although the complicated interpretations of these proxies depend on the rainfall amount, moisture transport and monsoon wind intensity changes [46][47][48][49] .
The leaf wax δ 2 H record at ODP Site 1146 in the northern South China Sea exhibits a dominant 100 ka cycle and relatively light isotopic values during glacial periods (Fig. 2e) 12 . Such depleted glacial δ 2 H values have also been observed in a previously conducted water isotope-enabled climate simulation 16 . Some studies proposed that changes in the ratio of local Pacific Ocean-to remote Indian Ocean-sourced moisture that reaches eastern China could probably regulate these rainfall isotopes 46,50 . However, during glacial periods, expansive ice sheets strengthened the pressure gradient between the East Asian continent and the Pacific Ocean, as well as enhance the Western North Pacific Subtropical High, and thus increased the amount of moisture Positive values indicate that the 23 ka cycle of annual rainfall variabilities (natural log ratio, unitless) was stronger than the 100 ka cycle. The star denotes IODP Site U1429 (this study); the northern China records have been shown in red, including Chinese loess sections collected in Gulang and Jingyuan 9,39 and 10 Be record 37 , sediment borehole of Chahanchi Lake in the Tengger Desert 25 and coastal sediment core Lz908; 42 the southern China records including ODP Site 1146 12 and 251PC 53 have been shown in blue; the southeast Asian record of MD98-2152 17 has been shown in purple. Note that the half red and half blue in symbol of IODP Site U1429 has been used to indicate both northern China rainfall based on silicate weathering proxy (this study) and southern China rainfall based on seawater δ 18 O 11 at this site. Besides, the grey triangles represent the locations of stalagmite records in China (Hulu, Sanbao and Dongge) 44   A previous interpretation also linked the depleted glacial δ 2 H at ODP Site 1146 to the increased rainout along the long moisture transport path over the emergent continental shelf during the glacial periods 11 . A major challenge to this perspective is other low-latitude precipitation isotope records including Borneo stalagmite δ 18 O records 43,51 and southern Sumatra δ 2 H record 16 in the western Maritime Continent with broad shelf regions all suggest strong glacial-interglacial cycles, but characterized by the heavier isotopes during glacial periods. These contradictions indicate that the shelf exposure effect on longer moisture transport path was probably not the main driver of depleted rainfall isotopes during glacial periods in southern China.
Hence, this leaf wax δ 2 H record can be largely interpreted as a rainfall proxy in southern China, rather than the results of moisture transport and shelf exposure effect. This can be supported by the iTraCE simulation results, which suggest positive annual, June-July-August (JJA), March-April-May (MAM) and September-October-November (SON) rainfall anomalies in most regions of southern China, accompanied by negative δ 18 O anomalies at the last glacial maximum (LGM) relative to 10 ka ( Supplementary Fig. 4).
The Yangtze River Valley rainfall and runoff proxy of seawater δ 18 O record at IODP Site U1429 is also dominated by glacialinterglacial cycles (Fig. 2e). This record has been reconstructed based on the high correlation between modern rainfall in Yangtze  44 . c Rainfall K/Al record at IODP Site U1429 indicates the rainfall in northern China (this study). d Simulated rainfall in northern China (this study, box NC in Fig. 1b). e The δ 18 O sw record obtained at IODP Site U1429 11 , the rainfall δ 2 H record at ODP Site 1146 12 (marked by circles), and the SSS record of core 251PC in the eastern South China Sea 53 (marked by rectangles). All of them indicate the rainfall in southern China and the South China Sea. f Simulated rainfall in southern China (this study, box SC in Fig. 1b). g δ 13 C wax record MD98-2152 from southern Sumatra 17 indicate the rainfall in the western regions of Southeast Asia. h Simulated Southeast Asian rainfall (this study, box SEA in Fig. 1b). i LR04 stacked benthic δ 18 O curve 72 . Spectral analysis has been performed with PAST software; 88 the window function is rectangle. The gray vertical bars mark the timing of glacial terminations and orbital periodicities, as shown in the left and right panels, respectively. River Valley and salinity at northeastern East China Sea 11,52 . Evaporation effects on the surface salinity in northeastern East China Sea is minimal because precipitation plus runoff dominates over evaporation 11 . This Yangtze River Valley rainfall record shows an identical trend as the leaf wax δ 2 H record at ODP Site 1146 and a previously published sea surface salinity (SSS) record in the South China Sea 53 (Fig. 2e). The high consistency among these records is not a coincidence but most likely implies their response to local rainfall variations.
All of these records suggest dominant glacial-interglacial cycles with relatively wet conditions in southern China and the South China Sea during glacial periods. They show large consistencies with our modeled annual rainfall results in southern China, except for some differences in the fluctuation amplitudes at~40-60 ka and 150-180 ka and some precession scale deviations (Fig. 2f). Such discrepancies probably result from the differences between the domain of the model box, and the region that can be covered by these rainfall proxies, as well as from the age bias between modeling and rainfall reconstructions. Besides, relatively wet glacial periods in southern China have also been observed in reconstructions based on pollen databases [54][55][56] and magnetic investigations 57 ( Supplementary  Fig. 5), as well as in LGM topsoil moisture and net precipitation simulations based on fully coupled climate models 58,59 .
Dominant glacial-interglacial cyclicity of rainfall variability in western Maritime Continent. Critical understanding of hydroclimate indication in rainfall isotopes of Southeast Asia remains lacking. Paleoclimate studies have mainly associated changes of rainfall isotopes in this region with rainfall amount 60 , moisture transport 61,62 , moisture convergence 16 and global ocean δ 18 O variance with glaciation 43 . In our iTraCE results, annual and JJA simulation between 20 ka and 10 ka suggest positive δ 18 O loading corresponding to decreased rainfall in Sumatra and southern Borneo; however, most regions of Sumatra and Borneo suggest decreased spring and autumn rainfall corresponding to negative δ 18 O (Supplementary Fig. 4). These complex seasonal patterns suggest large uncertainty in the interpretation of rainfall isotopes as rainfall amount in this region.
Here, a leaf wax δ 13 C record from MD98-2152 in southern Sumatra offshore has been used to indicate past rainfall changes of the western Maritime Continent (Fig. 2g). Variations in this leaf wax δ 13 C record mainly depend on the transition between C 3 and C 4 plant in Sumatra 17 . C 4 plants are more tolerant of high temperatures and arid climates than C 3 plants, so they are often utilized to reconstruct the local rainfall seasonality 17,63 . Generally, this leaf wax δ 13 C record has enriched δ 13 C values during glacial relative to the interglacial periods, suggesting that the 100 ka cycle is the dominant frequency. Enriched δ 13 C values during glacial periods indicate more C 4 plants and higher rainfall seasonality (enhanced dry season) in Sumatra 17 . This is consistent with our annual rainfall simulation results in this region (Fig. 2h), which suggests a relatively drier western Maritime Continent during glacial periods.
Forcing of northern China rainfall changes. Here we simulated the annual and seasonal rainfall differences between LGM and 10 ka with regards to the full forcing (TraCE-all) and single forcing (TraCE-ORB, TraCE-GHG and TraCE-ICE) simulations, in order to reveal rainfall changes during glacial-interglacial cycles and their forcing mechanisms. The annual, JJA, MAM and SON rainfall differences in TraCE-all all suggest a drier northern China during LGM than 10 ka (Fig. 3a, e, i). More prominent negative loading of JJA rainfall difference than that during MAM and SON indicates a dominant summer rainfall contribution in a year round in this region. This is consistent with the modern observation of seasonal cycle in northern China rainfall (Fig. 1d).
Coherent variations between the annual and summer rainfall in our 300-ka model simulation also support the leading role of summer rainfall (Supplementary Fig. 6a, b). The JJA rainfall differences over northern China of TraCE-all simulation is similar with that forced by TraCE-ORB and TraCE-GHG, indicating the combination effects of ORB and GHG forcings on local rainfall ( Fig. 3e-g). However, given the dominant glacialinterglacial mode of variability in atmospheric CO 2 concentrations, the dominance of precession cycle in northern China rainfall record implies that insolation plays a leading role in local hydroclimate changes.
Our reconstructed northern China rainfall record shows highly consistent trend with Chinese cave stalagmite δ 18 O record (Fig. 2b-d), which has been considered as an indicator of summer monsoon intensity and is driven by upstream depletion of moisture in an increasing number of works based on data-model approaches in recent years 23,28,32,64,65 . This suggests the linear response of rainfall in northern China to East Asian summer monsoon wind intensity across the glacial-interglacial cycles. Such relationship can also be observed in the modern hydroclimate, which shows highly correlation coefficients between northern China rainfall and summer southerly wind intensiy 9,28 . With increased insolation, the surface temperatures increase more over land than over the oceans, leading to an increased land-sea temperature contrast, which enhances the southerly wind transport of moisture to northern China and thus increasing monsoon rainfall in this region 66 . Our simulated rainfall and moisture flux differences between 56 ka (a high insolation stage in glacial period) and present show stronger JJA rainfall and moisture transport in northern China, providing a strong support for summer monsoon winds could still penetrate northern China and bring rainfall deep inland during high insolation stage in the glaciation (Fig. 4b).
Forcing of southern China rainfall changes. The annual, MAM and SON rainfall differences in TraCE-all suggest stronger rainfall occurred in the most regions of southern China and South China Sea during LGM than 10 ka (Fig. 3a, i). This simulated rainfall anomaly between LGM and 10 ka has also been found in a previous study 67 . The JJA rainfall difference suggests a relatively dry strip extends from southwestern China to northern South China Sea; other regions of southern China and South China Sea were wetter during LGM than 10 ka (Fig. 3e). Thus the relatively stronger local rainfall during glacial periods should result from the combined effects of JJA, MAM and SON rainfall. This is consistent with the modern observation of the seasonal cycle of rainfall changes (Fig. 1d). It has been proposed that rainfall amount in southern China are evenly split between summer and non-summer months 12 . A recent model simulation of hydroclimate conditions during the last deglacial period in southern China also highlighted the role of autumn rainfall 31 . This again emphasizes the importance of proxy interpretation, as summer and annual rainfall should both be carefully evaluated in this region 65 . The similar rainfall responses over this region reflected in the TraCE-all, TraCE-ORB and TraCE-ICE forcings suggest that rainfall was largely driven by a combination of ICE and ORB forcings (Fig. 3a, b, d-f, h-j, l). However, the observed dominant glacial-interglacial mode of local rainfall variations indicates that it should be firstly driven by ice-sheet forcing.
During glacial periods, the expansion of high-latitude ice sheets intensified the Siberian High, resulting in anomalous northwesterly winds 68 , which brought cold and dry air to southern China and the South China Sea and enhanced convergence with warm-wet air masses transported by southwesterly winds in the spring and autumn months, thus strengthening rainfall in these regions (Fig. 5c). Thus dominant glacial-interglacial cycle in annual rainfall probably largely depends on the intervention of spring and autumn rainfall. That's the reason why remarkable differences exist in cyclicity and variation trend between annual and summer rainfall over the southern China, which suggest dominant 100 ka cycle in annual rainfall but 23 ka cycle in summer rainfall, and relatively stronger annual rainfall during glacial periods (Supplementary Fig. 6c, d).
Forcing of Southeast Asian rainfall changes. The JJA, MAM and SON rainfall changes all accounted for the observed weakened glacial rainfall in Southeast Asia (Fig. 3a, e, i). Similar rainfall responses over this region reflected in the TraCE-all and singleforcing simulations suggest that rainfall changes over the western Maritime Continent were forced by a combination of GHG and ICE forcings, whereas ORB forcing was mainly present on the eastern side (Fig. 3a, c-e, g-i, k, l). In addition, ORB forcing also played a role in JJA rainfall of the western Maritime Continent (Fig. 3f).
Previous studies proposed the exposure of Sunda Shelf cooled the surface relative to the surrounding ocean and weakened the local ascending motions and air convection, thus resulting in less rainfall occurring in this region 17,69 . In our simulation, divergence over Sumatra and Borneo of JJA moisture supplied from eastern Indian Ocean during glacials may suppress local rainfall (Fig. 5b). Cooler sea surface temperature in eastern Indian Ocean response to the exposure of Sunda Shelf 70 probably favor such moisture transport. This nearby moisture supply may explain the positive rainfall δ 18 O in our iTraCE results (Supplementary Fig. 4e). During spring and autumn, increased meridional SST gradient associated with highlatitude ice sheet expansion produced and maintained anticyclone in western Pacific 31 . This anticyclone promoted the anomalous  moisture transport to Borneo and Sumatra, led to convergence with moisture derived from eastern Indian Ocean, thus resulted in wetter northern Sunda Shelf and parts of southern Sumatra and northwestern coast of Borneo, but drier other regions of Indonesia (Fig. 5c). Besides, such remote moisture source of western Pacific can explain the negative rainfall δ 18 O values of most regions of Indonesia observed in iTraCE simulation ( Supplementary Fig. 4f).

Conclusions
This study provides a quantitative reconstruction of rainfall in northern China from marine sediments throughout the last four glacial-interglacial cycles. Combined with multiple model simulations, the ideal precession cycle of rainfall variations in northern China has been uncovered. Our results also reveal a dominant glacial-interglacial cycle in southern China-South China Sea and Southeast Asian rainfall changes, with more abundant rainfall in the former region but less rainfall in the latter during glacial periods. In contrast with previous views, we propose a reversed latitude dependence of cyclicity in the East-Southeast Asia hydroclimate, which is mainly controlled by the different responses of seasonal rainfall to insolation, ice sheet and greenhouse gas forcings.
Our study explains many of the existing controversies regarding the orbital cycles associated with rainfall variations among proxies and model simulations. Of particular importance is the highlighting of seasonal contributions to orbital-scale of rainfall changes. The large proportions of spring and autumn rainfall in southern China indicate that rainfall proxies previously established in this region cannot be simply interpreted as summer rainfall indexes. This finding can be used to reexamine the proposed heterogeneous East Asian summer monsoon rainfall variability on multiple time scales, and sheds light on the Asian hydroclimate conditions associated with high-low-latitude climate interactions in the past and future.

Materials and methods
Site location and age model. IODP Site U1429 (31°37.04'N, 128°59.85'E) is located on the western continental slope of the northeastern East China Sea in the western North Pacific at a water depth of 732 m (Supplementary Fig. 1a). At this site, three holes were cored to a composite depth below the seafloor (CCSF-D) of 200 m. The lithology of the section recovered from the study site is mainly composed of calcareous nannofossil ooze and calcareous nannofossil-rich clay, with tephra layers and minor dispersed volcaniclastic materials throughout the succession 71 . The age model covering the last 400 ka at Site U1429 was developed by 11 by correlating the δ 18 O record obtained from benthic foraminifera with the global benthic stack 72 . In addition to the core samples, surface sediment samples were taken from the Yellow River drainage basin (from upper to lower reaches: Guide, Lanzhou, Qingtongxia, Zhaojunmu, Jundu, Huayuankou, Pingyin, Beizhen, Lijin, Kenli) in dry season of 1983 73 (Supplementary Fig. 1a), to exclude the effect of dam construction on chemical composition of sediments 74 .
Major and trace element analysis. Geochemical analyses were performed on the clay-sized (<2 μm) siliciclastic sediments collected from IODP Site U1429 and Yellow River surface sediments. The clay-sized fractions were separated based on Stoke's settling velocity principle and recovered by centrifugation 75 after organic matter and carbonates were removed with hydrogen peroxide (15%) and acetic acid (25%) treatments, respectively. The samples were then rinsed with deionized water three times and dried at 80°C before being ground into powder.
The pretreated sediments were digested by a concentrated HF + HNO 3 + HClO 4 mixture in Teflon vessels for the subsequent elemental analysis. The major and trace element concentrations were determined at the Institute of Oceanology, Chinese Academy of Sciences (IOCAS), using a Thermo Icap6300 inductively coupled plasmaatomic emission spectrometer (ICP-AES) and a Perkin-Elmer ELAN DRC II ICP-mass spectrometer (ICP-MS), respectively 76,77 . Several United States Geological Survey (USGS) and Chinese rock and sediment standards (BCR-2, BHVO-2, GBW07315 and GBW07316) and blanks were repeatedly digested and analyzed in parallel with the samples to monitor the quality of the ICP-MS and ICP-AES measurements. The results were generally within the range of ±6% of the certified values. The analytical precisions were generally better than 1% for major elements and 3% for trace elements.
Climate models and simulations. The fully coupled NCAR-CCSM3 model 26,27,78,79 was employed in this work; this model contains active atmosphere, land, ocean and sea ice components. The model has a spatial resolution of~3.75°f or the atmosphere (26 vertical layers) and land, a variable latitudinal resolution (~0.9°near the equator, becomes coarser poleward) and longitudinal resolution of 3.6°for the ocean (25 vertical layers) and sea ice.
To investigate rainfall changes across East to Southeast Asia in the recent glacialinterglacial cycles, we analyzed transient simulation results (hereafter ORB + GHG + ICE) driven by the response to orbital, greenhouse gas and continental ice-sheet forcings over the last 300 ka. The simulation was initialized from the pre-industrial climate conditions which represents an interglacial state. Then the simulation was forced by transient boundary conditions of the orbital parameters (precession, obliquity and eccentricity), greenhouse gases (CO 2 and CH 4 concentration), and continental ice-sheets (thickness and extent). Notably, the acceleration technique was used to save computational resources. Orbital parameters and greenhouse gas concentrations were multiplied by a factor of 100; therefore, 1 model year in the simulation represents real-time climate evolution over 100 years. The continental ice sheets (mainly over Antarctica, North America and Northern Europe) were updated temporally (every 1000 or 500 years) using the ICE-5G (VM2) reconstruction product 27 . More details of the model set-up and potential caveats can be found in 26 .
To elucidate the forcing mechanisms on different timescales, we also analyzed the transient climate conditions over the last 21,000 years (TraCE-21 ka) using the same fully coupled NCAR-CCSM3 model. In these full TraCE simulation, the model was forced with orbital insolation (ORB), greenhouse gas (CO 2 ), meltwater flux (MWF), and continental ice sheet (ICE) forcings. TraCE simulations were not accelerated. We also focused on three sensitivity experiments driven by ORB (TraCE-ORB), ICE (TraCE-ICE) and GHG (TraCE-GHG) variations to better understand the forcing signals at 23 ka and 10 ka 29 , while the millennial-timescale climate variability induced by fresh water discharge due to melting ice-sheets (reproduced in MWF) was not investigated. All TraCE simulations were conducted at the same low spatial resolution as the ORB + GHG + ICE data 79 . The modeled and post processed atmospheric data used to calculate decadal mean seasonal averages were derived from the Climate Data Gateway at NCAR (https://www. earthsystemgrid.org/). Many model works provide an overview of the model performance of the low-resolution CCSM3 in terms of orbital-timescale monsoon evolution in East Asia and India 11,64,80 . For instance, the model can reasonably reproduce low-level southeasterly penetrating deep into northern China 64 and the spatially heterogeneous of summer monsoon rainfall distribution between northern and southern China 11 .
To test the relationship between the rainfall amount and water isotope (δ 18 O) record in East-Southeast Asia, iTraCE simulations based on iCESM were also conducted. This modeling process was conducted largely following the modeling strategy of TraCE-21 ka but with water isotopes enabled 30,31,81 . The iCESM1.3 model is composed of a community atmosphere model (CAM5.3) 82 , the Parallel Ocean Program, version 2 (POP2) 83 , the Los Alamos Sea Ice Model, version 4  84 , and the Community Land Model, version 4 (CLM4) 85 . The resolution of the atmosphere and land data was approximately 2°(1.9°in latitude and 2.5°in longitude), with 30 vertical levels in the atmosphere; the resolution of the ocean and sea ice data was 1°horizontally (gx1v6), with 60 vertical levels in the ocean. The simulations performed herein covered the climate evolution process from 22 ka to 11 ka and considered the transient forcings associated with the ORB, ICE and GHG forcings (more details regarding these simulations can be found in ref. 32).