Introduction

The Indian summer monsoon (ISM) is a large-scale coupled land-ocean-atmosphere phenomenon that plays a dominant role in transporting water vapour to the Indian subcontinent1. Despite the strong bearing and societal relevance of the monsoon, reliable projections of ISM rainfall remain a challenge2, mainly because of the complex dynamics and the myriad of interactions with other tropical and extra-tropical processes across various timescales3.

On a year to year basis, the ISM exhibits considerable complex and time-varying interactions with the El Niño/Southern Oscillation (ENSO)4, the Indian Ocean Dipole5 and Tibetan snow cover6. Additionally, intraseasonal active or break events on short timescales, presumably associated with internal dynamics of the monsoon system itself, have large impacts on the interannual variability of monsoon rainfall7,8. Example includes the famous break of July 2002, where rainfall was reduced by 50%9. These large year-to-year variations in rainfall associated with internal and external interactions pose a considerable challenge for Indian monsoon forecasting10 and for attributing monsoon variability to large-scale forcings. This is particularly true for multi-decadal and longer time scales given that observational data only covers few cycles and paleoclimate records span longer time scales which fill this gap.

There is currently a considerable debate on the impact of different forcing agents versus the role of internal variability in determining past and future long-term trends of the ISM11. Global climate model experiments tend to suggest an intensification of the ISM in response to global warming partly because of the larger water-holding capacity of a warmer atmosphere, resulting in an estimated 10% increase by the end of the 21st century12; yet, there are large inter-model discrepancies (e.g., ref. 10). This is in striking contrast with the recent observational record, which indicates an overall decrease of the ISM rainfall during the period 1950–201213,14,15, followed by a weak revival between 2012–201416. Decadal-scale internal climate variability associated with ENSO and the Atlantic Multi-decadal Oscillation (AMO) can also play an important role in explaining the observed variations in NH monsoon precipitation during the last century17 and possibly beyond18,19. Thus a number of drivers, can contribute to modulating ISM rainfall, including internal variability of the climate system as well as interactions with natural and anthropogenic factors, either of global (e.g., greenhouse gases) or regional (e.g., anthropogenic aerosols and land-use) nature11,20,21,22.

More integrated studies of the wider Northern Hemisphere (NH) summer monsoon suggest coherent changes of monsoon precipitation and circulation at hemispheric scales23. The relatively coherent response of the various monsoons organised into the NH monsoon system to natural and anthropogenic forcing suggests that monsoon rainfall is sensitive to common hemispheric and/or planetary scale controls. Implicit in this view is the hypothesis that large-scale forcings may override regional-scale factors both natural and anthropogenic in determining future NH monsoon rainfall trends23,24. Conversely, several studies also suggest a varying response of the individual monsoon systems to the same forcing due to internal feedback dynamics25. Therefore, a better understanding of the large-scale controls of past decadal-scale variability of the monsoon is key to achieve more robust projections of its future changes. However, observations are spatially and temporally limited, mostly to the 20th century, which only allow identifying few long-term fluctuations17. Further extending these records back in time using paleoclimate reconstructions is essential to reduce current uncertainties, to identify its local and global scale forcing mechanisms, and to quantify the background monsoon variability on top of which shorter timescale fluctuations as well as longer-term anthropogenic-induced changes are superimposed24. Here we present a new record of ISM rainfall and SST variability from the Northern Indian Ocean spanning the last 2,000 years at multi-decadal resolution and investigate the large-scale forcing mechanisms of ISM rainfall.

Salinity Variations in the Bay of Bengal

The Bay of Bengal (hereafter BoB; see Fig. 1) is a natural laboratory for reconstructing past ISM rainfall variations for two main reasons: Firstly, the BoB receives an estimated annual freshwater discharge of 2950 km3 from major (Ganges-Bramhaputra and Irrawaddy) and minor (Mahanadi, Krishna and Godavari) rivers fed by the ISM monsoon precipitation and glacier melt water25. Thus the catchment of these rivers cover a large swath of the summer monsoon rainfall regions of North, Northeast and Central India. Secondly, the BoB directly receives heavy precipitation particularly in the northern parts. Precipitation plus runoff exceeds evaporation on an annual cycle. Fresh water from these sources is mixed into the upper layers of the BoB interacting with the surface circulation. After the summer monsoon some of the strongest lowering of seasonal salinity occurs along the Coromandel Coast in the northern and western BoB (Fig. 1). We exploited this prominent low salinity cap to reconstruct ISM rainfall and SST records from the Sediment Core MD161/17. The core site receives sediments from Mahanadi, Krishna and Godavari, and its tributaries, draining the northern portion of the Indian plateau, leading to enhanced sedimentation and expanded sediment sequence in the Krishna-Godavari Basin26,27.

Figure 1
figure 1

Location of MD161/17 Core and July through September salinity pattern in the northern Indian Ocean. Note MD161/17 Core was retrieved from the Krishna-Godavari Basin where a distinct low salinity prevails due the heavy over head precipitation and river discharge of Krishna, Godavari Rivers during Indian summer monsoon.

Full size image

Proxy Records of Indian Monsoon Rainfall

In this study, we use Mg/Ca (SST) and δ18O calcite (δ18Oc) from the calcite shells of planktonic foraminifera Globigerinoides ruber (hereafter G. ruber) to derive δ18Owater (δ18Ow) to reconstruct ISM rainfall variability in the BoB at muti-decadal to centennial time scales (see suppl. material for further details). δ18Ow exploits the higher O16/O18 isotopic ratios of precipitation associated with rainfall as opposed to sea water which has lower O16/O18 ratios after appropriate corrections. On land, δ18O of precipitation exhibits spatial variability and shows a negative relationship with precipitation which could vary over the monsoon season24. However, the effects of such spatial and temporal variability are minimised in our study as the BoB receives rainfall from a large area of the Indian monsoon fed region and salinity anomalies persist over the monsoon and post monsoon periods. Also, the δ18Ow and salinity show a strong positive relationship (r2 = 0.89 with a slope of 0.15) in the BoB28. Therefore, δ18Ow, at least in a qualitative sense, reflects salinity variations related to monsoon rainfall at the core site. Calcite tests of G. ruber are used for reconstructions because G. ruber fluxes from sediment trap time series from the BoB suggest that they record average surface conditions in the BoB29. Previous studies have successfully used G. ruber to reconstruct ISM rainfall variability at millennial scale in the BoB since the last glacial period30,31. High sedimentation rates (0.2 to 0.6 cmyr−1) in Core MD161/17 allow us to reconstruct ISM rainfall variability at multi-decadal temporal resolution for the first time with robust chronology derived from 11 radiocarbon dates of mixed planktonic foraminifera species (Fig. S1).

Results

The δ18Ow record shows a large variability (1.1‰) over the last 2000 years, with high and low δ18Ow values representing decrease and increase of ISM rainfall respectively (Fig. S2). Figure 2a presents δ18Ow as positive and negative anomalies over the mean of the record. Major portions of the Roman Warm Period (1950 to 1550 years before present (BP); RWP) are associated with negative δ18Ow anomalies, suggesting higher ISM rainfall. During Dark Age Cold Period (DACP) (1550 to 1300 year BP), high positive δ18Ow anomalies are indicative of weaker ISM rainfall. The Medieval Warm Period (MWP) from 1200 to 800 years BP is characterized by negative δ18Ow anomalies compared to the preceding period (i.e. 2000 to 1300 years BP), which reveals that BoB experienced highest ISM rainfall and strong freshening during this MWP over the last 2000 years. In contrast, during the Little Ice Age (LIA), positive δ18Ow anomalies indicate a decrease in ISM rainfall from 700 to 200 years BP. Our ISM rainfall record from the BoB shows similarity with the speleothem record from the Sahiya Cave (Fig. 2b) with statistically significant positive correlation at zero-lag (Fig. S3) reflecting continental rainfall over the Northern Himalayan foothills19 (Fig. 1). This confirms that our reconstructed variability in ISM rainfall is broadly representative of a large swath of Northern & Central India including the Indo-Gangetic plains.

Figure 2
figure 2

Variability of Indian summer monsoon rainfall over last two millennia; (a) δ18Ow anomaly derived from δ18Ow from δ18Oc of Globigerinoides ruber from MD161/17 Core; (b) δ18Ow anomaly of Sahiya Cave obtained from the δ18O of speleothem19; (c) sea surface temperature anomaly from MD161/17 core obatined from the Mg/Ca ratios in Globigerinoides ruber; (d) Northern Hemisphere Air Temperature anomaly32, (e) Multi-decadal variability of Atlantic multi-decadal oscillations51.

Full size image

Reconstructed SST from foraminiferal Mg/Ca ratios over the last 2000 years varies between 26.7°–30.3 °C (Fig. S2). SST anomalies are more than 1 °C warmer during LIA and DACP and about 1 °C cooler during the MWP as compared to the mean SST (Fig. 2c). Figure 3 shows the strong negative correlation between reconstructed ISM rainfall (represented by δ18Ow) and Mg/Ca based SST (r = −0.75; p < 0.001; n = 84). Cross correlation analyses of these two time series records also showed a significant negative correlation with zero-lag suggesting a strong temporal coupling between BoB SST and ISM rainfall over the last 2000 years at this core location (Fig. S4 and Table S2). In Fig. 2d, we compare the anomalies of ISM rainfall and BoB SST with the Northern Hemisphere surface air temperature (NHT) anomaly data of Mann and Jones32 for the last 1500 years. The composite NHT is derived from Atlantic and European regions but also include data from central and eastern Eurasian regions30. Increased ISM rainfall and cooler BoB SSTs occur during periods of warmer NHT and vice versa (Fig. 2d). Also increased and decreased ISM rainfall during MWP and LIA corresponds with increase and decrease of AMO respectively (Fig. 2e). Kernel cross-correlation (red line) analysis of δ18Ow anomaly (proxy of ISM rainfall) and NHT show significant positive correlation at zero year lag (Fig. S5). In summary, our reconstruction reveals the positive temporal relationship between ISM rainfall and Eurasian temperature, as well as the inverse relationship between ISM rainfall and BoB SST. In the following sections, we further explore these links by using recent (last 100 years) observations which, despite being temporally limited, provides extensive spatial coverage useful to understand the nature and mechanisms underpinning decadal ISM rainfall variability in our reconstructions of the last 2000 years.

Figure 3
figure 3

Relationship between the sea surface temperature and Indian Summer Monsoon Rainfall varibility over last two millennia. Here the depleted and enriched δ18Ow represent more and less ISM variability respectively. Note the significant negative correlation between the SST and ISM.

Full size image

Discussion

The Indian summer monsoon is a fully coupled ocean–land–atmosphere phenomenon whose fundamental driving mechanism can be identified in the differential heating between the land and the ocean to the south. During boreal spring, the land region over Asia and SE Asia warms faster than the ocean due to the lower heat capacity, setting up a low sea level pressure anomaly over Northern India and the Middle East. The resulting southward pressure gradient drives a cross-equatorial flow at the surface and a return flow aloft, forming a thermally-direct meridional cell8. A number of studies has found the Himalayas and Tibetan Plateau to also play an important thermal and mechanical role contributing to anchoring and intensifying the monsoon circulation by, respectively, heating the mid troposphere as well as by preventing cold air of mid-latitude origin to be advected to the south33,34. The large-scale inter-hemispheric thermal contrast that develops in late-spring and early summer leads to the northward seasonal migration of the Inter-Tropical Convergence Zone10. A large proportion of monsoon rainfall falls over northeastern peninsular India and in other orographic hotspots such as the Western Ghats and Arakan Range of Myanmar35. Once the monsoon is established over the Indian subcontinent, increased (land and ocean) evaporation and cloud cover lead to summertime cooling of the Indian Peninsula as well as of the north-equatorial Indian Ocean, with additional contribution from increased upwelling driven by the intense monsoon winds over its western sector1. As the Himalayas acts as a physical barrier for the southerly moisture-laden flow, precipitation and associated surface cooling are confined largely to the Indian Subcontinent, while the Tibetan Plateau is mainly dry and warm35. This maintains the meridional temperature gradient sustaining monsoon circulation and rainfall during the summer monsoon season36. Our reconstructions shown in Fig. 2a,c capture two key characteristics of the fundamental mechanics underpinning the annual cycle of the ISM1: salinity variations in the BoB reflecting monsoon rainfall and freshwater flux and2 lowering of the surface SSTs of the northern Indian Ocean during the monsoon period representing evaporative cooling associated with intensified monsoon winds and cloud cover. The inverse relationship (Fig. 3) and the coherent response (Fig. 2a,c) of these two reconstructed elements provide assurance that the changes we see are related to ISM variability. Furthermore, the close correspondence between these reconstructions and the NHT records (Fig. 2d,e) support the causal link between NH temperature and ISM rainfall via seasonal changes in the meridional surface/tropospheric thermal gradients23,24.

Drivers of Multi-decadal ISM Rainfall Variability

We further explore the link between large-scale temperature changes and ISM variability at multi-decadal scale by means of observational data during 1901–2012. To make a direct comparison with the multi-decadal proxy records we remove interannual and short time scale variability in the observational data by low-pass filtering. This isolates the decadal or longer time-scale components of the signal from higher-frequency variability which is known to be largely influenced by internal monsoon dynamics10. The purpose of using observational data here is twofold1: to evaluate whether seasonal features such as cooling of the BoB related to monsoon circulation are represented in decadal time scales and2 to identify large-scale forcing mechanisms of All-India Monsoon Rainfall (AIMR) fluctuations at these time scales. Figure 4 presents the spatial map of June-September (JJAS) surface temperature anomalies associated with decadal or longer variations in the AIMR. At the surface, excess AIMR is associated with significant positive temperature anomalies over Eurasia and in particular eastern Asia and the Tibetan Plateau (Fig. 4a) and negative temperature anomalies over the Indian Subcontinent and the Tropical Indian Ocean (Fig. 4b). The elevated heating over the Tibetan Plateau, concurrently with slight cooling of the north equatorial Indian Ocean, has been identified as a dominant factor driving the annual cycle and interannual variability of the monsoon system7,37. Further afield, a strong positive association between ISM rainfall and SST anomalies in the North Atlantic Ocean and surface temperature over western Eurasia is evident in Fig. 4b 38. Large basin-wide sea surface temperature anomalies are also identifiable across the Pacific, with cooling along the equatorial central-eastern Pacific, and an extensive warming in the extra tropics stretching from Japan to the western coast of the US. The tropical anomalies are somewhat reminiscent of a La Niña-like cooling in the eastern equatorial Pacific. However, the pattern also features substantial differences in two important and fundamental characteristics of the actual La Niña sea surface temperatures anomalies or of its longer-lived manifestation, the Pacific Decadal Oscillation (PDO)39. First, it does not display zonal opposite-sign SST anomalies over the western equatorial Pacific, and second, the anomalies do not extend over the northeastern Pacific. We will return to this point later. Overall, this corroborates the assembled paleo-reconstructions presented in Fig. 2 and allows us to place ISM rainfall variability over the last 2000 years in the wider context of mechanisms associated with NH climate variability.

Figure 4
figure 4

Regressions of observed low-pass filtered June-September (JJAS) mean (a) surface temperature (°C) and (b) sea-surface temperature (°C) on observed JJAS All-India summer monsoon rainfall during 1901-2012. All data were linearly detrended prior to the analysis. The analysis focuses on decadal or longer scale variability, isolated by applying a low-pass Lanczos filter with a cut-off frequency of 10 years to JJAS-mean data. The cross-hatching marks regions where the correlation exceeds the 90% (grey) and 95% (black) confidence levels, estimated by using a Monte Carlo approach with 1,000 random samples of the time series68. Land temperature data are from the Climate Research Unit (CRU), University of East Anglia, United Kingdom, dataset (CRU-TS3.21) at 0.5° × 0.5° resolution69, while sea surface temperature are from the Hadley Centre Sea Ice and Sea Surface Temperature data set (HadISST) dataset version 1 at 0.5° × 0.5° resolution70. Rainfall data are taken from the homogeneous rainfall data set of 306 raingauges in India, developed by the Indian Institute of Tropical Meteorology71.

Full size image

The AMO has been suggested to have played a major role in modulating hemispheric-scale temperature gradients at multi-decadal time scales during the last 2000 years40,41 via its link with Atlantic Thermohaline Circulation (THC) variations and associated inter-hemispheric heat transport fluctuations42,43. The AMO has also been suggested to play a major role in modulating the 20th century multi-decadal variations of Indian and African (Sahel) summer rainfall, as well as the tropical Atlantic atmospheric circulation (e.g., refs. 17,44). Figure 5 shows regression/correlation patterns of JJAS precipitation (Fig. 5a), three-dimensional circulation (Fig. 5b), and upper-tropospheric temperature anomalies (Fig. 5c) on the AMO. A positive AMO phase is associated with overall increased precipitation over the Indian subcontinent, with a large positive anomaly over central India (peak values statistically significant at the 95% level) and reduced precipitation to the south and over the central Himalayas (Fig. 5a). This pattern is associated with stronger low-tropospheric southwesterly winds over the Arabian Sea bringing more moisture toward India (Fig. 5b). The mid- and upper-tropospheric heating anomaly (Fig. 5c) features significant warming from the Atlantic across northern Africa to the Middle East and northern India, which intensifies the meridional thermal contrast between the continent and the Indian Ocean- a pattern consistent with several other modelling studies45,46. The regional circulation anomalies associated with the southwesterlies in the northern Indian Ocean are part of a hemispheric-wide response pattern to AMO variability featuring enhanced easterly trade winds across the northern equatorial Pacific and enhanced Walker cell north of the equator, leading to larger ascent (precipitation) over Indo-china and the Maritime Continent and subsequent changes in the meridional direction (regional Hadley circulation). These patterns, including the Indo-Pacific sea surface temperature anomalies, bear remarkable similarity to those described in previous studies47 and related to variations in the Atlantic oceanic circulation. Figure 5b also displays an interesting wave pattern in the ascending and descending air masses in the northern hemisphere extratropics which some have indicated as responsible for the propagation of the AMO signal from the Atlantic across Eurasia (e.g., ref. 48). This analysis, albeit not exhaustive, is suggestive of an AMO-ISM teleconnection with signatures of both an atmospheric pathway via modulation of the upper-tropospheric temperature across Eurasia (e.g.45,48,49,50), and one that also involves atmosphere-ocean feedbacks in the Indo-Pacific region (e.g.47,51). Overall, this is indicative of an AMO imprint on the global interhemispheric temperature gradient with a positive phase of the AMO associated with a warmer NH than the SH (Fig. 4b), with ensuing modulation of the NH monsoon systems. We explore this further using palaeoclimate and observational time series.

Figure 5
figure 5

Regressions of observed low-pass filtered June-September (JJAS) mean (a) precipitation (mm day-1), (b) 500-hPa vertical velocity (hPa day-1) and 850-hPa winds (m s-1), where blue and red areas represent ascending and descending air masses respectively, (c) mean upper-tropospheric (200 – 500 hPa) air temperature (°C) on the AMO index shown in Fig. 6 during 1901-2012. All data were linearly detrended prior to the analysis. Precipitation data are from the CRU-TS3.21 dataset, while winds and air temperature are from NOAA-CIRES 20th Century Reanalysis dataset. Anomalies are given per °C change of the AMO. The cross-hatching marks regions where the correlation exceeds the 90% (grey) and 95% (black) confidence levels, estimated by using a Monte Carlo approach with 1,000 random samples of the time series.

Full size image

The time series displayed in Fig. 6a–c show that, in general, decadal-scale wet periods over India (e.g., 1926–1965) are in phase with a positive AMO (warmer North Atlantic SSTs) and warmer Asian surface temperature anomalies; conversely, dry periods of the ISM (i.e., 1901–1926 and 1965–1995) occur simultaneously with a negative phase of the AMO and cooler Asian surface temperature anomalies. Such forcing mechanism and distal teleconnections appear to be at play also during the last 2000 years based on our paleo-monsoon reconstructions (Fig. 2) revealing a coherent response between AMO and ISM, further substantiating previous studies based on sparse monsoon records24,52. On millennial to century time scales, ISM strength has been correlated with Dansgaard–Oeschger cycles in Greenland ice cores, where weak monsoon events are associated with cold events in the NH during the last glaciation and the Holocene52,53,54,55,56. Such weakening of the ISM has been related to Northern hemispheric cooling resulting from ice-sheet instability inducing changes in the Atlantic THC32,57. Therefore, a strong relationship between NHT and Indian monsoon variability appears to be a persistent and coherent feature of the past and observed climate at decadal or longer time scales.

Figure 6
figure 6

Time series of observed low-pass filtered (a) Atlantic Multidecadal Oscillation (AMO) index (°C), (b) June-September (JJAS) All-India summer monsoon rainfall anomalies (mm day-1), (c) JJAS surface temperature anomalies (°C) over Eurasia (20° - 140°E, 30° - 75°N), (d) December-February Niño3.4 index (°C, right y-axis), and (e) standardised annual mean Pacific Decadal Oscillation index during 1901–2012. All data were linearly detrended prior to the analysis. The analysis focuses on decadal or longer scale variability, isolated by applying a low-pass Lanczos filter with a cut-off frequency of 10 years to the detrended data. Land temperature data are from the Climate Research Unit (CRU), University of East Anglia, United Kingdom, dataset (CRU-TS3.21) at 0.5° × 0.5° resolution68, while sea surface temperature are from the Hadley Centre Sea Ice and Sea Surface Temperature data set (HadISST) dataset version 1 at 0.5° × 0.5° resolution69. Rainfall data are taken from the homogeneous rainfall data set of 306 rain gauges in India, developed by the Indian Institute of Tropical Meteorology70. The AMO index is defined as the area average over the North Atlantic (80°W - 0°E, 0° - 60°N) of annual mean sea surface temperatures17. The Niño3.4 index is calculated by averaging the time series of sea surface temperature (SST) anomalies over (5°S-5°N, 170°-120°W), with SST coming from the HadISST dataset as in Fig. 4. The unfiltered (yearly) values of the Niño3.4 index are also plotted in (a) as the continuous black line (left y-axis). The PDO index (obtained from http://research.jisao.washington.edu/pdo/PDO.latest) is derived as the leading principal component of monthly sea surface temperature anomalies in the North Pacific Ocean, poleward of 20°N. The monthly mean global average SST anomalies are removed to separate this pattern of variability from any global warming signal that may be present in the data.

Full size image

In contrast, the association of ISM rainfall with ENSO/Pacific Decadal Oscillation (PDO) is poor in both observational data (Fig. 6d,e) and palaeoclimate records on multi-decadal time scales (Figs. S7 and S8). This is in contrast with one recent study of the last 500 years, which revealed that AMO, ISM rainfall and ENSO exhibit a common 50–80 year variability suggesting that this mode is an integral part of global multi-decadal oscillations arising from large scale coupled ocean-atmosphere-land interactions58. Collectively, such discordance suggests that ENSO forcing on ISM rainfall is time varying and complex as suggested by several recent studies59,60. This is perhaps because ENSO forcing mechanism operates mainly on higher frequencies on which internal dynamics of monsoon itself dominates and hence not always detected in longer time scales. High frequency variability in rainfall associated with the ENSO4, the Indian Ocean Dipole5 and Tibetan snow cover6 are relatively small, the interannual standard deviation of summer-mean rainfall being around 10% of the climatological mean10. Observational data also suggest that the decadal pattern of ocean-atmospheric oscillations of the Pacific Ocean do not conform with the typical ENSO or PDO modes (e.g., Figs. 4b and 5b). Further investigations are required to understand fully the nature of this oscillation pattern.

In summary, we suggest that multi-decadal ISM rainfall variability during the last 2000 years was modulated by AMO related fluctuations in NH temperatures whereby warmer (cooler) Asian landmasses along with a cooler (warmer) Indian Ocean set up stronger (weaker) surface and tropospheric meridional temperature gradients strengthening (weakening) monsoon circulation and rainfall. Therefore, the reconstructed temporal relationships between ISM rainfall, BoB SST and NHT over the last 2000 years allow to extend and strengthen the findings based on the rather limited observational record for the 20th century shedding light on distal teleconnections underpinning multi-decadal ISM variations. It is interesting to note that Fig. 6 suggests a weakening of the AMO-ISMR relationship from the mid-1990s, when the AMO entered a positive phase while the ISM rainfall remained below average (ref. 17 Fig. 6). This recent decoupling is discussed further in the context of anthropogenic forcing.

Implications for ISM Variations Under Anthropogenic Forcing

Our results suggesting a link between the AMO and the associated NH differential warming and ISM rainfall variability at multi-decadal timescales during the last 2 millennia have important implications for the current debate on the future evolution of ISM rainfall in response to global warming. Our study strongly supports the notion that multi-decadal variations of ISMR are sensitive to extra tropical forcing that occurs beyond regional scales. It is also on these timescales that anthropogenic forcing can interact with internal variability of the monsoon. Broadly, AMO forcing on the monsoon acts by changing inter-hemispheric thermal gradients (Fig. 4). NH warming associated with positive phases of this oscillation enhances meridional pressure gradients that drive cross-equatorial flows and intensify NH monsoon circulation and rainfall17 (Fig. 5). Global warming trends since the early seventies have been asymmetric across the hemispheres, with the NH warming faster than the SH23,61. The interaction between this global warming pattern and positive phases of the AMO would further accentuate the inter-hemispheric thermal contrast, reinforcing NH temperature positive anomalies and intensifying NH monsoon circulation and rainfall. The post mid-1990s decoupling of ISM from AMO (Fig. 6) has been widely attributed to the influence of aerosols and land use changes8,11,19. In view of our results, such regional influences dampening ISMR are likely to be temporary depending partly on future levels of these aerosol emission and land use change. On longer time scales, NH temperature anomalies will have an overriding influence on enhancing monsoon. In this regard, we note that ISM rainfall registered a recovery since 2002 with above average rainfall after 20122,23,43,62. Finally, we suggest that interaction between AMO and global warming is likely to be a crucial factor in ISM rainfall trajectory into the future. Therefore, large-scale interactions involving extratropical factors in response greenhouse gas warming cannot be ignored in future ISM rainfall prediction.

Methods Summary

Core MD161/17 was collected at a water depth of 790 m from the KGB. Age model was based on linear interpolation between 11 AMS radiocarbon dates which covers a time span of 2 kyr (Supplementary Table 2). Globigerinoides ruber (150 to 250 µm) was analyzed for δ18O by reacting with 100% orthoposphoric acid at a constant 75 °C in a Kiel Carbonate III preparation device and evolved CO2 gas was analyzed by using Thermo Electron delta plus advantage stable isotope ratio mass spectrometer. Cleaning protocol of G. ruber for Mg/Ca analyses were followed by Barker et al.63, in addition a reductive cleaning step was used according to Martin and Lea64. Mg/Ca ratios were measured by using an axial viewing Varian 720 ICP-OES (for details see the supplementary material). Mg/Ca values were then used to estimate SSTs using the equation Mg/Ca = 0.449exp(0.09*T)65. δ18Osw were computed by applying the following equation of Bemis et al.66: δ18Osw = 0.27 ((T − 16.5 + 4.8*d18O)/4.8). The derived δ18Osw estimates were corrected for continental ice volume using Shackleton’s data set67 and presented here as δ18Ow. SST and δ18Ow anomalies are computed by subtracting the mean values of all data points within two millennia from the SST and δ18Ow values respectively.