Intra-interglacial climate variability : model simulations of Marine Isotope Stages 1 , 5 , 11 , 13 , and 15

Using the Community Climate System Model version 3 (CCSM3) including a dynamic global vegetation model, a set of 13 time slice experiments was carried out to study global climate variability between and within the Quaternary interglacials of Marine Isotope Stages (MISs) 1, 5, 11, 13, and 15. The selection of interglacial time slices was based on different aspects of interand intra-interglacial variability and associated astronomical forcing. The different effects of obliquity, precession, and greenhouse gas (GHG) forcing on global surface temperature and precipitation fields are illuminated. In most regions seasonal surface temperature anomalies can largely be explained by local insolation anomalies induced by the astronomical forcing. Climate feedbacks, however, may modify the surface temperature response in specific regions, most pronounced in the monsoon domains and the polar oceans. GHG forcing may also play an important role for seasonal temperature anomalies, especially at high latitudes and early Brunhes interglacials (MIS 13 and 15) when GHG concentrations were much lower than during the later interglacials. Highversus low-obliquity climates are generally characterized by strong warming over the Northern Hemisphere extratropics and slight cooling in the tropics during boreal summer. During boreal winter, a moderate cooling over large portions of the Northern Hemisphere continents and a strong warming at high southern latitudes is found. Beside the well-known role of precession, a significant role of obliquity in forcing the West African monsoon is identified. Other regional monsoon systems are less sensitive or not sensitive at all to obliquity variations during interglacials. Moreover, based on two specific time slices (394 and 615 ka), it is explicitly shown that the West African and Indian monsoon systems do not always vary in concert, challenging the concept of a global monsoon system on astronomical timescales. High obliquity can also explain relatively warm Northern Hemisphere high-latitude summer temperatures despite maximum precession around 495 ka (MIS 13). It is hypothesized that this obliquity-induced high-latitude warming may have prevented a glacial inception at that time.


Introduction
The Quaternary period is characterized by the cyclic growth and decay of continental ice sheets associated with global environmental changes (see, e.g., Lisiecki and Raymo, 2005;Tzedakis et al., 2006;Jouzel et al., 2007;Lang and Wolff, 2011).While it is commonly accepted that the transitions between glacial and interglacial stages are ultimately triggered by varying astronomical insolation forcings (Hays et al., 1976), climate research is just beginning to understand the internal climate feedbacks that are required to shift the Earth system from one state to the other (see, e.g., van Nes et al., 2015).The astronomical forcing, with its characteristic periods of ca.400 and 100 kyr (eccentricity), 41 kyr (obliquity), and ca.19 and 23 kyr (precession) as in Berger (1978), also acts as an external driver for long-term climate change within the interglacials (i.e. the long-term intra-interglacial climate variability) and likely contributes to interglacial diversity since the evolution of astronomical parameters differs between all Quaternary interglacial stages (cf.Tzedakis et al., 2009;Yin and Berger, 2015).Understanding both interglacial climate diversity and intra-interglacial variability helps to estimate the sensitivity of the Earth system to differ-Published by Copernicus Publications on behalf of the European Geosciences Union.
ent forcings and to assess the rate and magnitude of current climate change relative to natural variability.
Numerous interglacial climate simulations have been performed in previous studies using Earth system models of intermediate complexity (e.g.Kubatzki et al., 2000;Crucifix and Loutre, 2002;Loutre and Berger, 2003;Yin andBerger, 2012, 2015).While the present and the last interglacial have also been extensively investigated with fully coupled atmosphere-ocean general circulation models (e.g.Braconnot et al., 2007;Lunt et al., 2013), earlier interglacial periods have received much less attention by climate modellers.Coupled general circulation model (CGCM) studies of earlier interglacial climates have recently been performed for Marine Isotope Stage (MIS) 11 (Milker et al., 2013;Kleinen et al., 2014) and MIS 13 (Muri et al., 2013).Using the CGCM CCSM3 (Community Climate System Model version 3), Herold et al. (2012) presented a set of interglacial climate simulations comprising the interglaciations of MIS 1,5,9,11,and 19.Their study, however, focussed on peak interglacial forcing (i.e.Northern Hemisphere summer occurring at perihelion) and intercomparison of interglacials (i.e.interglacial diversity) only.In particular, they found that, compared to the other interglacials, MIS 11 exhibits the closest resemblance to the present interglacial, especially during boreal summer.
Here, we present a different and complementary CGCM (CCSM3) study which takes intra-interglacial climate variability into account by simulating two or more time slices for each interglacial stage of MIS 1,5,11,13,and 15.For the interglacial of MIS 5 (last interglacial, MIS 5e; ca.130-115 kyr ago), proxy data suggest a peak global mean temperature of about 1 • C higher than during the pre-industrial period (see, e.g., Otto-Bliesner et al., 2013;Dutton et al., 2015).The maximum global mean sea level has been estimated to 6-9 m above the present-day level (Kopp et al., 2009;Dutton and Lambeck, 2012;Dutton et al., 2015).The interglacial of MIS 11 was unusually long, about 30 000 years (ca.425-395 kyr ago).Global average temperatures of MIS 11 are highly uncertain, but a peak global mean temperature anomaly of up to 2 • C relative to pre-industrial values cannot be ruled out (Lang and Wolff, 2011;Dutton et al., 2015).Maximum global mean sea level may have been 6-13 m higher than today (Raymo and Mitrovica, 2012;Dutton et al., 2015).Interglacials before MIS 11 (early Brunhes interglacials), such as MIS 13 and 15, are generally characterized by lower global mean temperatures, larger continental ice sheets, lower global sea level, and lower atmospheric greenhouse gas (GHG) concentrations relative to the more recent interglacials (see, e.g., Yin and Berger, 2010;Lang and Wolff, 2011;Dutton et al., 2015).
The goal of this study is to disentangle the effects of obliquity, precession, and GHG on global surface climate.Our selection of interglacial time slices takes into account different aspects of inter-and intra-interglacial variability and associated astronomical forcing.As such, our approach dif-fers from and complements previous model studies that focussed on peak interglacial forcing and intercomparison of interglacials (Yin and Berger, 2012;Herold et al., 2012).The selection of the time slices is described in detail in Sect.2.3.
In contrast to previously performed climate model experiments with idealized astronomical forcing, in which obliquity and precession have usually been set to extreme values (see, e.g., Tuenter et al., 2003;Mantsis et al., 2011Mantsis et al., , 2014;;Erb et al., 2013;Bosmans et al., 2015), our analyses are based on realistic astronomical configurations.We note that realistic and idealized forcing experiments are equally important and complementary.Idealized experiments provide important insight into the climate system's response to astronomical forcing.However, since this response may be nonlinear, using extreme values of astronomical parameters in idealized experiments may hide important aspects of astronomical forcing.Obviously, realistically forced experiments have a stronger potential for model-data comparison.
Special focus is on the sensitivity of the West African and Indian monsoon systems to obliquity and precession forcing.In particular, the applicability of the global monsoon concept (Trenberth et al., 2000;Wang et al., 2014) will be tested for astronomical timescales.

Model description
We use the fully coupled climate model CCSM3 with the atmosphere, ocean, sea-ice and land-surface components interactively connected by a flux coupler (Collins et al., 2006).We apply the low-resolution version of the model (Yeager et al., 2006), which enables us to simulate a large set of time slices.In this version, the resolution of the atmosphere is given by T31 spectral truncation (3.75 • transform grid) with 26 layers, while the ocean model has a nominal horizontal resolution of 3 • (as has the sea-ice component) with 25 levels in the vertical.The land model shares the same horizontal grid with the atmosphere and includes components for biogeophysics, biogeochemistry, the hydrological cycle as well as a dynamic global vegetation model (DGVM) based on the Lund-Potsdam-Jena (LPJ) DGVM (Sitch et al., 2003;Levis et al., 2004;Bonan and Levis, 2006).The DGVM predicts the distribution of 10 plant functional types (PFTs) which are differentiated by physiological, morphological, phenological, bioclimatic, and fire-response attributes (Levis et al., 2004).In order to improve the simulation of land-surface hydrology and hence the vegetation cover, new parameterizations for canopy interception and soil evaporation were implemented in the land component (Oleson et al., 2008;Handiani et al., 2013;Rachmayani et al., 2015).PFT population densities are restored annually, while the land and atmosphere models are integrated with a 30 min time step.

Set-up of experiments
To serve as a reference climatic state, a standard preindustrial (PI) control simulation was carried out following PMIP (Paleoclimate Modelling Intercomparison Project) guidelines with respect to the forcing (see, e.g., Braconnot et al., 2007).The PI boundary conditions include astronomical parameters of 1950 AD, atmospheric trace gas concentrations from the 18th century (Table 1) as well as pre-industrial distributions of atmospheric ozone, sulfate aerosols, and carbonaceous aerosols (Otto-Bliesner et al., 2006).The solar constant is set to 1365 W m −2 .The PI control run was integrated for 1000 years starting from modern initial conditions, except for the vegetation which starts from bare soil.
In total, 13 interglacial time slice experiments were carried out, all branching off from year 600 of the PI spin-up run and running for 400 years each.Note that the present study only focusses on the surface climate, for which this spin-up time should be sufficient, whereas the deep ocean usually needs more time to adjust to changes in forcing (Renssen et al., 2006).
Boundary conditions for the selected time slices, which span the last 615 kyr, comprise astronomical parameters (Berger, 1978) and GHG concentrations as given in Table 1, while other forcings (ice-sheet configuration, ozone distribution, sulfate aerosols, carbonaceous aerosols, solar constant) were kept as in the PI control run.The mean of the last 100 simulation years of each experiment was used for analysis.
We note that a fixed calendar based on a 365-day year is used for all experiments (Joussaume and Braconnot, 1997;Chen et al., 2011).The greatest calendar biases are known to occur in boreal fall, whereas the effects in boreal summer and winter (the seasons discussed in the present study) are generally small (see, e.g., Timm et al., 2008).

Selection of interglacial time slices
For MIS 1, the mid-Holocene time slice of 6 ka using standard PMIP forcing (Braconnot et al., 2007) was complemented by an early-Holocene 9 ka simulation when Northern Hemisphere summer insolation was close to maximum (Fig. 1).Two time slices, 125 and 115 ka, were also chosen for the last interglacial (MIS 5e).Similar to 9 ka, the 125 ka time slice is also characterized by nearly peak interglacial forcing, although the MIS 5 insolation forcing is stronger due to a greater eccentricity of the Earth's orbit.Moreover, the global benthic δ 18 O stack is at a minimum around 125 ka (Lisiecki and Raymo, 2005).By contrast, boreal summer insolation is close to a minimum at 115 ka, which marked the end of MIS 5e (Fig. 1).GHG concentrations for the MIS 5 time slices were taken as specified by PMIP3 (Lunt et al., 2013).
For the unusually long interglacial of MIS 11 (see, e.g., Milker et al., 2013), three time slices were chosen: 394, 405, and 416 ka.The middle time slice (405 ka) coincides with the δ 18 O minimum of MIS 11 (Lisiecki and Raymo, 2005;Milker et al., 2013).The time slices of 394 and 416 ka are characterized by almost identical precession and similar GHG concentrations (Table 1) but opposite extremes of obliquity (maximum at 416 ka, minimum at 394 ka; Fig. 1).This allows to study the quasi-isolated effect of obliquity forcing (Berger, 1978) during MIS 11 by directly comparing the results of these two time slices.As opposed to idealized simulations of obliquity forcing (see, e.g., Tuenter et al., 2003;Mantsis et al., 2011Mantsis et al., , 2014;;Erb et al., 2013), our approach considers quasi-realistic climate states of the past using realistic forcings.In the same vein, time slices for MIS 13 have been chosen.Obliquity is at a maximum at 495 ka and at a minimum at 516 ka, while precession is almost identical.Unlike the 394 and 416 ka time slices of MIS 11, which are characterized by intermediate precession values, precession is at a maximum at 495 and 516 ka, i.e.Northern Hemisphere summer occurs at aphelion causing weak insolation forcing (Yin et al., 2009).In addition, the 504 ka time slice was picked because of peak Northern Hemisphere summer insolation forcing, while obliquity has an intermediate value (Fig. 1).Finally, two time slice experiments were performed for MIS 15 to assess the climatic response to minimum (579 ka) and maximum (609 ka) precession.Accordingly, Northern Hemisphere summer insolation is near maximum and minimum at 579 and 609 ka, respectively.In addition, a third MIS 15 experiment was carried out (615 ka) with insolation forcing in between the two others (Fig. 1).Moreover, the 615 ka time slice has a very special seasonal insolation pattern as we will see in the next section.All three MIS 15 time slices coincide with minimum δ 18 O values (Lisiecki and Raymo, 2005).(Lisiecki and Raymo, 2005), climatic precession, obliquity, and insolation in July at 65 • N (Berger, 1978) for the different interglacials.The dots mark the time slices simulated in this study.
Table 1 summarizes the GHG forcing of all experiments with values based on Lüthi et al. (2008), Loulergue et al. (2008), andSchilt et al. (2010) using the European Project for Ice Coring in Antarctica (EPICA) Dome C timescale EDC3, except for the MIS 1 and MIS 5 experiments, where GHG values were chosen following the PMIP guidelines (see above).We note that due to the uneven distribution of methane sources and sinks over the latitudes, values of atmospheric CH 4 concentration derived from Antarctic ice cores present a lower estimate of global CH 4 concentration.We further note that some results from the MIS 1 (6 and 9 ka), MIS 5 (125 ka), and MIS 11 (394, 405, and 416 ka) experiments were previously published (Lunt et al., 2013;Milker et al., 2013;Kleinen et al., 2014;Rachmayani et al., 2015).

Insolation anomalies
Annual cycles of the latitudinal distribution of insolation at the top of the atmosphere (as anomalies relative to PI) are shown in Fig. 2 for each experiment.The insolation patterns can be divided into three groups which differ in their seasonal distribution of incoming energy.Group I is characterized by high Northern Hemisphere summer insolation as exhibited for the 6 and 9 ka (MIS 1), 125 ka (MIS 5), 405 and 416 ka (MIS 11), 504 ka (MIS 13), and 579 ka (MIS 15) time slices.In most (but not all, see below) cases this is due to an astronomical configuration with northern summer solstice at or close to perihelion.Group II comprises anomalies with low boreal summer insolation as shown for 115 ka (MIS 5), 495 and 516 ka (MIS 13), and 609 ka (MIS 15).In these cases, northern winter solstice is near perihelion.Group III is characterized by changes in the sign of the Northern Hemisphere insolation anomalies from spring to summer and consists of two dates (394 and 615 ka).At 394 (615 ka) the insolation anomaly spring-to-summer change is from positive (negative) to negative (positive).In these cases, spring equinox (394 ka) or fall equinox (615 ka) are close to perihelion.

JJAS surface temperature anomalies
The response of boreal summer (June-July-August-September, JJAS) surface temperature to the combined effect of insolation and GHG in all individual climates (Fig. 3) shows warm conditions (relative to PI) over most parts of the continents in Group I (6, 9, 125, 405, 416, 504, and 579 ka), with the three warmest anomalies at 9, 125, and 579 ka.The warm surface conditions can largely be explained by the immediate effect of high summer insolation and a reduction in the Northern Hemisphere sea-ice area by about 15-20 % (not shown) relative to PI.The large thermal capacity of the ocean explains a larger temperature response over land than over the ocean (Herold et al., 2012;Nikolova et al., 2013).Simulated cooling over northern Africa (10-25 • N) and India in the Group I experiments is caused by enhanced monsoonal rainfall in these regions, which is associated with increased cloud cover, i.e. reduced shortwave fluxes, and enhanced land-surface evapotranspiration, i.e. greater latent cooling (see, e.g., Braconnot et al., 2002Braconnot et al., , 2004;;Zheng and Braconnot, 2013).Cooling in some parts of the Southern Ocean in most Group I experiments is likely attributable to an austral summer remnant effect of local insolation (see below) as in Yin and Berger (2012).The 416 ka time slice, however, differs from the other Group I members by anomalously cold conditions over the Southern Hemisphere continents.Again, this behaviour can be explained by the immediate effect of the insolation, which shows negative anomalies in the Southern Hemisphere during the JJAS season (Fig. 2).As such, the 416 ka time slice must be considered a special case in Group I.While high Northern Hemisphere summer insolation is related to low precession in most Group I members, positive anomalies of Northern Hemisphere summer insolation at 416 ka are attributable to a maximum in obliquity (Fig. 1), yielding the Northern versus Southern Hemisphere insolation contrast.
In contrast to Group I, Group II climates exhibit anomalously cold JJAS surface temperatures globally with the three coldest anomalies at 115, 516, and 609 ka.Again, the temperature response can largely be explained by the direct response to insolation forcing, amplified at high latitudes by an increase in the sea-ice cover (about 5 % in the Arctic compared to PI).Due to a particular combination of high precession and eccentricity with low obliquity, the insolation forcing and surface temperature response is strongest for the 115 ka time slice.Group II warming in the northern African and Indian monsoon regions is associated with increased aridity and reduced cloudiness.
Group III climates (394 and 615 ka) show rather complex temperature anomaly patterns, especially in the tropics.In the 394 ka time slice, however, northern continental regions show a distinct cooling, whereas continental regions exhibit an overall warming in the Southern Hemisphere (except for Antarctica).To a large extent, the 394 ka time slice shows a reversed JJAS temperature anomaly pattern compared to the 416 ka simulation over the continental regions, except for Antarctica.

DJF surface temperature anomalies
Boreal winter (December-January-February, DJF) surface temperature anomalies are presented in Fig. 4. Generally low DJF insolation in Group I time slices (Fig. 2) results in anomalously cold surface conditions over most of the globe, particularly strong in the 579 ka (MIS 15) time slice.However, anomalously warm conditions in the Arctic stand in contrast to the global DJF cooling at 6, 9, 125, 405, and 416 ka.The Arctic warming is due to the remnant effect of the polar summer insolation through ocean-sea-ice feedbacks (Fischer and Jungclaus, 2010;Herold et al., 2012;Yin and Berger, 2012;Kleinen et al., 2014).Anomalous shortwave radiation during the Arctic summer leads to enhanced melting of sea ice and warming of the upper polar ocean.
The additional heat received by the upper ocean delays the formation of winter sea ice, reduces its thickness, and finally leads to a warming of the winter surface atmospheric layer by enhanced ocean heat release (Yin and Berger, 2012).Arctic winter warming is not present in the 504 ka (MIS 13) and 579 ka (MIS 15) time slices in Group I, where the summer remnant effect in the Arctic is probably masked by a global cooling that is induced by low GHG concentrations typical for early Brunhes (MIS 13 and before) interglacials.
To a large extent, DJF surface temperature anomaly patterns are reversed in Group II, with warming over most continental regions.Moreover, the summer remnant effect reverses to a substantial cooling in the Arctic region.Temperature anomaly patterns in Group III are, again, rather complex.Interestingly, most Northern Hemisphere continental regions remain relatively cold during boreal winter (as in summer) in the 394 ka simulation.Relatively low GHG concentrations, especially CH 4 , contribute to the year-round extratropical cooling in this time slice.

JJAS precipitation anomalies
Boreal summer precipitation shown in Fig. 5 exhibits intensified rainfall in the monsoon belt from northern Africa to India, via the Arabian Peninsula, in all Group I simulations in response to high summer insolation (Prell and Kutzbach, 1987;de Noblet et al., 1996;Tuenter et al., 2003;Braconnot et al., 2007).By contrast, the same monsoon regions experience anomalously dry conditions in the Group II (low boreal summer insolation) experiments.The most interesting results regarding the tropical rainfall response to astronomical forcing appear in Group III, where the monsoonal precipitation anomalies show opposite signs in northern Africa (Sahel region) and India.

Net primary production (NPP) anomalies
Vegetation responds to changes in surface temperature and precipitation and, in certain regions, may feed back to the climate (cf.Rachmayani et al., 2015).Figure 6 shows the simulated changes in NPP, reflecting increase or decrease in and expansion or retreat of vegetation covers, relative to PI.At high Arctic latitudes, NPP increases in the Group I simulations, except for 405 ka, where temperature changes are probably too small to substantially affect the vegetation.By contrast, Arctic NPP declines in the Group II experiments, albeit only in the easternmost part of Siberia in the 495 ka experiment.A substantial decline in Arctic NPP is also simulated for 394 ka (Group III).In the tropical regions, vegetation changes are mostly governed by precipitation.Consequently, enhanced rainfall results in increased NPP over northern Africa, the Arabian Peninsula and India in all Group I experiments.In northern Africa increased NPP is associated with a northward shift of the Sahel-Sahara boundary.The largest shifts are simulated for 125 and 579 ka in accordance with maximum northern African rainfall anomalies.In these experiments, a complete greening of the Arabian Desert is simulated.Opposite NPP anomalies in the tropical monsoon regions are simulated in the Group II experiments.In Group III, NPP increases result from anomalously high rainfall in northern Africa (615 ka) or India (394 ka).

Climatic effects of obliquity variations during MIS 11 and MIS 13
The MIS 11 time slices 394 and 416 ka show opposite obliquity extremes (at similar precession), as do the MIS 13 time slices 495 and 516 ka (Fig. 1).Insolation differences between the high-obliquity (416, 495 ka) and low-obliquity (394, 516 ka) cases (i.e.416−394 and 495−516 ka) are dis-played in Fig. 7.The effect of high obliquity is to strengthen the seasonal insolation cycle.At low latitudes, the effect of obliquity on insolation is small.For the maximum obliquity time slices (416 and 495 ka) relatively high boreal summer insolation directly translates into positive surface temperature anomalies over Northern Hemisphere continents, except for the low latitudes where reduced local insolation (especially in the MIS 13 case) and higher monsoon rainfall (especially in the MIS 11 case, see below) lead to surface cooling (Fig. 8a, b).By contrast, receiving anomalously low insolation during austral winter, Southern Hemisphere continents exhibit anomalously cold surface temperatures.For the 416-394 ka case, however, the Antarctic continent and the Southern Ocean show large-scale warming during the JJAS season, which can be attributed to a south polar summer remnant effect as the austral summer insolation anomaly is extremely high in this experiment (Fig. 7a).Higher GHG concentrations at 416 compared to 394 ka may add to this warming.Owing to a smaller south polar summer insolation anomaly (Fig. 7), the summer remnant effect is smaller in the 495-516 ka case and even surpassed by anomalously low GHG forcing in the 495 ka time slice, leading to negative austral winter temperature anomalies in the Southern Ocean and Antarctica (Fig. 7b).
During boreal winter, Northern Hemisphere continents show large-scale cooling in response to high obliquity (and hence relatively low insolation), except for the Arctic realm where the summer remnant effect results in substantial positive surface temperature anomalies (Fig. 8c, d).During the same season (DJF), anomalously high insolation causes surface warming in the Southern Hemisphere in response to high obliquity.As a general pattern in the annual mean, maximum-minus-minimum obliquity forcing causes anomalous surface warming at high latitudes and surface cooling at low latitudes caused by seasonal and annual insolation anomalies in combination with climate feedbacks like the polar summer remnant effect and monsoon rainfall.
Despite the weak insolation signal at low latitudes, substantial obliquity-induced changes in tropical precipitation are simulated (Fig. 8e, f).The strongest signal is found in the northern African monsoon region in the MIS 11 experiments, where greater JJAS precipitation occurs during maximum obliquity at 416 ka than during the obliquity minimum at 394 ka.A positive Sahel rainfall anomaly is also found in the MIS 13 experiments (495-516 ka), but it is much weaker than in the MIS 11 case (416-394 ka).We suppose that the obliquity-induced increase in northern African monsoonal rainfall is counteracted by the high precession at 495 ka that tends to weaken the monsoon.Considering the spatiotemporal insolation patterns (Fig. 7), the Northern Hemisphere tropical summer insolation anomaly is less negative and the meridional summer insolation gradient anomalies are generally greater in the 416-394 ka case compared to the 495-516 ka case.Both features of the insolation anomaly favour a strong northern African monsoon (see the "Discussion" section).

Evaluating the climatic effects of astronomical and GHG forcings through correlation maps
In order to evaluate the climatic effects of obliquity, precession and GHG concentrations, linear correlations between the individual forcing parameters and climatic fields (surface temperature, precipitation) were calculated from the 14 time slice experiments (13 interglacial time slices plus PI).To this end, each climate variable (temperature, precipitation) was averaged over the last 100 years of each experiment.Linear correlation coefficients between a climatological variable and a forcing parameter (obliquity, precession, GHG radiative forcing) were calculated at each grid point.Significance of correlations was tested by a two-sided Student's t test with 95% confidence level.Total radiative forcing from CO 2 , CH 4 , and N 2 O in each experiment was calculated based on a simplified expression given in Table 3 (IPCC, 2001).
Figure 9 shows the corresponding correlation maps for annual mean, boreal summer, and boreal winter surface temperature.As expected, GHG forcing is positively correlated with surface temperature over most regions of the globe (Fig. 9a), which is particularly pronounced in the annual mean.For the seasonal correlation maps (boreal summer and winter) the correlation coefficients are smaller because of the dominant impact of obliquity and precession forcing.
As already described in the previous subsection, the general surface temperature pattern of high-obliquity forcing is warming at high latitudes and cooling at low latitudes  (Fig. 9b).High precession (northern solstice near aphelion) leads to boreal summer surface cooling over most extratropical regions (Fig. 9c).However, surface warming occurs in some tropical regions as a response to weaker monsoons.During boreal winter, anomalously high insolation causes anomalous surface warming except in the Arctic (due to the summer remnant effect) and northern Australia (due to a stronger regional monsoon).
Correlation maps for annual mean, boreal summer, and boreal winter precipitation are shown in Fig. 10.GHG radiative forcing exhibits no clear response in precipitation except for the high latitudes where the hydrologic cycle accelerates with higher GHG concentrations (Fig. 10a).Since the GHG variations are relatively small, the effects of astronomical forcing on the monsoons are way larger than the effects of GHG variations during the interglacials.Arctic precipitation is also amplified by high obliquity during summer (Fig. 10b).Obliquity also strengthens the monsoonal rainfall in northern Africa (Sahel region), whereas no effect of obliquity can be detected for the Australian monsoon.The sensitivity of other monsoon systems to obliquity changes is also weak or even absent in our experiments.The most robust response of the hydrologic cycle is found for precession (Fig. 10c).In particular, high precession reduces summer rainfall in the monsoon belt from northern Africa to India as well as in the Arctic realm.East Asian rainfall shows a somewhat heterogeneous pattern and is, in general, only weakly coupled with the Indian and African monsoons.This finding is consistent with a recent model intercomparison study by Dallmeyer et al. (2015).During boreal winter, the hydrologic cycle strengthens in the Arctic and Antarctic regions, while Southern Hemisphere monsoon systems amplify, resulting in enhanced rainfall over South America, southern Africa, and northern Australia in response to high precession.We note that these monsoonal rainfall changes go along with distinct surface temperature signals in the annual mean (Fig. 9c).

Discussion
While most time slices presented in this study were simulated for the first time using a comprehensive CGCM, the 6, 115, and 125 ka time slices have been studied extensively in previous model studies.In general, the CCSM3 results are in line with these previous studies in terms of large-scale temperature and precipitation patterns.Warm boreal summer conditions (relative to PI) over most parts of the continents and the Arctic are a general feature in paleoclimatic simulations of the mid-Holocene (6 ka), while the northern African and South Asian monsoon regions are anomalously cold due to enhanced rainfall (Braconnot et al., 2007).Though evidenced by proxy records (see, e.g., McClure, 1976;Hoelzmann et al., 1998;Fleitmann et al., 2003), several models fail to simulate wetter mid-Holocene conditions over the Arabian Peninsula (cf.https://pmip3.lsce.ipsl.fr/database/maps/),while CCSM3 simulates not only enhanced rainfall but also greening of the Arabian Desert.The 125 ka surface temperature pattern shows similar features to the 6 ka pattern, but they are much more pronounced due to the larger eccentricity and hence stronger precessional forcing.However, compared to other simulations of the last interglaciation, our CCSM3 simulation produces a relatively cold MIS 5e surface climate as shown by Lunt et al. (2013).At 115 ka, surface tempera- ture anomalies show the opposite sign with dramatic cooling over the Arctic and the northern continental regions providing ideal conditions for glacial inception (see, e.g., Khodri et al., 2005;Kaspar and Cubasch, 2007;Jochum et al., 2012).A retreat of the vegetation at high northern latitudes tends to amplify the insolation-induced cooling (cf.Gallimore and Kutzbach, 1996;Meissner et al., 2003).
A recent simulation of the MIS 13 time slice at 506 ka using the CGCM Hadley Centre Coupled Model version 3 (HadCM3) (Muri et al., 2013) can be compared to our 504 ka time slice using CCSM3.Global patterns of surface temperature anomalies (relative to PI) are remarkably similar in the two different simulations with warm anomalies over all continents (except for the northern African and South Asian monsoon regions) in boreal summer and worldwide cold anomalies during boreal winter.Moreover, both simulations show anomalously high boreal summer precipitation over northern South America, North and central Africa as well as the South Asian monsoon region.
Although our CCSM3 results show general agreement with other model studies, the validation of model results with data is usually not straightforward.The reader is referred to previous work where our CCSM3 simulation of 125 ka (Lunt et al., 2013) as well as the MIS 11 simulations have been extensively compared to proxy data (Milker et al., 2013;Kleinen et al., 2014).Taken together, these and other stud- ies (e.g.Lohmann et al., 2013) indicate that CGCMs tend to produce generally smaller interglacial temperature anomalies than suggested by the proxy records.So far, the reason for these discrepancies is unsolved (cf.Liu et al., 2014), but Hessler et al. (2014) pointed out that uncertainties associated with sea surface temperature reconstructions are generally larger than interglacial temperature anomalies.Thus, currently available surface temperature proxy data cannot serve as a target for benchmarking interglacial model simulations.
Two time slices of MIS 11 (394 vs. 416 ka) and two time slices of MIS 13 (495 vs. 516 ka) allow the investigation of (almost pure) obliquity effects on global climate, although the GHG and precession are not exactly the same in the different time slices.The results from these simulations can be compared to previously performed idealized model experiments in which obliquity has been changed from maximum to minimum values (Tuenter et al., 2003;Mantsis et al., 2011;Erb et al., 2013;Bosmans et al., 2015).The common results of those idealized model experiments and our experiments can be summarized as follows.High-versus low-obliquity climates are characterized by strong warming over the Northern Hemisphere extratropics and slight cool-ing in the tropics during boreal summer.During boreal winter, a moderate cooling over large portions of the Northern Hemisphere continents and a strong warming at high southern latitudes is found.The obliquity-induced Northern Hemisphere summer warming appears to be of particular interest for the MIS 13 climate evolution.At 495 ka, precession is at a maximum, but the global benthic δ 18 O stack by Lisiecki and Raymo (2005) does not show the expected increase towards heavier values which would indicate colder conditions and Northern Hemisphere cryosphere expansion (Fig. 1).In fact, despite high precession, the 495 ka simulation exhibits the warmest Northern Hemisphere summer temperatures from all Group II experiments (Fig. 3), which can be attributed to concomitant high obliquity.We therefore hypothesize that the Northern Hemisphere summer climate at 495 ka was not cold enough for ice sheets to grow and global ocean δ 18 O to increase.We note, however, that the benthic δ 18 O stack is subject to age model uncertainties of a few thousand years.
Moreover, our CCSM3 results as well as the studies by Tuenter et al. (2003) and Bosmans et al. (2015) suggest a significant effect of obliquity on West African monsoon rainfall despite the weak insolation signal at low latitudes.Bosmans et al. (2015) have shown that obliquity-induced R. 5,11,13,and  changes in moisture transport towards northern Africa result from changes in the meridional insolation gradient (Davis and Brewer, 2009).However, the impact of obliquity on the monsoon also depends on precession.In the 495-516 ka experiment the obliquity effect on the West African monsoon is minor, as both time slices (495 and 516 ka) are characterized by precession maxima leading to extremely weak monsoonal circulation and rainfall in both cases.The existence of a ∼ 41 kyr cyclicity (in addition to astronomically related ∼ 100 and 19-23 kyr cycles) in reconstructions of northern African aridity during the Quaternary has usually been attributed to obliquity-forced Northern Hemisphere cryosphere effects on the monsoon climate (see, e.g., Bloemendal and deMenocal, 1989;deMenocal et al., 1993;Tiedemann et al., 1994;deMenocal, 1995;Kroon et al., 1998).Our model results along with the studies by Tuenter et al. (2003) and Bosmans et al. (2015) complement this picture, showing that the direct insolation gradient forcing associated with obliquity can contribute to West African monsoon changes without involving high-latitude remote climate forcing associated with Northern Hemisphere ice sheets.
According to the CCSM3 results, the Indian monsoon is less sensitive to direct obliquity (insolation gradient) forcing than the West African monsoon.This finding is consis-tent with proxy records from the Arabian Sea that show substantial 41 kyr (obliquity) periodicity only after the onset of Quaternary glacial cycles when the waxing and waning of northern ice sheets could have worked as an agent for the transfer of obliquity forcing to the Indian monsoon region (Bloemendal and deMenocal, 1989).In general, it is found that the two monsoon systems do not always vary in concert.This is particularly evident in the Group III experiments (394 and 615 ka), where the precipitation anomalies over northern Africa and India have opposite signs (Table 2).Considering the annual insolation maps of the 394 and 615 ka experiments (Fig. 2), West African monsoon rainfall turns out to be most sensitive to changes in summer insolation, whereas spring or early summer insolation is more important for monsoon rainfall over India.Similar results have been found by Braconnot et al. (2008).It has been argued that the reason is a resonant response of the Indian monsoon to the insolation forcing when maximum insolation anomalies occur near the summer solstice and a resonant response of the African monsoon -which has its rainfall maximum 1 month later in the annual cycle than the Indian monsoon -when the maximum insolation change is delayed after the summer solstice.The different responses to specific forcings and the sometimes out-of-phase behaviour of the African and Indian monsoon systems challenge the global monsoon concept -according to which all regional monsoon systems are part of one seasonally varying global-scale atmospheric overturning circulation in the tropics (Trenberth et al., 2000;Wang et al., 2014) -on astronomical timescales.
Another important result of our study is associated with obliquity forcing of high-latitude precipitation anomalies.As obliquity increases, high latitudes become warmer and the gradient in solar heating between high and low latitudes decreases, while precipitation over high-latitude continental regions increases (Fig. 10b).This result clearly contradicts the "gradient hypothesis" by Raymo and Nisancioglu (2003), according to which low obliquity would favour polar ice-sheet growth through enhanced delivery of moisture owing to an increased meridional solar heating gradient.
Since CO 2 and other GHG variations are relatively small during the interglacials, the effects of astronomical forcing on the monsoons are substantially larger.Hence, GHG forcing shows a clear response in precipitation only for the high latitudes where the hydrologic cycle accelerates with higher GHG concentrations.In the monsoon regions, interglacial rainfall variations are almost entirely controlled by astronomical forcing.
The use of a modern ice-sheet configuration for all interglacial time slice experiments, however, must be considered a limitation of the present study.Future studies should include the effects of changing ice sheets and associated meltwater fluxes in shaping interglacial climates.Large Northern Hemisphere ice sheets might have played an important role for regional and global climates, especially during early Brunhes interglacials (MIS 13 and before) as suggested by, e.g., Qiuzhen Yin et al. (2008) and Muri et al. (2013).Model studies suggest an influence of changing land ice on the interglacial climate evolution also during late Brunhes interglacial stages, like the Holocene (Renssen et al., 2009;Marzin et al., 2013).The tremendous uncertainties regarding ice-sheet reconstructions beyond the present interglacial could be taken into account by performing sensitivity experiments.

Conclusions
Using CCSM3-DGVM, 13 interglacial time slice experiments were carried out to study global climate variability between and within Quaternary interglacials.The selection of interglacial time slices was based on different aspects of inter-and intra-interglacial variability and associated astronomical forcing.As such, our approach is complementary to both idealized astronomical forcing experiments (e.g.Tuenter et al., 2003;Mantsis et al., 2011Mantsis et al., , 2014;;Erb et al., 2013;Bosmans et al., 2015) and climate simulations that focussed on peak interglacial forcing (Herold et al., 2012;Yin and Berger, 2012).
In this study, the different roles of obliquity, precession, and GHG forcing on surface temperature and precipitation patterns have been disentangled.In most regions seasonal surface temperature anomalies could largely be explained by local insolation anomalies induced by the astronomical forcing.Climate feedbacks modify the surface temperature response in specific regions, particularly in the monsoon domains and the polar oceans.GHG forcing may also play a role for seasonal temperature anomalies, especially at high latitudes and in the early Brunhes interglacials MIS 13 and 15, when GHG concentrations were much lower than during the later interglacials.
A significant role of obliquity in forcing the West African monsoon was found, whereas the Indian monsoon -as well as the other regional monsoon systems -appear to be less sensitive (or not sensitive at all) to obliquity changes during interglacials.Despite this important role of obliquity in West African monsoon variability, the response to precession is still stronger.Different responses to specific forcings and the obvious anti-phase behaviour of the African and Indian monsoon systems in the 394 and 615 ka experiments, where the northern African rainfall anomaly has the opposite sign compared to the Indian anomaly, clearly point to the fact that the two regional monsoon systems do not always vary in concert and challenge the global monsoon concept on the astronomical timescale.
As a general pattern in the annual mean, maximum-minusminimum obliquity forcing causes anomalous surface warming at high latitudes and surface cooling at low latitudes caused by seasonal and annual insolation anomalies in combination with climate feedbacks such as the polar summer remnant effect and monsoon rainfall.High obliquity may also explain relatively warm Northern Hemisphere highlatitude summer temperatures despite maximum precession around 495 ka (MIS 13).We hypothesize that this obliquityinduced high-latitude warming may have prevented a glacial inception at that time.Moreover, our results suggest highlatitude precipitation increase with increasing obliquity, contradicting the "gradient hypothesis" by Raymo and Nisancioglu (2003), according to which low obliquity would favour polar ice-sheet growth through enhanced delivery of moisture owing to an increased meridional solar heating gradient.
Future studies should include the effects of changing ice sheets and associated meltwater fluxes in shaping interglacial climates.With increasing computer power, long-term transient simulations of interglacial climates will become more common.So far, transient CGCM simulations have been performed for the present (see, e.g., Lorenz and Lohmann, 2004;Varma et al., 2012;Liu et al., 2014) and the last interglacial (see, e.g., Bakker et al., 2013;Govin et al., 2014).More transient simulations of earlier interglacials, ideally with coupled interactive ice-sheet models, will help to develop a significantly deeper understanding of interglacial climate dynamics.

Figure 2 .
Figure 2. Insolation anomalies (relative to PI) for the time slices simulated in this study.Patterns of insolation anomaly are classified as Groups I, II, and III (see text).The calculation assumes a fixed present-day calendar with vernal equinox on 21 March.

Figure 3 .
Figure 3. Boreal summer surface temperature anomalies (relative to PI) for the different interglacial time slices.Classification as Groups I, II, and III (see text) is indicated.

Figure 5 .
Figure 5.As in Fig. 3, but for boreal summer precipitation.

Figure 6 .
Figure6.As in Fig.3, but for annual net primary production.

Figure 9 .
Figure 9. Linear correlation maps between surface temperature and GHG radiative forcing (a), obliquity (b), and climatic precession (c) as calculated from the entire set of experiments.Summer refers to JJAS, winter to DJF.Only significant values are shown, according to a two-sided Student's t test at 95 % confidence level.
As in Fig.9, but for precipitation.

Table 1 .
Atmospheric GHG concentrations used in the interglacial experiments.

Table 2 .
Summer (JJAS) precipitation over northern Africa (20 • W-30 • E and 10-25 • N) and over India (70-100 • E and 10-30 • N) along with anomalies relative to PI. Absolute precipitation values are given with standard error (2σ ) based on 100 simulation years of each experiment.

Table 3 .
Simplified expressions for calculation of radiative forcing due to CO 2 , CH 4 , N 2 O. C is CO 2 in ppmv, M is CH 4 in ppbv, and N is N 2 O in ppbv.The subscript 0 denotes the unperturbed GHG concentration of PI.