Impact of Climate Change on Past Indian Monsoon and Circulation: A Perspective Based on Radiogenic and Trace Metal Geochemistry

: The Indian summer monsoon (ISM), one of the dramatic illustrations of seasonal hydrological variability in the climate system, affects billions of lives. The ISM dominantly controls the northern Indian Ocean sea-surface salinity, mostly in the Bay of Bengal and the Andaman Sea, by the Ganga-Brahmaputra-Meghna and Irrawaddy-Salween rivers outﬂow and direct rainfall. In the past decade, numerous studies have used radiogenic neodymium ( ε Nd ) isotopes of seawater to link Indian subcontinent erosion and the ensuing increase in discharge that results in changes in the north Indian Ocean sea surface. Here we synthesized the state of the ISM and ocean circulation using the neodymium and hafnium isotopes from north Indian Ocean deep-sea sediments. Our data suggest that the Bay of Bengal and north Indian Ocean sea-surface conditions were most likely modulated by changes in the ISM strength during the last glacial-interglacial cycle. These ﬁndings contrast to the hypothesis that suggests that the bottom water neodymium isotopes of the northern Indian Ocean were modulated by switching between two distant sources, namely North Atlantic Deep Water and Antarctic bottom water. Furthermore, the consistency between the neodymium and hafnium isotopes during the last glacial maximum and Holocene suggests a weak and dry ISM and strong and wet conditions, respectively. These data also indicate that the primary source of these isotopes was the Himalayas. Our results support the previously published paleo-proxy records, indicating weak and strong monsoons for the same periods. Moreover, our data further support the hypothesis that the northern Indian Ocean neodymium isotopes were decoupled from the global ocean neodymium budget due to the greater regional inﬂuence by the great Ganga-Brahmaputra-Meghna and Irrawaddy-Salween discharge draining the Indian subcontinent to the Bay of Bengal and the Andaman Sea.


Introduction
Monsoons represent one of the Earth's most dynamic interactions between the atmosphere, ocean, and continent [1][2][3]. The livelihoods of approximately one-half of the Earth's population depend on the monsoon's predictable annual return, particularly in West Africa, Central and East Asia, and Australia [4,5]. The Indian summer monsoon (ISM), part of the broader Asian monsoon system [6,7], is one of the Earth's most notable demonstrations of seasonal hydrological variability of the climate system. The asymmetric heating between the warmer Indian subcontinent and the cooler Indian Ocean during the summer (May-September) results in a steep pressure gradient, allowing the transfer of a large amount of moisture and causing severe precipitation [8,9]. In contrast to the ISM, the weaker winter monsoon brings little rainfall to the Indian subcontinent. The ISM predominantly controls the sea-surface salinity, especially in the Bay of Bengal and the Andaman Sea, by the direct rainfall and Ganga-Brahmaputra-Meghna and Irrawaddy-Salween rivers' outflow [10]. With the potential impacts of increasing greenhouse gases and the ensuing global temperature rise of 1.5-5.8 • C by 2100, most scientific investigations, including the 2019 Intergovernmental Panel for Climate Change 6th Assessment Report [11], predict that up to a quarter of the global mountain glacier mass could disappear by 2050 and up to a half or more could be lost by 2100 [12,13]. Furthermore, the Indian Himalayan glaciers could suffer similar or worse consequences as a result of the atmospheric warming [14] caused by increasing greenhouse gases [15], augmented by the absorption of heat by atmospheric black carbon [16]. If the 2019 Intergovernmental Panel for Climate Change predictions of global temperature rise are correct, the freshwater availability resulting from the melting of these glaciers and monsoonal rainfall becomes a significant concern for sustaining the livelihood of a billion people. In this context, past climate studies that refer to the present issues of global warming and glacial melt could shed light on the monsoon's future aspects for the Indian subcontinent.
Historical records of the ISM that use organic and inorganic proxies diverge in reconstructing the monsoon's strength due to complications inherent in paleo-proxies [17]. For example, Arbuszewski et al. [18] suggested that the seawater oxygen isotopes, used as a proxy for the past strength of the continental runoff, potentially overestimate surface hydrology due to the impact of high salinity on the oxygen isotopes in planktonic foraminifers [19]. In contrast to the inorganic proxies, paleo-proxy records using organic biomarkers may also suffer from complications arising from lateral transport (i.e., long-distance transport), bottom current reorganization, influence of mass-transport, etc., especially on the continental margins like the Bay of Bengal, Andaman Sea or Arabian Sea [20]. Radiogenic isotopes in marine sediments, particularly the neodymium (Nd), were also used to reconstruct past ISM strength through assessing erosion of the Indian subcontinent and circulation [21][22][23][24][25][26]. The isotopic distribution of Nd in the ocean reflects local source provenance (related to the intensity of continental weathering and orogenic and volcanic processes) and changes in paleo-circulation due to its short residence time (500-1000 years) [27][28][29][30][31]. The Hafnium (Hf) radiogenic isotopes behave similarly to the Nd isotopes but have a shorter residence time [32][33][34], and the elemental concentration of Hafnium is two to three times lower than that of the Nd [33]. Application of the ε Hf to document continental erosion or circulation is low in the global ocean [35] but even rarer in the Indian Ocean. Moreover, the combined use of Lu-Hf and Sm-Nd isotope systems, two radiogenic isotopes widely used to trace the long-term fractionation processes occurring in the Earth system because these isotopes are robust and little affected by alteration processes, also offer a unique perspective on silicate weathering when the detrital fraction is analyzed. When marine sediments are plotted on the "terrestrial array" (representing the correlation between Nd and Hf isotopes; [36]), it defines a large field in the Nd versus Hf isotopic plot with a wide range of Hf isotopes at any given Nd isotopic composition, compared to the tight array defined by both oceanic and continental rocks [37,38]. In short, marine sediments composed of clays and muds have more variable Nd isotopic values, and their Hf isotopes are systematically lower than those of the Fe-Mn crusts and nodules. The decoupling between Lu-Hf and Sm-Nd systems is closely linked to the so-called "zircon-effect": Hf is sequestered in erosionresistant zircons that remain on the continents or are deposited on the continental margin. Simultaneously, the rare-earth elements (REEs) are entrained in clays and particulates traveling to the deep ocean [39]. Therefore, the study of the detrital fraction of marine sediments provides additional information on sources and past silicate weathering intensity on continents. In any case, both Nd and Hf isotopic signatures are expressed in epsilon (ε) notation, as the deviation from the Chondritic Uniform Reservoir, by ε Nd or ε Hf = (R sample /R CHUR − 1) × 10 4 , where R represents 143 Nd/ 144 Nd and 176 Hf/ 177 Hf ratios. The deviation of Chondritic Uniform Reservoir values for 143 Nd/ 144 Nd and 176 Hf/ 177 Hf is 0.512638 [40] and 0.282785 [41], respectively.
We present combined seawater and detrital ε Nd and ε Hf data from one of the tropical northern Indian Ocean sites-758 of the Ocean Drilling Program (ODP)-to reconstruct the past strength of the ISM. Furthermore, we make an attempt, by integrating published radiogenic and oxygen isotopes from the northern Indian Ocean, to provide answers to the outstanding questions for understanding the ISM and ocean circulation. Our current contribution seeks to (1) provide an insight into the past 145,000 years, i.e., last glacial cycle (ε Nd and ε Hf ), and (2) offer perspectives for the last glacial maximum (LGM) to the present, due to the rapid proliferation of ε Nd data from this region during the previous decade. In so doing, we address one fundamental question: whether ε Nd and ε Hf isotopes can concomitantly be used to assess the strength of past erosion of the Indian subcontinent and, by extension, the strength of the ISM.

Modern Oceanographic Setting of the Indian Ocean
The northern Indian Ocean surface circulation seasonally reverses due to the monsoon [42]. During the ISM, the East Indian Coastal Current flows northward along the eastern Indian coast in the Bay of Bengal, whereas the West Indian Coastal Current flows equatorward (Figure 1), following the west Indian coast in the Arabian Sea [2,43]. However, both the East and West Indian Coastal currents reverse direction during the winter monsoon. The other notable surface circulation that changes direction in the southwest and northeast monsoon currents links the western Java coastal waters to the Arabian Sea [2]. The westward flowing South Equatorial Current between 5 • and 15 • S brings western Pacific waters through the Indonesian throughflow and around 5 • N of the Somali margin; these waters turn eastward and thus feed the thermocline [44]. The westward flowing North Equatorial Current in winter also brings surface waters through the Strait of Malacca, reaching 7.5 • N in the Bay of Bengal.
The deep Indian Ocean water masses originate mainly from the south. In the southeast, deep waters comprise equal amounts of lower Circumpolar Deep Water and Antarctic bottom water [45,46] and enter into the southeast Indian Ocean through the Perth Basin, providing about one-quarter of the global Antarctic bottom water ( Figure 1). This Antarctic bottom water originates from the Ross Sea and travels along the Wilkes-Adélie Coast [47][48][49], where it spills into the central Indian Basin through numerous deep gaps of the 90 • E ridge [46,50]. In the southwest, lower Circumpolar Deep Water (>3.8 km) enters into the Crozet-Madagascar and Mozambique Basins and is overlain by the upper Circumpolar Deep Water (between 2 and 3.8 km; [48]). These waters flow into the Somali and Arabian basins as the western boundary current after passing through the Amirante Passage [46,51]. In the northern Indian Ocean, the nutrient content of the northwardflowing deep waters progressively ages due to lack of ventilation and is not modified by mixing with any other deep-water mass [48,[51][52][53]. The only deep-water mass exiting from the Indian Ocean in the 2-3.5 km depth range is the Indian Deep Water [48,51,54], which is a mixture of upwelled bottom water and the upper Circumpolar Deep Water; it flows into the western Indian Ocean as the southward boundary current [55].
The Antarctic bottom water forms mainly at the Weddell Sea and Wilkes-Adélie coast of the Antarctic margin through a combination of surface water cooling, wind stress, and salt gain in the course of seasonal sea ice formation [48,54]. This cold and dense water floods the deepest part of the Antarctic Ocean and moves northward to fill the global ocean. The sinking of the Antarctic bottom water and Antarctic Intermediate Water balances the upwelling of Circumpolar Deep Water, which is fed by injection of the North Atlantic Deep Water between the Antarctic Circumpolar Current and Antarctic Polar Front.

Materials and Methods
Marine sediments from the Arabian Sea and Bay of Bengal basin are archives representing a complete record of the erosional history of the Himalayas over time, allowing reconstruction of past ISMs. Records are available in these basins; however, most of these records [57][58][59][60] are unavailable, have an insufficient temporal resolution, or have a length of records inadequate for this study. Consequently, we have focused our study mainly on cores ODP 758 (this study; [23], RC12-343 [21], and CR2 [22,24] and use other minor data, such as Hein et al. [59], Liu et al. [61], Naik et al. [25], and Yu et al. [62], for the discussion.

Samples and Data
The ODP site 758 was drilled on the 90 • E ridge in 1989 at 2925 m water depth ([63]; Figure 1). Published ε Nd data of cores RC12-343 (15.2 • N; 90.6 • E) and SK129-CR2 (3 • N; 76 • E) collected at 2666 and 3800 m water depths, respectively, detailed in Stoll et al. [21], Piotrowski et al. [22], and Wilson et al. [24], were also used in the study and are introduced in the respective sections. Data from piston cores RC12-344 [26] and MD77-176 [62], collected on the northern Andaman Sea shelf and eastern Bay of Bengal (exit pathway of the Andaman Sea water) at 2100 and 1375 m water depths, respectively, were also taken into consideration.

Revision of the Stratigraphy
The SPECtral MAPping Project orbital chronology [64] was used to construct the stratigraphy of ODP site 758 and core RC12-343 [21]. Piotrowski et al. [22] used the 14 C-AMS dates and the earlier version of the benthic oxygen isotope stack of Lisiecki and Raymo [65] to build the age model of core SK129-CR2. We revised the age models due to the availability of the latest version of the radiocarbon calibration program CALIB 8.1 [66] and Marine20 reservoir ages [67] for the 14 C-AMS dates used in cores SK129-CR2 and the deep Indian Ocean benthic stack of Lisiecki and Stern [68]. A brief description of the age model in each core is given below.
ODP Site 758: Farell and Janecek [63] reported stratigraphy of the ODP Site 758, measuring oxygen and carbon isotopes in the planktonic foraminifer Globigerinoides sacculifer and benthic foraminifer Cibicidoides wuellerstorfi, magnetic stratigraphy, and astronomical tuning. We updated the age model of site 758 using the stratigraphic framework proposed by Lisiecki and Stern [68], applying the deep-Indian regional stack.
SK129-CR2 (hereafter CR2): Thirteen 14 C-AMS dates were converted to calendar years using the Fairbanks calibration curve 01.07 [69] reported by Piotrowski et al. [22] and Wilson et al. [24]. Moreover, the authors applied a uniform 350-year reservoir age following the reports by Butzin et al. [70] and Cao et al. [71] to these 14 C-AMS dates. The 14 C-AMS dates were converted to calendar years using the CALIB 8.1 [66] and Marine20 reservoir ages [67] to maintain consistency and harmonize the age model with other records used in this study.
RC12-343: Planktonic foraminifera Globigerinoides ruber (white) oxygen isotopes were used to construct the stratigraphy, which was tuned to the SPECtral MAPping Project orbital stratigraphy [64] to build the age model of core RC12-343 [21]. The outdated SPECtral MAPping Project chronology is no longer operational, and hence we updated the stratigraphy of RC12-343, applying the deep-Indian regional stack of Lisiecki and Stern [68].
RC12-344: Nine 14 C-AMS dates were converted to calendar years using the older versions of the CALIB program [26,72] to construct the age model. The 14 C-AMS dates were converted to calendar years using the updated CALIB 8.1 [66] and Marine20 reservoir ages of Heaton et al. [67].

ODP Site 758
The analytical protocol for ε Nd seawater extraction and measurement was detailed in Gourlan et al. [23] and Rashid et al. [26]. For ε Nd and ε Hf of detrital fraction analyzed on the same sediments, we briefly describe both protocols below. Approximately 300-400 mg of pulverized sediments were treated by CH 3 COOH (1N) in excess, which allowed a total dissolution of the carbonate fraction and leaching of most of the Mn-oxides around the microfossils in which the seawater ε Nd was trapped as well as the intratest organic matter [73]. As the seawater signal was already analyzed, the supernatant fluid carrying the seawater signal was removed. The solid residue was rinsed twice with water, ultrasonically agitated in 4 mL of 1N HBr, and centrifuged to remove the remaining oxides to extract the detrital signal. Then, the residue was digested for one week using a classic mixture of HF-HNO 3 -HClO 4 (3mL-3mL-15 drops) at 135 • C in Parr Bombs. Nd and Hf chemical separation were based on the method published by Chauvel et al. [74] and performed at the Université de Savoie Mont-Blanc in Grenoble. Total procedural blanks for Nd and Hf were 33 pg (n = 11) and 12 pg (n = 6), respectively, accounting for less than 0.5% of the sample in total. The ε Nd and ε Hf were measured on a Nu Plasma MC-ICP-MS in ENS-Lyon and corrected for mass fractionation bias using 146 Nd/ 144 Nd = 0.7219 and 179 Hf/ 177 Hf = 0.7325, respectively. The Ames-Rennes-Nd and Ames-Grenoble Hf reference standards were run every two or three samples and yielded an average value of 143 Nd/ 144 Nd = 0.511964 ± 0.000019 (2σ, n = 34) and 176 Hf/ 177 Hf = 0.282161 ± 0.000025, respectively, for all the sessions; they were corrected after every session to the recommended values for these two standards [74,75]. Several complete duplicate analyses were also performed. For Nd and Hf isotopic ratios, the measurements were reproduced within analytical errors (0.25 ≤ ε Nd or ε Hf ).

Sample Preparation Methods of the Published Records
It would have been ideal to follow a uniform sample preparation protocol to determine the seawater ε Nd data of the three major sites (RC12-343, ODP Site 758, and SK129-CR2) used in this study. However, Stoll et al. [21], Piotrowski et al. [22], and Gourlan et al. [23] used two different sample preparation methods, which are (i) leaching of authigenic sedimentary Fe-Mn oxyhydroxides from bulk sediments and (ii) extraction of Nd by dissolving planktonic foraminifers. It is beyond the scope of this paper to provide a detailed description of these methods; however, in short, Piotrowski et al. [22] followed a technique analogous to Gourlan et al. [73] to prepare samples for ε Nd determination and Stoll et al. [21] acquired ε Nd data from planktonic foraminifers, building on the work of Burton and Vance [57] and following the Cd/Ca cleaning protocol of Boyle and Keigwin [76]. The method comprises repeated ultrasonication in ultra-pure water and methanol to dislodge the clay and other fine-grained sediments adhering to the foraminiferal tests subjected to oxidative and reductive cleaning to remove the Fe-Mn oxide coatings. The cleaned tests were dissolved using nitric acid to determine the ε Nd [21].

Paleoclimate Records of the Last Glacial Cycle in the Northern Indian Ocean
The composite data of the three records selected (presented in Figure 2) provide a latitudinal transect that allows us to evaluate changes in the surface and deep-water circulation in the northern Indian Ocean. As previously described, we revised the age models by tuning records to the regional deep-Indian benthic δ 18 O stack [68]. Hence, any minor leads and lags between the water mass properties were not considered; instead, we focused our discussion on the broader aspects of these data and their implications for understanding the surface and deep ocean circulation.
At site 758, cold glacial ε Nd were more radiogenic (~−7.4 to −8.4), whereas warm interglacial periods were less radiogenic (~−11 to −9.5). This variability in ε Nd is in tandem with the benthic oxygen isotopic (δ 18 O) stack of the glacial-interglacial cycle [68], suggesting an overriding link between the Northern Hemisphere ice-sheet growth and decay and the tropical hydrological cycle. To explain these fluctuations, Gourlan et al. [23] considered the primary sources of Nd contributing to ε Nd values in the northern Indian Ocean: (i) the Himalayan rivers (less radiogenic input; [59,62,77]), (ii) the Indonesian Throughflow (where radiogenic input is related to Pacific waters passing through the complex Indonesian volcanic islands system; [23,53,[78][79][80][81]), and (iii) the North Atlantic Deep Water (less radiogenic input, which is transported into the southern Indian Ocean by the Antarctic Circumpolar Current [79,82,83]. There are other local Nd sources, like the local wind [84,85], minor inputs from Crozet and Kerguelen plateaus [86], and southern African terrains [87,88]. However, the impact of these sources on the Nd budget is considered limited compared to the three primary sources (i-iii) and thus is inadequate to explain the large glacial-interglacial ε Nd fluctuations of the northern Indian Ocean. Considering the small ε Nd variations recorded at ODP Site 757, located in front of the Indonesian throughflow, Gourlan et al. [23] considered this primary source negligible to modify the ε Nd values in the Indian Ocean. Moreover, previous studies of salinity distributions, Nd isotopic composition of marine sediments, and δ 18 O foraminifera from the Bay of Bengal [58,62,80,[89][90][91] have suggested a mixing between the mean Indian Ocean and fresh water and sediments derived by the Ganga-Brahmaputra-Meghna outflow. Hence, we hypothesize that the ε Nd values at site 758 were modulated mostly by changes in the Ganga-Brahmaputra-Meghna outflow [23] rather than a deep-water circulation. In this context, higher radiogenic values resulted from the dramatic reduction of the Ganga-Brahmaputra-Meghna outflow during the glacial period, while during the interglacial periods, these outflows were significantly increased, resulting in less radiogenic ε Nd.
In contrast to site 758, Stoll et al. [21] reported seawater ε Nd data using surface-living foraminifers, following an identical method used at site 758 by Burton and Vance [57] for the last 195 ka from core RC12-343. The glacial-interglacial ε Nd varies from radiogenic (−6.95) to less radiogenic (−11.12), consistent with the findings at site 758 [23]. Moreover, using a suite of cores, Stoll et al. [21] reconstructed the meridional gradients between 20 • and 5 • N for the Holocene and LGM periods ( Figure 3). The authors demonstrated that meridional gradients in Holocene ε Nd , where the southern sites are more radiogenic compared to those of the northern sites, are consistent with the findings of modern seawater ε Nd [53,80,92]. Padmakumari et al. [93] also showed a north-south ε Nd gradient in the glacial surface waters of this region. In this regard, less radiogenic ε Nd was attributed to dilution by the Ganga-Brahmaputra-Meghna outflow at the northern sites, while higher radiogenic values in the southern sites were assumed to be associated with less sea-surface dilution.   These modern fingerprinting data were then used to interpret the LGM εNd, where the north-south gradients appear to be significantly weakened, suggesting a reduction in dilution by the Ganga-Brahmaputra-Meghna outflow at the northern sites [21]. By extension, this suggests that the ISM was weak, and a dry climate prevailed during the LGM. The findings of Stoll et al. [21] are consistent with other independent paleo-proxy records of the ISM, such as ancient seawater oxygen isotopic composition (δ 18 Osw), which also suggests a dry climate in the Bay of Bengal and the Andaman Sea during the LGM [72,95,96].
Compared to the data of Bay of Bengal sites, Piotrowski et al. [22] reported seawater εNd for the last 170 ka from the northern tropical Indian Ocean core CR2 (Figure 1). It was collected at 3800 m water depth located at 12° and 2° S, i.e., 2094 and 1596 km from the Bay of Bengal sites RC12-343 and 758, at 2666 km and 2925 km water depths, respectively. The authors demonstrated that the glacial εNd was more radiogenic and averaged ~−7.5, whereas interglacial εNd was less radiogenic and varied from −8.4 to −10 ( Figure 2). In a subsequent study using the same core, Wilson et al. [24] extended the εNd record up to 250 ka, following a method similar to Gourlan et al. [23], and tested various sample preparation methods (see above). The authors did not find any changes in the εNd values of the interval between 0 and 145 ka, the interval of interest of Piotrowski et al. [22]. In any event, there are similarities between the pattern (Figure 2) of ancient seawater εNd and δ 13 C These modern fingerprinting data were then used to interpret the LGM ε Nd , where the north-south gradients appear to be significantly weakened, suggesting a reduction in dilution by the Ganga-Brahmaputra-Meghna outflow at the northern sites [21]. By extension, this suggests that the ISM was weak, and a dry climate prevailed during the LGM. The findings of Stoll et al. [21] are consistent with other independent paleo-proxy records of the ISM, such as ancient seawater oxygen isotopic composition (δ 18 O sw ), which also suggests a dry climate in the Bay of Bengal and the Andaman Sea during the LGM [72,95,96].
Compared to the data of Bay of Bengal sites, Piotrowski et al. [22] reported seawater ε Nd for the last 170 ka from the northern tropical Indian Ocean core CR2 (Figure 1). It was collected at 3800 m water depth located at 12 • and 2 • S, i.e., 2094 and 1596 km from the Bay of Bengal sites RC12-343 and 758, at 2666 km and 2925 km water depths, respectively. The authors demonstrated that the glacial ε Nd was more radiogenic and averaged~−7.5, whereas interglacial ε Nd was less radiogenic and varied from −8.4 to −10 ( Figure 2). In a subsequent study using the same core, Wilson et al. [24] extended the ε Nd record up to 250 ka, following a method similar to Gourlan et al. [23], and tested various sample preparation methods (see above). The authors did not find any changes in the ε Nd values of the interval between 0 and 145 ka, the interval of interest of Piotrowski et al. [22]. In any event, there are similarities between the pattern (Figure 2) of ancient seawater ε Nd and δ 13 C (measured in the epifaunal benthic foraminifera C. wuellerstorfi, believed to faithfully record the bottom-water dissolved inorganic carbon, if not influenced by the high productivity in the surface water). However, a divergence between these water mass properties from 90 to 120 ka is prominent. Piotrowski et al. [22] interpreted that the fluctuation of ε Nd was due to changes in the proportion between the North Atlantic Deep Water and Circumpolar Deep Water, i.e., changes due to the bottom water circulation. Using this reasoning, the North Atlantic Deep Water carries less radiogenic ε Nd . Simultaneously, the Circumpolar Deep Water contains more radiogenic ε Nd, and any changes in the ratio between these water masses control the bottom water ε Nd values at site CR2 [22]. Wilson et al. [24] put forward a new hypothesis to explain the ε Nd data of Piotrowski et al. [22], suggesting that the ε Nd data reflect changes in the global meridional overturning circulation rather than the simple on and off mechanism of the North Atlantic Deep Water formation.
In summary, the Bay of Bengal and Andaman Sea records [21,23,26] emphasize that the Indian and East Asian rivers draining the greater Himalayas dominantly influenced the ε Nd values in the northern Indian Ocean, whereas the northern tropical Indian Ocean records [22,24] mostly favor deep ocean circulation, mainly by the Antarctic bottom water and North Atlantic Deep Water. Therefore, the question is whether the regional sediment discharge or a distant source modulates the ε Nd values in the northern Indian Ocean, which we explore below.

Factors Controlling the ε Nd Signals in the Northern Indian Ocean
The question of what controls the northern Indian Ocean ε Nd values is documented in many sedimentary records. The two main factors that can provide information about the ε Nd values are briefly discussed below.
The first factor could be the differences in the sample preparation techniques. Wilson et al. [24] tested three methods of Nd extraction [57], Stoll et al. [21], Gourlan et al. [73], and Piotrowski et al. [22]) using samples from the same tropical northern Indian Ocean core-CR2used by Piotrowski et al. [22]. To summarize, (a) bulk sediments were used to extract Nd from ferromanganese coatings on detrital and biogenic grains, i.e., sediment leaching; (b) Nd was extracted from the planktonic foraminifers following the reductive-oxidative cleaning protocol of Boyle and Keigwin [76]; (c) the cleaning protocols of Boyle and Keigwin [76] were applied to fish teeth and fish debris, omitting the oxidative step to extract Nd [24,97]. Surprisingly, Wilson et al. [24] reported identical ε Nd values using all three methods. Therefore, it appears that the heterogeneity in ε Nd of the northern Indian Ocean is not a result of the sample preparation techniques. Moreover, the surprisingly similar ε Nd values during the late Holocene and MIS4 at three and two sites (Figure 2), respectively, raise questions about the utility of sophisticated sample preparation techniques and the sea-surface versus deep ocean circulation hypotheses.
The second factor concerning the Northern Indian Ocean surface circulation could be an important issue; a debate exists whether the ε Nd values integrate the upper few hundred meters of the surface ocean or simply the bottom water signature. Piotrowski et al. [22] hypothesized that the North Atlantic Deep Water, or a modified version of the North Atlantic Deep Water, dominantly influenced the northern Indian Ocean deep water. However, the conductivity, temperature, and depth data in the Bay of Bengal [98], modern seawater [80], the GEOTRACERS program [33,92], and sediment core-top ε Nd data appear to contradict such a hypothesis, which was also questioned by Gourlan et al. [23]. It is plausible that the Ganga-Brahmaputra-Meghna outflow predominantly influenced the surface waters at site 758 during the summer, with an ε Nd of −12 to −12.5 [58,62]. However, using the modern monthly averaged sea-surface salinity, Rashid et al. [72,95] demonstrated that site 758 could also be affected by Irrawaddy and Salween river outflow. These waters might contain ε Nd (~11) signatures of eastern and northern Indo-Burman Ranges, as well as the eastern Himalayas, drained by the Irrawaddy and Salween river [99,100]. Thus, the most likely explanation for the ε Nd data at site 758 would be the integrated signature among the Ganga-Brahmaputra-Meghna and Irrawaddy-Salween river outflows and the bottom waters. The annual discharge of the Ganga-Brahmaputra-Meghna, Irrawaddy, and Salween rivers is 1090 km 3 [21], 450 km 3 , and~200 km 3 [26,101], respectively. Thus, these discharges provide clues to the relative contribution of ε Nd from these sources, which requires further documentation. Dia et al. [86] reported ε Nd data from 18 cores from various physiographic provenances (abyssal plains, ridges, aseismic ridges, a fracture zone, and a submarine fan) of the Indian Ocean. The authors then suggested that the northern and central Indian Ocean shows two isotopically homogeneous regions. However, an evaluation of these core tops data reveals no 14 C-AMS dates, which raises questions of whether those core tops represent modern bottom water ε Nd , given the extremely low sedimentation rates in the vast tropical Indian Ocean [102]. Further, Stoll et al. [21] suggested a change in the source input during the glacial period, with a dominance of more radiogenic Nd sources from Arakan coastal rivers (estimated at −6, but this value is uncertain). The authors evaluated a contribution of 80-90% for the Burman Arc/Arakan coastal river sources, with 20-10% for Ganga-Brahmaputra-Meghna sources related to higher precipitation in the Arakan region. Nevertheless, the drainage area and discharge of Arakan coastal rivers are seven and five times smaller, respectively, than those of the Ganga-Brahmaputra-Meghna rivers. This implies that the Arakan coastal rivers' impact on the Nd isotopic value of the Bay of Bengal seawater was most likely limited, as was also recently inferred by Naik et al. [25]. Therefore, the ε Nd seawater value at site 758 is probably a mixture of the Bay of Bengal seawater and the Himalayan river system, which only includes the Ganga, Brahmaputra, Meghna, Salween, and Irrawaddy rivers.
The topical Indian Ocean site CR2 is bathed by a combination of lower Circumpolar Deep Water and Indian Deep Water, with a distinct density, oxygen, silica content, and ε Nd signature [46,82]. A complex array of surface currents, with different surface and thermocline waters depending on the season, exits at site CR2 (Figure 1). During the ISM, the eastward flowing Southwest Monsoon Current delivers the West Indian thermocline water (ε Nd = −7.5) [82]. These waters most likely have a portion from the westward flowing south equatorial current, between 5 • and 15 • S, which brings western Pacific waters through the Indonesian throughflow with an ε Nd signature of −3 to −4 [80,103]. Around 5 • N of the Somali margin, these waters turn eastward and feed the thermocline [44]. In addition, the north-flowing East Indian Coastal Current off India reverses during the ISM, carrying low salinity freshwater with the dissolved Nd (ε Nd values of these waters are presently unknown) below the tip of Sri Lanka, near-site CR2 [2,26,72,95]. In contrast, the Northeast monsoon current in winter also carries dissolved particulates that again exit below the tip of Sri Lanka [2,44].
The westward flowing north equatorial current in winter also brings less radiogenic surface water with ε Nd values of −10.2 to −10.9 [80] (Figure 1) through the Strait of Malacca, which reaches the Bay of Bengal at 7.5 • N, close to site CR2. In summary, the ε Nd data at site CR2 most likely represent integrated water column values. Therefore, it would be insightful to differentiate between the various water masses using the ε Nd when most of these waters have a near-identical signature. Any effort to fingerprint these waters using the ε Nd values would be like trying to find a needle in a haystack. Thus, it appears that the initial intention of using ε Nd to trace deep ocean circulation accurately runs against the grain of not being complicated by various factors, at least in the Indian Ocean, as is the case for the paired use of Cd/Ca and δ 13 C in benthic foraminifers in the Southern Ocean [104].

The Questions Regarding the Northern Indian Ocean ε Nd Data
The three ε Nd curves (Figure 2) co-vary with one another despite their variations in data resolution and the absolute ε Nd values of the northern Indian Ocean. Our ε Nd data at site 758 are consistently less radiogenic throughout the last glacial cycle than ε Nd values at site CR2 (Figure 2d). It should be acknowledged that the revised age model for the ODP site 758 indicates a gap in ε Nd data, not only in the MIS4 but also much wider between 50 and 100 ka (Figure 2). The penultimate glacial and LGM periods show the highest radiogenic values, regardless of whether the records were generated close to the Ganga-Brahmaputra-Meghna outflow or in deeper water such as site CR2. Furthermore, the penultimate glacial, MIS4 and LGM ε Nd data are identically more radiogenic at sites RC12-343 and CR2. One of the most striking features is that the lowest radiogenic ε Nd values are at the northern site of RC12-343 for the penultimate interglacial (MIS5e), compared to higher radiogenic values at sites 758 and CR2. However, there are particularities in the ε Nd values at site RC12-343, showing crossovers such that (i) radiogenic ε Nd values are higher compared to the records of CR2, from 78 to 116 ka, which we label "Zone-I"; (ii) identical ε Nd values are found in cores RC12-343 and CR2 during the MIS4; (iii) intermediate ε Nd values during the MIS6 and from MIS3 to Holocene are found between cores CR2 and site 758 (Figure 2).
The lowest ε Nd values during the MIS5e could be explained by the wettest period at RC12-343, and by extension, the strongest ISM, resulting in the largest Ganga-Brahmaputra-Meghna outflow and East Asian rivers' discharge [21,99]. This hypothesis would be consistent with the large body of geochemistry data from the Bay of Bengal and Andaman Sea [105,106] and the δ 18 O records from the greater Asian monsoon system [107][108][109]. By extension, the lack of such light ε Nd values during the MIS5e at CR2 would suggest that it would be decoupled from the dynamics of the ISM if the site were to be influenced by the monsoons. Surprisingly, the late Holocene ε Nd values are statistically identical for the three sites (Figure 2), consistent with the new northern Bay of Bengal records from Achyuthan et al. [110] and Naik et al. [25]. The data of Naik et al. [25] indicate a dominant monsoonal impact in the surface and deep central BoB since 8 ka, increasing riverine particulates from the Ganga-Brahmaputra-Meghna system with the release of unradiogenic Nd.
The other striking feature is the more radiogenic ε Nd in the northern site at RC12-343 than the southern site at CR2 in Zone-I. It would have been easier to explain that the Ganga-Brahmaputra-Meghna outflow was weaker, and hence the relative flux of ε Nd compared to the southern site was smaller for this interval. The other possibility is that the position of the cores in terms of deep-water settings are such that the site RC12-343 is bathed by the upper Circumpolar Deep Water, while the southern site is bathed by lower Circumpolar Deep Water/Antarctic bottom water. Bottom water δ 13 C data at site RC12-343 are required but currently unavailable to understand the discrepancy of the ε Nd values further. The bottom water δ 13 C data between 0 and 78 ka at site 758 are almost identical to those of the CR2, even though the latter is located~1000 m deeper than the former, which may suggest that both sites were bathed by the same water mass [111], and similar surface water productivity prevailed.

A New Paradigm for the Northern Indian Ocean Neodymium Isotopes during the Last Deglaciation
To better understand the global marine Nd cycle, Du et al. [31] provide a synthesis of authigenic ε Nd records using the "bottom-up" concept for Nd budget in the deep ocean [112,113]. The authors used the modern global seawater ( [37,85,114,115]; and references therein) and paleo ε Nd data in a box model to evaluate the sensitivity of ε Nd as an ocean circulation tracer. In so doing, Du et al. [31] partitioned the global historical ε Nd data into four categories: (i) a deglacial time series, (ii) a Holocene time slice (0-6 ka); (iii) an LGM time slice (19-23 ka); (iv) a Heinrich Stadial 1 (HS1) time slice (15-18 ka). We used the deglacial time series to evaluate further our ε Nd , including the published northern Indian Ocean ε Nd data. Incorporating more than 270 seawater ε Nd records, Du et al. [31] performed principal component (PC) analysis in which~86% and~6% of the total variance of ε Nd values correspond to PC1 and PC2, respectively. In this context, PC1 reflects smooth ε Nd changes from the LGM to Holocene, with a rapid shift in mid-late Holocene.
A considerable degree of variability appears between the tropical northern Indian Ocean and Bay of Bengal ε Nd values ( Figure 4). However, a similar overall ε Nd trend, i.e., glacial values are more radiogenic than those of the Holocene, despite these data were generated using different methods (see Section 3.3) in different carriers of the ε Nd . Padmakumari et al. [93] also reported such a trend (i.e., a significant ε Nd shift) in the planktonic foraminifera from the LGM to Holocene from the Bay of Bengal. This lighter ε Nd trend from glacial to Holocene most likely suggests a transition from a weaker to stronger ISM, resulting in more river outflow and sediment discharge, consistent with other paleo-proxy records of the ISM, namely the seawater oxygen isotopes (δ 18 O sw ) [72,95,116,117]. Moreover, the northern Bay of Bengal and Andaman Sea ε Nd data differ from those of the tropical Indian Ocean (CR2). The former reflects millennial-scale variability, including the large changes in ε Nd . Using ε Nd data from the Andaman Sea core RC12-344, Rashid et al. [26] hypothesized that the contribution from the Irrawaddy-Salween river was due primarily to the erosion of the Indian subcontinent reflecting an alternation between strengthening and weakening of the ISM [94,108].
Our new ε Hf and ε Nd data at site 758 show that (i) the detrital fraction has a composition falling in the field of clays and biogenic muds in the Nd-Hf space without significant influence from sands, (ii) the Nd isotopic values obtained on the detrital fraction at site 758 fall in the range of the Himalayan values and are similar to "Higher Himalayan Crystalline" with a contribution of "Tethyan Sedimentary Sequence" [118]; (iii) Hf isotopic values of the detrital fraction do not correspond to the Himalayan values. This latest finding can be explained by the more complex behavior of Hf in sediments. In contrast to Nd, Hf isotopic ratios of sediments vary as a function of both source changes and mineralogy of the sediments. Here, the sediments include only fine-grain fraction without zircon; this fraction is not representative of the source rock for the Lu-Hf system. Part of the Hf source signal remained close to the continent within zircon-rich sands and is not sampled in this study. Moreover, a small shift is observed between the LGM and Holocene period for both systems, from higher to lower radiogenic values, indicating an increase of the riverine input or/and a change in the detrital source. Using authigenic FeMn oxy-hydroxides, Gutjahr et al. [35] provide similar reasoning to explain higher to lower radiogenic ε Hf values from the northern North Atlantic. This shift independently confirms the trend of ε Nd and δ 18 O sw , thus suggesting a reduction in the sediment discharge during the last deglacial from the Holocene (Figure 4). Recently, Naik et al. [25] indicated that the northern Bay of Bengal ε Nd values during the deglacial period were mostly controlled by the Antarctic bottom water, whereas the late Holocene values were influenced by the Ganga-Brahmaputra-Meghna discharge, inconsistent with the findings at sites RC12-344 and 758.  [22,24]; (b) Nd seawater data of cores SK129-CR2 [22], ODP 758 [23], and RC12-343 [21], in which various sample preparation methods were used; (c) Nd data of cores MD77-176 [62], SK157-20 [25], and RC12-344 [26]; (d) hafnium isotopes Hfdetr (black) and Ndetr (blue) from the detrital fraction of core ODP 758 (this study); (e) principal component 1 (PC1) of the global seawater and paleo-seawater Nd data [31]; (f) oxygen isotopes of North Greenland Ice Core Project [119].
The PC1 scores of the global Nd data appear to be a mirror image of both Bay of Bengal and tropical Indian Ocean Nd data (Figure 4). Du et al. [31] estimated a four to five epsilon unit offset between the north Indian Ocean Nd and the global PC1 scores. This large deviation cannot be explained by alternating sources between the North Atlantic Deep Water and Antarctic bottom water, as the modern Nd values of North Atlantic Deep Water and Antarctic bottom water are ~−13.5 and ~−4, respectively. Furthermore, the modern Bay of Bengal Nd values vary from −10.5 to −11 [53], in which the particulate phases of Ganga-Brahmaputra-Meghna rivers mostly impact the dissolved Nd. Hence, it appears that the Nd isotopes from the Bay of Bengal cannot be used as a tracer for the mixing of water masses in the global ocean. This implies that the Nd isotopes can be used as a robust  [22,24]; (b) ε Nd seawater data of cores SK129-CR2 [22], ODP 758 [23], and RC12-343 [21], in which various sample preparation methods were used; (c) ε Nd data of cores MD77-176 [62], SK157-20 [25], and RC12-344 [26]; (d) hafnium isotopes ε Hfdetr (black) and ε Ndetr (blue) from the detrital fraction of core ODP 758 (this study); (e) principal component 1 (PC1) of the global seawater and paleo-seawater ε Nd data [31]; (f) oxygen isotopes of North Greenland Ice Core Project [119].
The PC1 scores of the global ε Nd data appear to be a mirror image of both Bay of Bengal and tropical Indian Ocean ε Nd data ( Figure 4). Du et al. [31] estimated a four to five epsilon unit offset between the north Indian Ocean ε Nd and the global PC1 scores. This large deviation cannot be explained by alternating sources between the North Atlantic Deep Water and Antarctic bottom water, as the modern ε Nd values of North Atlantic Deep Water and Antarctic bottom water are~−13.5 and~−4, respectively. Furthermore, the modern Bay of Bengal ε Nd values vary from −10.5 to −11 [53], in which the particulate phases of Ganga-Brahmaputra-Meghna rivers mostly impact the dissolved Nd. Hence, it appears that the Nd isotopes from the Bay of Bengal cannot be used as a tracer for the mixing of water masses in the global ocean. This implies that the Nd isotopes can be used as a robust proxy for reconstructing terrestrial discharge adjacent to continents, consistent with the Nd isotopes of the modern Amazon discharge [120]. Using a sediment core from the central Bay of Bengal, Naik et al. [25] completely discredited the role of the dissolved or particulate Nd from the Ganga-Brahmaputra-Meghna rivers by suggesting that the effects of surface particles on seawater Nd isotopes are mostly confined to the surface layers of the oceans, with little to no impact on bottom water [121]. Assuming that the assumptions of Naik et al. [25] are correct, and if the contributions from the Ganga-Brahmaputra-Meghna and Irrawaddy-Salween rivers were held constant or no contribution from these rivers to the northern Indian Ocean was considered, the sole supply from the North Atlantic Deep Water, modified North Atlantic Deep Water, or Antarctic bottom water simply fail to explain the observed ε Nd values. Realizing the shortcomings of the ε Nd values, Du et al. [31] hypothesized that the ε Nd values in the northern Indian Ocean are decoupled from the global ocean. To explain the ε Nd values at site RC12-344, Rashid et al. [26] argued that the Andaman Sea ε Nd values were mostly modulated by the Irrawaddy and Salween discharge with stiff (and unprofessional) resistance from Nd community. However, Du et al. [31] appear to confirm the hypothesis that the discharge proximal to continental margin must be the dominant factor in determining Nd isotopes.

Summary
Nd isotopes have been used in the northern Indian Ocean, Bay of Bengal, and the Andaman Sea during the last decade to reconstruct the past strength of the ISM and deep ocean circulation. One of the factors that made the Nd isotopes a good proxy is its apparent conservative nature, undiluted by the regional influences such as the complications associated with the Mg/Ca or B/Ca ratios in biogenic calcite. The signature of the Nd isotopes in the global ocean was believed to be governed by binary sources, i.e., less radiogenic old continental crust and more radiogenic relatively young volcanic rocks, which were distributed by the North Atlantic Deep Water and Antarctic bottom water. However, in this broader understanding, the boundary exchange or regional sources of Nd were considered insignificant to the global Nd budget. This long-held paradigm appears to be weak [31,113], which allows our northern Indian Ocean Nd isotopes to offer insight.
By compiling a large set of ε Nd data from northern Indian Ocean marine sediments, including our ε Nd data, this study highlights the significant contribution to the northern Indian Ocean Nd budget of the Himalayan sources drained by the Ganga-Brahmaputra-Meghna and Irrawaddy-Salween rivers and suggests that the ISM was the strongest and wettest during the last interglacial and Holocene periods, as inferred from less radiogenic ε Nd values. In contrast to the interglacial periods, glacial ε Nd and ε Hf values are more radiogenic, which suggests a reduction in the Ganga-Brahmaputra-Meghna and Irrawaddy-Salween discharge and implies less erosion of the Indian subcontinent, which could be tied to the weaker ISM. Our hypothesis is consistent with the global analysis of the modern seawater Nd budget and modeling results of the past seawater ε Nd data (e.g., [31]). However, high-resolution seawater and detrital ε Nd and ε Hf data from the locus of the Ganga-Brahmaputra-Meghna outflow are needed to test the veracity of the hypothesis.