The South American Monsoon Variability over the Last Millennium in Climate Models

In this paper we assess South American monsoon system (SAMS) variability in the last millennium as depicted by global coupled climate model simulations. High-resolution proxy records for the South American monsoon over this period show a coherent regional picture of a weak monsoon during the Medieval Climate Anomaly and a stronger monsoon during the Little Ice Age (LIA). Due to the small external forcing during the past 1000 years, model simulations do not show very strong temperature anomalies over these two specific periods, which in turn do not translate into clear precipitation anomalies, in contrast with the rainfall reconstructions in South America. Therefore, we used an ad hoc definition of these two periods for each model simulation in order to account for model-specific signals. Thereby, several coherent large-scale atmospheric circulation anomalies are identified. The models feature a stronger monsoon during the LIA associated with (i) an enhancement of the rising motion in the SAMS domain in austral summer; (ii) a stronger monsoon-related upper-tropospheric anticyclone; (iii) activation of the South American dipole, which results in a pole-ward shift of the South Atlantic Convergence Zone; and (iv) a weaker upper-level subtropical jet over South America. The diagnosed changes provide important insights into the mechanisms of these climate anomalies over South America during the past millennium.


Introduction
It is well established that monsoon systems respond to orbital forcing (Kutzbach and Liu, 1997;Kutzbach et al., 2007;Bosmans et al., 2012). At orbital timescales (especially related to the precessional cycle of ca. 19 and 21 kyr), changes in the latitudinal insolation gradient, and hence temperatures, force the monsoon circulation globally (e.g. Bosmans et al., 2012). In the precession frequency band the summer insolation is in anti-phase between hemispheres (for example, when Northern Hemisphere (NH) summer insolation is at its maximum, summertime insolation in the Southern Hemisphere (SH) is at its minimum). This results in a weakening of the monsoon circulation and precipitation in one hemisphere while in the other the monsoon is strengthened. The mechanism for the orbitally induced monsoon variability is therefore mainly related to meridional temperature gradients. Thus, it is not surprising that other phenomena that produce important changes in hemispheric temperature gradients are also responsible for monsoon variability. Examples of these are abrupt Dansgaard-Oeschger events during the last glacial Cheng et al., 2013) and Heinrich events, including the Heinrich 1 event, during the last deglaciation (ca. 17 ka; e.g. Griffiths et al., 2013;Deplazes et al., 2014;Cruz et al., 2006;Strikis et al., 2015).
In recent years, similar variability has also been observed for shorter timescales, in particular between the two most prominent climate anomalies over the last millennium (LM), the Medieval Climate Anomaly (MCA, ca. 950-1250 CE) and the Little Ice Age (LIA, ca. 1450-1850e.g. Masson-Delmotte et al., 2013a). Recent high-resolution records from the area of the South American monsoon system (SAMS) domain have been used to reconstruct precipitation over this region. Records include speleothems (Novello et al., , 2016Kanner et al., 2013;Apaestegui et al., 2014), pollen (Ledru et al., 2013), lake sediments (Bird et al., 2011), and tree-ring reconstructions (Morales et al., 2012).  reviewed currently available proxy records for the SAMS region. Most reconstructions show good correlations with NH temperature and Intertropical Convergence Zone (ITCZ) reconstructions. According to these palaeoclimate studies, the LIA was characterized by a cool north equatorial Atlantic and a warm south equatorial Atlantic (Haug et al., 2001;Polissar et al., 2006), whereas an opposite pattern was present during the MCA. This meridional temperature gradient led to a southward (northward) migration of the Atlantic ITCZ during the LIA (MCA; Haug et al., 2001). Indeed, SAMS reconstructions during the last millennium show a weaker monsoon during the MCA period and a relatively stronger monsoon during the LIA period (e.g. Bird et al., 2011;Vuille et al., 2012;Ledru et al., 2013;Apaestegui et al., 2014), indicating an anti-correlation with reconstructions of the Southeast Asian monsoon (Zhang et al., 2008;Shi et al., 2014;Polanski et al., 2014), as well as with the North African and North American monsoons (Asmerom et al., 2013), for those periods.
Moreover, modelling studies support a southward (northward) shift of the Atlantic ITCZ during LIA (MCA). For instance, model simulations by Vellinga and Wu (2004) suggest that anomalous northward ocean heat transport during the MCA was linked to an enhanced cross-equatorial temperature gradient in the Atlantic and a northward movement of the ITCZ. Kageyama et al. (2013) analysed freshwater hosing simulations over the North Atlantic to force fluctuations in the strength of the Atlantic Meridional Overturning Circulation (AMOC). Their analyses suggest that the model response to an enhanced high-latitude freshwater flux is characterized by a general cooling of the North Atlantic, a southward shift of the Atlantic ITCZ, and a weakening of the African and Indian monsoons. Furthermore, modelling experiments discussed by Broccoli et al. (2006) and Lee et al. (2011) indicate that cooler-than-normal temperatures imposed in the North Atlantic domain, as occurred during the LIA, shift the Atlantic ITCZ southward. In their experiments, this shift is related to a strengthening of the northern Hadley cell in austral summer and a slight shift in its rising branch to the south. Thus, a number of palaeoclimate reconstructions and modelling studies suggest that the particular temperature anomalies observed during the MCA and LIA periods, especially in the North Atlantic, were large enough to modify the location of the ITCZ over the tropical Atlantic, thereby affecting the strength of the summer SAMS throughout the past millennium (see also a review by Schneider et al., 2014).
Recent climate modelling experiments for the LM (850-1850 CE) have been incorporated in the third phase of the Palaeoclimate Modelling Intercomparison Project (PMIP3). About a dozen models included in the Climate Model Intercomparison Project Phase Five (CMIP5) ran this experiment, which considers solar, volcanic, greenhouse gases, and landuse scenarios during the LM (Schmidt et al., 2011(Schmidt et al., , 2012. In this paper, we explore if and how these coupled general circulation model (GCM) simulations capture the variability of the SAMS associated with LIA and MCA temperature anomalies, as suggested by rainfall reconstructions and diverse modelling studies in the region. This evaluation provides further insights regarding the response of the current generation of GCMs to external forcing during the LM. We focus on the models' ability to simulate the variability of the main characteristics of the South American climate during the MCA and LIA. These characteristics are analysed by concentrating on three main features: precipitation, temperature, and atmospheric circulation.
This paper is organized as follows: Sect. 2 presents a short description of the model simulations considered and the method used to identify the MCA and LIA periods; Sect. 3 presents the main results from the climate simulations of the SAMS during both periods; and Sect. 4 presents a discussion and the main conclusions from this study.

Methods and model simulations
We use nine available coupled GCMs, eight of which correspond to CMIP5/PMIP3 LM simulations and one that does not follow the exact PMIP3 experimental setup. The model simulations are listed in Table 1. These simulations cover the period 850-1850 CE, although some of them continue up to the present. But since not all modelling groups have continuous runs to the present (including the period 1850-2000) available, the analysis in this paper covers only the period until 1850 CE. The LM simulations were forced with orbital variations (mainly shifts in the perihelion date), common solar irradiance, two different volcanic eruption reconstructions, land-use change, and greenhouse gas (GHG) concentrations. A full description of the exact forcings used in these LM simulations is given by Schmidt et al. (2011Schmidt et al. ( , 2012. Furthermore, a detailed list of individual forcings applied in each simulation is given in Annex 2 of Masson-Delmotte et al. (2013a).

Definition of periods
The Intergovernmental Panel on Climate Change (IPCC) Fifth Assessment Report (AR5; IPCC, 2013) defined the two periods of most prominent climate anomalies over the past millennium as the MCA (950-1250 CE) and the LIA (1450-1850 CE). This report also concluded that the MCA was a period of relative global warmth, albeit in general less homogenous than the current warmth, whereas the LIA  was a much more globally uniform cold period (Masson-Delmotte et al., 2013a). A recent analysis of the consistency of the CMIP5/PMIP3 LM temperature simulations indicates that these simulations often differ from available temperature reconstructions in their long-term multi-centennial trends, which is related to the transition from the MCA to the LIA period (Bothe et al., 2013;Fernández-Donado et al., 2013). Figure 1a shows the NH temperature anomaly time series for each of the nine models considered, as well as its ensemble mean. For comparison, reconstructions used in Fig. 5.7 of the fifth IPCC report are shown (Masson-Delmotte et al, 2013a, b). From Fig. 1a it is clear that the temperature anomalies over the last millennium are small, and clearly common MCA and LIA periods are not easily identifiable across models. This is consistent with the notion that both periods are partially a result of internal climate variability (PAGES 2k Consortium, 2013), particularly the MCA (Neukom et al., 2014). The hypothesis that guides our analysis used to asses the SAMS variability in models is that both periods are primarily a result of internal (non-forced) variability In addition, given that not all GCM simulations used the exact same forcing, we cannot expect the models to exactly reproduce the temporal variability as indicated by the reconstructions. Therefore, we identify these two periods individually in each model, using two criteria. First, for each model, the warmest period during 950-1250 CE (MCA) and coldest period during 1450-1950 CE (LIA) are defined by calculating the annual temperature anomaly over the NH (north of 30 • N) with respect to the 1000-1850 mean (the longest common period in the simulations) and lying above and below the mean for the MCA and LIA, respectively. Second, given the evidence for Atlantic southward/northward shifts of the ITCZ related to altered meridional sea surface temperature gradients between the tropical North and South Atlantic, we also verify that the periods identified with the first criterion correspond to periods when the surface temperature difference between the boxes (5-20 • N) and (20-5 • S) in the Atlantic were negative (positive) for LIA (MCA). We then verify that both criteria coincide. For example, for the LIA, the period with cold NH temperature anomalies coincide with temperature anomalies in the North Atlantic box colder than that its South Atlantic box counterpart (negative gradient, not shown).The MCA and LIA periods identified in each model are shown in Table 1. Note that in general the periods are of the order of 80-110 years long: shorter than the more gen-  Figure 1b shows the Gaussian fit of the frequency distribution of NH temperatures of all the years defined as LIA years (red curve) and MCA years (blue curve). The difference between the two periods is statistically significant (bootstrap test, 5 % significance level). Even though the anomalies are rather weak during both periods, a comparison with the values from their respective control simulation (piControl) shows that both periods are also significantly different at the 5 % significance level from the long-term mean. In addition, Fig. 2 shows the maps of the annual mean temperature anomalies during LIA and MCA, as well as their difference, for the ensemble mean. Temperature anomalies in the models are largest over the NH and in particular over the North Atlantic domain. Importantly, however, the LIA and MCA periods identified in the models are not synchronous, as shown in Table 1.

Variables used
To identify the main differences in LM simulations of the SAMS, particularly during the LIA and MCA periods, we analyse monthly CMIP5/PMIP3 output for rain rate, and 850 and 200 hPa horizontal winds. We also analysed vertical wind (omega) at 500 hPa in order to evaluate regions of ascending motion. However, not all modelling groups have saved this variable, so this analysis was done with five models only, and a figure of the vertical wind changes is not included in the paper. All variables have been re-gridded using a simple linear interpolation to a common 2 × 2 • grid.
The oceanic ITCZ is identified following the method proposed by Frierson and Hwang (2012). They define the ITCZ location by the precipitation centroid, as the tropical latitude with the maximum precipitation, at all longitudes over the ocean. Following their method, the precipitation was first interpolated onto a grid of 0.1 • to allow the precipitation centroid to vary. We explicitly do not consider the precipitation maxima over continents due to known problems in the correct definition of the ITCZ (e.g. see Laderrach and Raible, 2013;Nicholson, 2009).
The next section examines the performance of the models and whether they simulate a stronger SAMS during the LIA, in comparison to the MCA, as suggested by precipitation proxies and previous modelling experiments. In addition, since the SAMS is a dominant feature of the South American climate during austral summer (e.g. Vera et al., 2006), we focused on its mature phase, the December-January-February (DJF) season. Figure 3a shows the annual mean precipitation difference between the LIA and MCA periods. Blue and red curves correspond to the annual mean position of the oceanic ITCZ during LIA and MCA periods, respectively. The ensemble mean shows that the precipitation differences are small and statistically significant only in some regions (bootstrap test, p < 0.05). There is more precipitation during the LIA compared with the MCA in northeastern Brazil and across the tropical Atlantic, which are regions directly affected by the ITCZ position in the current climate. The mean position of the ITCZ between the two periods does not show any significant shifts (see Fig. 3b), but a small southward shift in the Atlantic during the LIA is found, in accordance with the precipitation signal. Individually, models do show that during the LIA the ITCZ was shifted further southward at some lon- gitudes (Pacific and Atlantic Oceans) when compared with the MCA (not shown). Figure 4a shows climatological precipitation and 850 hPa atmospheric circulation over the SAMS region during austral summer. In general, models are able to reproduce the main summer circulation and precipitation characteristics over South America observed in present-day climate. A narrow oceanic ITCZ, a broad area of maxima rainfall over the continent (SAMS), and a southeast-northwest-oriented South Atlantic Convergence Zone (SACZ) are observed in LM simulations, consistent with present-day observations (e.g. Garreaud et al., 2009). However, some models exhibit a double ITCZ over the eastern Pacific. This bias has been previously identified in CMIP3 and CMIP5 simulations, especially during austral summer and autumn seasons (Hirota and Takayabu, 2013;Sierra et al., 2015). Despite the limitations of model resolution, austral summer lower-tropospheric circulation simulated by the ensemble mean reproduces a cyclonic circulation over southeastern Bolivia (aka "Chaco low") as seen in observations and its associated northerly low-level jet, which is channelled by the Andes topography, transporting moisture to southern South America (Marengo et al., 2004).

Precipitation
When LIA and MCA composites for DJF are compared (Fig. 4b), the models exhibit an increased easterly flow at approximately 5 • S over the Atlantic and a weaker northerly low-level jet north of the Chaco low region, consistent with less precipitation over the SACZ during the LIA. Models also simulate less summer SAMS precipitation during LIA over the Amazon and the SACZ but more in the Nordeste. When analysing the associated vertical motion in models for which this variable is available, during the LIA, compared to MCA, there is stronger ascending motion in most of the SAMS domain (not shown). This precipitation pattern is in opposition to rainfall reconstructions over the western Amazon, the SACZ, and the Nordeste (e.g. Vuille et al., 2012;Novello et al., 2012Novello et al., , 2016Apaestegui et al., 2014). By contrast, with regard to annual mean simulations (Fig. 3), most models show a southward migration of sections of the Atlantic ITCZ (not very visible in the ensemble mean) and enhanced precipitation over the SAMS domain during the LIA, particularly over the eastern and southern Amazon, in agreement with palaeoclimatological records for this period. This indicates that the LM simulations are not able to reproduce the expected changes of the austral summer Atlantic ITCZ location and SAMS rainfall during LIA and MCA periods. The positive changes in the annual mean seen in Fig. 3 are due to the spring and autumn transition seasons.

Bolivian high and subtropical jet
The well-documented southward migration of the Hadley cell and its rising centre from 10 • N in JJA to 10 • S in DJF is only a part of the monsoon rainfall seasonal migration over the Americas, which reaches a more southward location in austral summer (Dima and Wallace, 2003). Furthermore, this wide area of continental convection, although related to local convergence zones, is not only a result of the shift of the ITCZ into subtropical latitudes. The establishment of the Bolivian high, the characteristic monsoon upper-level anticyclone located over the central Andes during austral summer, and the position and strength of the SH subtropical jet (SHSJ) in South America are also related to this monsoonal convective activity (Lenters and Cook, 1997;Garreaud et al., 2003;Yin et al., 2014).
To identify changes in the Bolivian high during the LM, we analyse the austral summer upper-troposphere circulation during the LIA and MCA (Fig. 5). Results indicate a stronger and more southeastward location of the SAMS anticyclone during the LIA. This strengthening of the Bolivian high is consistent with a stronger SAMS circulation. The southward shift of this upper-level anticyclone is related to an enhanced summer easterly flow over the central Andes, as suggested by previous studies (Lenters and Cook, 1999), and in turn would favour moisture transport and rainfall over the region (Garreaud et al., 2003). Moreover, the upper-tropospheric wind anomalies strikingly resemble the South American dipole (e.g. Robertson and Mechoso, 2000), a primary mode of variability over this region. An anticyclonic anomaly is associated with a diffuse SACZ, enhancing moisture convergence and precipitation on its southwestern flank (i.e. leading to a poleward shift in the location of the SACZ). Again, model simulations do not show this enhanced austral summer rainfall in the Amazon and central Andes during the LIA; they feature only marginally more precipitation to the southwest (Fig. 3).
On the other hand, recent studies have identified that the strength and location of the SHSJ, which corresponds to the southward extent of the Hadley cell, is a key factor for triggering convection during the dry-to-wet season transition in the Amazon (Yin et al., 2014). Particularly, when the SHSJ is weaker and/or reaches a more equatorward location, it promotes the incursion of synoptic disturbances to subtropical South America (e.g. Garreaud, 2000), enhancing lowertroposphere convergence and triggering the wet-season onset over the region (e.g. . To identify simulated changes of the SHSJ during the LIA and the MCA, Fig. 6 shows the 30 m s −1 isotach of the climatological September-November 200 hPa zonal wind as well as the difference between LIA and MCA periods. In general, the ensemble mean does not exhibit significant changes in the SHSJ location over South America during either period, as also indicated by Fig. 5b; however, the models simulate a weaker SHSJ during the LIA, not only in austral spring but also for the annual mean and summer seasons (not shown). This weaker SHSJ, particularly during austral spring (i.e. the transition season from dry to wet conditions in the SAMS), would allow a stronger influence of cold air incursions to trigger SAMS convection and probably maintain a stronger monsoon during the LIA.

Discussion and conclusions
According to our analysis, LM simulations are able to identify circulation features coherent with a stronger SAMS during the LIA: (i) an enhancement of the rising motion in the SAMS domain in austral summer; (ii) a stronger monsoon-related upper-troposphere anticyclone; (iii) activation of the South American dipole, which results to a certain extent in a poleward shift in the SACZ; and (iv) a weaker spring SHSJ over South America. However, austral summer simulations do not exhibit the expected increase in precipitation in this region during this cold period, as suggested by proxy evidence, except over the Nordeste, where it is not expected based on proxy data . Furthermore, LM simulations only reproduce a slight, but insignificant, southward (northward) shift of the austral summer Atlantic ITCZ during the LIA (MCA), unlike results found in other modelling studies (Vellinga and Wu, 2004;Lee et al., 2011;Kageyama et al., 2013). This disagreement might be partially related to the fact that the above-mentioned modelling studies impose much stronger external forcing than the forcing used in the LM simulations. This meridional shift of the Atlantic ITCZ is commonly considered a key aspect to explain the changes in SAMS rainfall observed during these periods (e.g. Vuille et al., 2012).
Recent studies indicate that the new generation of models included in the CMIP5 still tend to perform poorly in simulating precipitation in South America, especially over the Amazon Basin and the Atlantic ITCZ (Yin et al., 2013;Siongco et al., 2014;Sierra et al., 2015). However, CMIP5 models have shown further improvement in simulating precipitation over the region, in comparison to the CMIP3 generation (Jones and Carvalho, 2013;Yin et al., 2013;Hirota and Takayabu, 2013).
What could bias the simulated austral summer SAMS rainfall response of the CMIP5 models during the past millennium? Recent studies indicate that CMIP5 simulations tend to overestimate rainfall over the Atlantic ITCZ (Yin et al., 2013) and exhibit either an east or west Atlantic bias, in association with overestimated rainfall along the African (Gulf of Guinea) or South American (Brazil) coast, respectively (Siongco et al., 2014). Such a misinterpretation of the local ITCZ has been shown to bias rainfall simulations in the core of the SAMS (Bombardi and Carvalho, 2011). A stronger Atlantic ITCZ, for example, may contribute to enhanced surface divergence over tropical South America, inducing drier conditions in the region (e.g. , as observed in CMIP5 historical simulations (Yin et al., 2013;Sierra et al., 2015). However, a stronger local ITCZ does not necessarily translate into reduced SAMS rainfall since moisture convergence in this region is mainly influenced by the SACZ (Vera et al., 2009). Thus, the weaker SACZ during the LIA simulated by these models (Fig. 3) could reduce moisture convergence and rainfall over the SAMS. Furthermore, positive feedbacks between land surface latent heat flux, rainfall, surface net radiation, and large-scale circulation are also found to contribute to the dry biases over the Amazon and SAMS in most of the CMIP5 historical simulations (Yin et al., 2013).
Another circulation feature related to SAMS rainfall is the intensity and location of the South Atlantic subtropical high. The eastward displacement of this anticyclone and its inwww.clim-past.net/12/1681/2016/ teraction with the SACZ provide favourable conditions for monsoon precipitation (Raia and Cavalcanti, 2008). Recent analysis of CMIP5 projections under different scenarios suggests that this surface anticyclone is likely to strengthen in association with globally warmer conditions . Thus, a detailed examination of the response of this subtropical high to LM forcing is necessary in order to provide further explanations for the inadequate CMIP5/PMIP3 simulations of the SAMS rainfall variability throughout the past millennium.
The previous generation of LM model simulations reproduced warmer temperatures during the MCA when compared with the LIA but generally underestimated the regional changes detected from available reconstructions or failed to simulate a synchronous response in accordance with these reconstructions (e.g. Gonzalez-Rouco et al., 2011). The latter has been mainly related to uncertainties in the forcing estimates, as well as reduced sensitivity to external perturbations, underestimated internal variability, or incorrect representation of important feedbacks in GCMs (e.g. Goosse et al., 2005;Braconnot et al., 2012). Some of these problems still persist in the PMIP3 LM simulations (PAGES 2k-PMIP3 group, 2015). Furthermore, a recent model simulation of the global monsoon during the LM, performed in a non-PMIP3 experiment, indicates that the NH summer monsoon responds more sensitively to GHG forcing than the SH monsoon rainfall, which appears to be more strongly influenced by solar and volcanic forcing (Liu et al., 2012;Colose et al., 2016;Novello et al., 2016). Hence, a stronger sensitivity of SAMS rainfall to LM forcing estimations and the inadequate response of current GCMs to such forcings may also bias the CMIP5/PMIP3 simulations of the summer SAMS rainfall during the past millennium. Therefore, weak temperature response seen in these models during the MCA (Figs. 1 and 2) could contribute to the inadequate changes of austral summer rainfall in South America between LIA and MCA ( Figs. 3 and 4).
This evaluation of the SAMS throughout the past 1000 years in the latest generation of LM simulations confirms previous findings regarding the ability of the current generation of GCMs to reproduce large-scale circulation features in South America and their lack of an adequate representation of precipitation over the region. However, the weak or absent temperature and precipitation response to the imposed forcing in climate models provides a formidable challenge for proxy-model comparisons. To better compare and eventually reconcile model reconstructions with proxy evidence will require a more detailed analysis of precipitationgenerating mechanisms in climate models. Our results indicate that the CMIP5/PMIP3 models quite accurately reproduce changes in the large-scale circulation that in turn are consistent with proxy evidence of precipitation changes over the past millennium. These changes, however, do not translate into corresponding precipitation changes. This implies that the models may lack relevant feedbacks or that precipita-tion in the models may be too dependent on the microphysics and convective parameterization schemes but not sufficiently sensitive to large-scale circulation mechanisms. On the proxy side, despite recent multi-proxy reconstructions of temperature in South America (e.g. PAGES 2k Consortium, 2013;, Neukom et al., 2011), a stronger effort to not only reconstruct surface climate at individual locations but also focus on reconstructions of modes of variability or entire climate components such as the SAMS, which implicitly include circulation changes, is needed. Proxies such as pollen or stable hydrogen and oxygen isotopes from lakes, speleothems, and ice cores have shown potential to record larger-scale climate signals and changes in the tropical hydrological cycle over South America (Vuille and Werner, 2005;Vimeux et al., 2009;Bird et al., 2011;Vuille et al., 2012;Ledru et al., 2013;Flantua et al., 2016;Hurley et al., 2015). Multi-proxy reconstructions from such networks, which implicitly incorporate remote and large-scale circulation aspects, may therefore provide a better tool to assess the performance of climate models than reconstructions that are based solely on local precipitation estimates.