Carbon isotopic ratios of modern C₃ and C₄ vegetation on the Indian peninsula and changes along the plant–soil–river continuum – implications for vegetation reconstructions

The large difference in the fractionation of stable carbon isotopes between C₃ and C₄ plants is widely used in vegetation reconstructions, where the predominance of C₃ plants suggests wetter and that of C₄ plants drier conditions. The stable carbon isotopic composition of organic carbon (OC) preserved in soils or sediments may be a valuable (paleo-)environmental indicator, based on the assumption that plant-derived material retains the stable carbon isotopic value of its photosynthetic pathway during transfer from plant to sediment. In this study, we investigated the bulk carbon isotopic values of C₃ and C₄ plants (δ¹³C) and of organic carbon (δ¹³C_org) in soils, river suspended particulate matter (SPM) and riverbed sediments to gain insight into the control of precipitation on C₃ and C₄ plant δ¹³C values and to assess changes in δ¹³C_org values along the plant–soil–river continuum. This information allows us to elucidate the implications of different δ¹³C end-members on C₃ / C₄ vegetation reconstructions. Our analysis was performed in the Godavari River basin, located in the core monsoon zone in peninsular India, a region that integrates the hydroclimatic and vegetation changes caused by variation in monsoonal strength. The basin has distinct wet and dry seasons and is characterised by natural gradients in soil type (from clay-rich to sandy), precipitation (∼ 500 to 1500 mm yr⁻¹) and vegetation type (from mixed C₃ / C₄ to primarily C₃) from the upper to the lower basin. The δ¹³C values of Godavari C₃ plants were strongly controlled by mean annual precipitation (MAP), showing an isotopic enrichment of ∼ 2.2 ‰ from ∼ 1500 to 500 mm yr⁻¹. Tracing δ¹³C_org values from plant to soils and rivers revealed that soils and riverbed sediments reflected the transition from mixed C₃ and C₄ vegetation in the dry upper basin to more C₃ vegetation in the humid lower basin. Soil degradation and stabilisation processes and hydrodynamic sorting within the river altered the plant-derived δ¹³C signal. Phytoplankton dominated the δ¹³C_org signal carried by SPM in the dry season and year-round in the upper basin. Application of a linear mixing model showed that the %C₄ plants in the different subbasins was ∼ 7 %–15 % higher using plant end-members based on measurement of the Godavari vegetation and tailored to local moisture availability than using those derived from data compilations of global vegetation. Including a correction for the ¹³C enrichment in Godavari C₃ plants due to drought resulted in maximally 6 % lower estimated C₄ plant cover. Our results from the Godavari basin underline the importance of making informed choices about the plant δ¹³C end-members for vegetation reconstructions, considering characteristics of the regional vegetation and environmental factors such as MAP in monsoonal regions.

1 Introduction

Vegetation reconstruction of the coverage of C₃ and C₄ plants uses the distinct δ¹³C values of both vegetation types with the assumption that organic carbon (OC) retains the stable carbon isotopic composition during transfer from plant to soils and sediments (e.g. Koch, 1998; Dawson et al., 2002; Wynn and Bird, 2007). The C₃ and C₄ photosynthetic pathways fractionate carbon isotopes to a different extent; this is reflected in bulk δ¹³C values of ∼ −20 ‰ to −37 ‰ in C₃ and ∼ −10 ‰ to −16 ‰ in C₄ plants (e.g. Bender, 1971; Farquhar et al., 1989; Kohn, 2010). C₃ plants fix CO₂ using the Calvin–Benson cycle, and they are prevalent in relatively cold and humid environments. C₄ plants add an initial CO₂ fixation step using phosphoenolpyruvate (PEP) carboxylase to concentrate CO₂ via bundle sheath cells inside the leaf. This additional step in the photosynthesis pathway allows for C₄ plants to maintain relatively high photosynthesis rates under low stomatal conductance with limited water loss, which enables them to thrive in high temperature and (semi-)arid environments (e.g. Farquhar, 1983; Sage and Monson, 1999; Sage, 2004).

Given the sensitivity of vegetation type to water availability, shifts in the relative contribution of C₃ and C₄ plants can be used to infer hydroclimatic changes, with a shift to more dominant C₄ vegetation in drier periods or areas and more dominant C₃ vegetation in wetter conditions (Koch, 1998; Sage, 2004; Still et al., 2003, and references therein). The bulk carbon isotopic composition of organic matter (δ¹³C_org) preserved in soils or in river-dominated sediments in the marine realm is often used as a proxy for vegetation reconstructions, as it is considered to represent an integrated signal of the vegetation (e.g. Galy et al., 2007; Sarangi et al., 2021). For vegetation reconstructions, different techniques can be employed, including bulk or compound-specific isotope analyses. Although less source specific, bulk isotope analyses provide a low-cost, high-throughput approach that can be applied at high resolution and/or for large geographic areas, also in (sub)tropical regions where carbon- and vegetation-specific compound concentrations are generally low and may undergo compound-specific degradation patterns and/or differential settling into sediments (e.g. Hou et al., 2020; Li et al., 2020). Subsequently, isotope mixing models provide a means to infer the distribution of C₃ and C₄ plants and changes therein at spatial or temporal scales, which is particularly relevant in the context of changing climatic conditions affecting C₃ / C₄ vegetation patterns. For instance, gradual aridification in a basin can result in drought-stressed C₃ plants, as well as an increased abundance of aridity-adapted C₄ plants. Both changes result in an increase in bulk δ¹³C_org values. A shift to more dominant C₄ vegetation is important to identify, as a shift to a dry ecosystem indicates a reduced resilience to changes in moisture availability (e.g. Cui et al., 2017; Ghosh et al., 2017).

The monsoon-influenced Indian subcontinent is particularly sensitive to changes in hydroclimate on both short (seasonal) and long (orbital) timescales (Turner and Annamalai, 2012; Sinha et al., 2011, 2015; Banerji et al., 2020; Dutt et al., 2021), resulting in changes in the C₃ and C₄ vegetation distributions over the Neogene and Quaternary period (e.g. Agrawal et al., 2012; Ghosh et al., 2017; Basu et al., 2018, 2019a; Roy et al., 2020). For example, changes in bulk δ¹³C_org and leaf-wax-specific δ¹³C values captured the late Miocene (7.4–7.2 Ma) expansion of C₄ plants recorded in Himalayan-derived Indus Fan sediments (Feakins et al., 2020) as well as in Gangetic Plain alluvial sediments and paleosols (Ghosh et al., 2017; Roy et al., 2020). C₄ plants spread southward over the Indian peninsula during the mid-Pliocene to mid-Pleistocene (3.5–1.5 Ma), linked to reduced rainfall (Dunlea et al., 2020), whereas an increase in monsoon strength over the last deglaciation led to a shift to a more C₃-dominated ecosystem (e.g. Galy et al., 2008a; Contreras-Rosales et al., 2014). Finally, an increase in C₄ plants across India from the mid- to late Holocene, establishing the modern-day vegetation, was linked to aridification (Ponton et al., 2012; Contreras-Rosales et al., 2014; Sarkar et al., 2015; Usman et al., 2018; Basu et al., 2019a). At the seasonal to decadal scale, changes in monsoon intensity and distribution can affect plant δ¹³C as well as preservation and provenance of the bulk δ¹³C_org signal in soils and sediments (e.g. Ittekkot et al., 1985; Galy et al., 2008b). However, the effect of monsoon variation on the bulk δ¹³C_(org) signal along the plant–soil–river continuum has not been thoroughly tested for peninsular India.

One approach to reconstructing vegetation is the application of a straightforward, linear isotope mixing model, which requires an informed choice about the C₃ and C₄ plant δ¹³C end-members and their variability, which are used as input. These end-members can be based on averages determined in compilations of vegetation occurring around the globe or be based on (modern) vegetation samples in a region. In addition, their date of sampling with respect to OC turnover rates in soils or sediments and controls by environmental factors need to be considered. Meta-analyses of global C₃ vegetation revealed that hydroclimatic conditions such as rainfall amount and seasonality affect the plant δ¹³C (Diefendorf et al., 2010; Kohn, 2010; Basu et al., 2019b, 2021). In particular drought stress results in a less negative δ¹³C in C₃ plants, where mean annual precipitation (MAP) has a much stronger control than other environmental factors such as temperature or altitude (Stewart et al., 1995; Diefendorf et al., 2010). Observed changes in C₃ plant δ¹³C due to drought stress are the result of changes in the ratio of leaf interior to atmospheric CO₂ concentrations, which can be the result of changes in stomatal conductance, photosynthetic capacity and photosynthetic rate or a combination thereof (Farquhar et al., 1989; Diefendorf et al., 2010; Liu et al., 2013). For C₄ vegetation, carbon isotope fractionation is generally unaffected by drought stress, although a few field and experimental studies have shown that C₄ plant δ¹³C becomes more negative under water-limiting conditions as a result of a less efficient CO₂-concentrating mechanism referred to as “bundle sheath leakiness” (e.g. Buchmann et al., 1996; Yoneyama et al., 2010; Basu et al., 2015; Ellsworth and Cousins, 2016). Plants using an alternative photosynthetic pathway, i.e. crassulacean acid metabolism (CAM), as adaption to aridity photosynthesise during the day and respire at night (i.e. temporal CO₂-concentrating mechanism), where moisture availability determines the expression of the C₃ or C₄ fixation pattern (Sankhla et al., 1975). Under water-stressed conditions CAM plants have stable carbon isotopic values similar to C₄ plants, but they are relatively rare in India (Sankhla et al., 1975; Ziegler et al., 1981) and are, therefore, not further considered here. Regardless, drought stress has the largest impact on C₃ plants and is recognised to cause high intraspecies variability in water-limited ecosystems (<1000 mm yr⁻¹; e.g. Ma et al., 2012; Liu et al., 2013, 2014; Luo et al., 2021). Hence, the existing plant community in a region and the impact of water availability on those plant species may vary locally and result in a C₃ plant δ¹³C value that differs depending on the regional conditions (Liu et al., 2014; Basu et al., 2019b).

Although it is possible to determine region-specific plant δ¹³C end-members that are representative of the regional conditions, including vegetation species, structure, density, agriculture/land use and important environmental controls such as MAP, this approach requires detailed knowledge of the δ¹³C values of the regional C₃ and C₄ vegetation as well as of rainfall distributions. Problematically, such detailed information is often unavailable and yet to be established for Indian plants in the core monsoon zone (CMZ). Furthermore, correction of the plant end-members for drought conditions that (seasonally) prevail in peninsular India requires details on regional rainfall distributions and depends on the δ¹³C value that was derived from measurement of regional vegetation or estimated based on data compilations of global vegetation. Recently, a study of δ¹³C values in region-specific vegetation along a precipitation gradient on the Gangetic Plain prompted a recalculation of the abundance C₃ and C₄ plants in sedimentary deposits accounting for drought-stress-induced ¹³C enrichment in C₃ plants (Basu et al., 2015, 2019b). They showed that earlier investigations likely underestimated the abundance of C₄ plants (∼ 20 %). Recalculation using an end-member and mixing-model approach revealed that C₄ plants existed in this region at an earlier date than anticipated, changing the timing of (Miocene) C₄ grassland expansion on the Gangetic Plain to ∼ 17 Ma (Basu et al., 2015, 2019b). This shift highlights the effect of plant δ¹³C end-member values on vegetation reconstructions.

Next to precipitation controls on vegetation δ¹³C, the initial plant δ¹³C signal may be altered during transit from the plant to the sedimentary archive, depending on physical and biogeochemical processes that determine the stability, i.e. protection against degradation and transport efficiency of this plant-derived OC (e.g. Battin et al., 2009; Ward et al., 2017). Soil and sedimentary deposits integrate a temporal signal, depending on OC turnover rates which are estimated to range from ∼ 10 years in tropical forest soils to ∼ 25–40 years in savanna soils (Martin et al., 1990; Bird et al., 1996). Comparison of the older δ¹³C value of soils and sediment with that of the modern vegetation requires consideration of the Suess effect that describes the rapid decline over the last few decades in the δ¹³C value of atmospheric CO₂ as a result of fossil fuel burning, causing a change in the OC stable carbon isotopic composition of vegetation over time. Furthermore, it is well-established that soil degradation processes lead to a ¹³C enrichment in the remaining soil OC, which is usually estimated to be ∼ 1 ‰–3 ‰ but can be as high as 6 ‰ in tropical and semi-arid regions (e.g. Krull et al., 2005). Possible factors that contribute to this ¹³C enrichment are preferred uptake and degradation to CO₂ of ¹³C-depleted OC by microbes, incorporation of ¹³C-enriched microbial and fungal biomass in the soil, and/or preferential adsorption of ¹³C by fine mineral particles (Krull et al., 2005; Wynn, 2007; Wynn and Bird, 2007). Soil OC thus comprises a mixture of plant-derived, fungal and bacterial biomass, and microbially processed carbon. The different compounds (e.g. lipids, proteins, carbohydrates) differ in their degradability but may also be associated with mineral surfaces, which protects them from degradation. Compound-specific degradation rates or preservation in soils via microbial processing or association with mineral particles may influence the reconstructed C₃ / C₄ vegetation balance depending on the targeted compound, as shown for vegetation and soils in the Gangetic Plain (Sarangi et al., 2021; Roy and Sanyal, 2022). This complex interplay between different inputs and microbial processing may challenge the use of stable carbon isotope ratios for vegetation reconstructions.

Furthermore, the marked hydrological changes in Indian monsoonal rivers can change the source and thereby the δ¹³C value of the OC that it contains at a seasonal scale, from mainly soil-derived OC in the wet season to aquatic produced OC in the dry season, or change its provenance by sourcing from particular parts of the basin with a different vegetation cover in response to the rainfall distribution (Gupta et al., 1997; Balakrishna and Probst, 2005; Aucour et al., 2006; Galy et al., 2008b, 2011; Kirkels et al., 2020a; Menges et al., 2020). For example, Galy et al. (2008b) showed that δ¹³C_org values of suspended particulate matter (SPM) in the Ganges–Brahmaputra River reflected dominant C₃ input in the Himalayan tributaries, but after in-river degradation, this signal was replaced by C₄ inputs in the Gangetic Plain. Finally, hydrodynamic sorting within the river may result in depth-specific OC distributions and thereby influence the δ¹³C signal that is transported downriver (Galy et al., 2008b; Bouchez et al., 2014; Feng et al., 2016; Repasch et al., 2022). Hence, interpretation of δ¹³C-based vegetation reconstructions needs to consider potential alterations during transit from plant source to sedimentary deposits.

In this study, we examine the δ¹³C values of C₃ and C₄ vegetation in the modern-day Godavari River basin, the largest monsoonal river of peninsular India (Fig. 1a), to examine links between monsoon-driven hydroclimate and plant δ¹³C values. In addition, we analyse the bulk δ¹³C_org values in soils, river SPM and riverbed sediments collected in the wet and dry season to explore the evolution of the initial plant-derived δ¹³C signal along the plant–soil–river continuum. Finally, we use our insights into the modern system to assess the influence of drought stress and the use of region-specific plant δ¹³C end-members on C₃ / C₄ vegetation mixing-model estimates and thus the uncertainty in δ¹³C-based vegetation reconstructions.

2 Materials and methods

2.1 Regional setting

The Godavari is the largest peninsular river of India (catchment area: 3.1×10⁵ km², length: 1465 km), with an annual discharge of 110 km³ and sediment load of 170 Mt, of which ∼ 2.8 Mt is OC (Biksham and Subramanian, 1988a, b; Gupta et al., 1997). The Godavari River starts in the Western Ghats mountains and flows across peninsular India before emptying in the Bay of Bengal and is situated in the core monsoon zone (Ponton et al., 2012; Sarkar et al., 2015; Giosan et al., 2017) (Fig. 1a, b), which dictates the seasonality of the Godavari River, with 75 %–85 % of the annual rainfall and 98 % of the sediment transport in the monsoon/wet season between June and September (Biksham and Subramanian, 1988a, b).

The basin is characterised by several natural gradients, where the upper basin developed on Deccan flood basalts which weathered into clay-rich soils, while the lower basin formed on felsic rock formations with sandy to loamy textured soils (Giosan et al., 2017). Petrogenic OC is absent in the upper basin (Reddy et al., 2021) and very sporadic (i.e. coal deposits) in the lower basin (Usman et al., 2018). In addition, the precipitation gradient ranges from ∼ 430 mm yr⁻¹ in the interior upper basin that is in the rain shadow of the Western Ghats mountain range to ∼ 2300 mm yr⁻¹ near the coast of the Bay of Bengal (Fig. 2). The natural vegetation reflects this gradient and varies from (C₄) grasses, dry deciduous forests and thorny shrublands in the upper basin to moist and evergreen deciduous forests with mostly C₃ flora in the lower basin (Olson et al., 2001; Asouti and Fuller, 2008; Fig. 1b). Agriculture covers ∼ 60 % of the basin, with dominant C₄ crops (sorghum, millet, maize, sugar cane) in the upper basin and rice fields (C₃) in the lower basin (CWC, 2014; Pradhan et al., 2014). Hence, the upper and lower Godavari basin are distinctly different in terms of C₃ / C₄ vegetation distributions, moisture conditions and bedrock geology.

The Godavari basin is divided into five subbasins: the upper (∼ 37 % of the total basin area), middle (6 %) and lower (2 %) Godavari cover the main stem river and are joined by the north (35 %; Wainganga, Penganga, Wardha and Pranhita rivers) and east tributaries (20 %; Indravati and Sabari rivers) (Babar and Kaplay, 2018) (Figs. 1b, S1 in the Supplement). Abundant dams in the upper basin limit the river flow, while a large dam with a reservoir lake at Rajahmundry (Rajamahendravaram) controls the flow into the tidally influenced delta (Pradhan et al., 2014).

Figure 1 (a) Location of the Godavari River basin in peninsular India. (b) Godavari River basin and sampling sites, with the major vegetation zones (Olson et al., 2001; Asouti and Fuller, 2008). Names of the subbasins (grey) and major rivers (blue) are indicated; a zoom for the Godavari delta is available in Fig. S1.

2.2 Sample collection

Samples of aboveground plant material were collected in February–March 2015 (dry season) across the Godavari basin, selecting the three to five most dominant species at each site and spanning the full range of plant life-forms (i.e. trees, shrubs, herbaceous plants and grasses). For shrubs and trees the leaves were collected, and for herbs and grasses the leaves and stems were combined. Given that the deciduous trees and shrubs shed their leaves annually in the dry season (Kushwaha and Singh, 2005; Elliott et al., 2006), leaves were considered the main contributor to soil OC rather than woody biomass. Depending on plant size, each sample consisted of approximately 10 to 50 g of leaves or aboveground plant parts of three to five individuals of the same species or multiple “sun” and “shade” leaves of the same individual. Dominant agricultural crops (e.g. sorghum, maize, millet, sugarcane) were also sampled. After collection, plant samples were air-dried and subsequently frozen upon arrival in the laboratory. In total, 77 samples of C₃ plants and 16 samples of C₄ plants were prepared for analysis.

Topsoils (1–10 cm) were collected during the same campaign (n=47), after removal of the litter layer and by combining three to five spatial replicates. This topsoil layer receives most plant input and gets most likely eroded and transported into the river. SPM (n=40) and riverbed sediments (n=37) were collected in the dry (February–March 2015) and wet (July–August 2015) season. For SPM, surface river water (10–80 L) was collected at a mid-channel position from a bridge or boat or 2–3 m out of the riverbank and filtered on pre-combusted (450 ^∘C, 6 h) GF/F filters (0.7 µm, Whatman) using pressurised steel filtration units (after Galy et al., 2007). Additional river depth profiles (two to three depths, one to three sites across the river) were sampled in the Godavari delta and in the middle Godavari (Fig. 1b; site 10 and 28). At these sites, river water was collected at equal increments to the riverbed with a custom-built depth sampler (after Lupker et al., 2011). Riverbed sediments were dredged at each location with a sediment grabber (Van Veen grab 04.30.01, Eijkelkamp) or with a shovel when the water level was low. The fine fraction (≤ 63 µm) was isolated by sieving for a selection of sites for soils (n=10) and riverbed sediments collected in the wet season (n=25). All samples were frozen upon arrival in the laboratory.

2.3 Elemental and bulk isotopic analysis

Prior to analysis, C₃ and C₄ plants, bulk soils, and riverbed sediments were freeze-dried, homogenised and ground into powder using an agate mortar and pestle or a steel ball mill. The bulk soils and sediments were decalcified by overnight treatment with 1 M HCl and then rinsed twice with deionised water and left to dry at 60 ^∘C, following van Helmond et al. (2017). SPM was decalcified by vapour acidification (Komada et al., 2008; van der Voort et al., 2016). In short, randomly selected, small pieces of GF/F filters containing the SPM were placed in pre-combusted (450 ^∘C, 6 h) Ag capsules and put in a desiccator at 70 ^∘C with 37 % HCl for 72 h and subsequently dried for minimal 120 h with NaOH. Fine-fraction (≤ 63 µm) soils and sediments were placed in pre-combusted Ag capsules and decalcified by addition of 100 µL 1 M HCl and then left to dry overnight at 60 ^∘C, following Vonk et al. (2008, 2010).

Total organic carbon (TOC) content (wt %) and bulk stable carbon isotopic composition (δ¹³C) of plants, (bulk) soils, riverbed sediments and SPM was measured with a Flash 2000 organic element analyser connected to a Thermo Delta V Advantage isotope ratio mass spectrometer (Thermo Scientific, Italy), at NIOZ (Royal Netherlands Institute for Sea Research; Texel, the Netherlands). Total nitrogen (TN) (wt %) was measured for plants and in non-decalcified bulk soils in the same way. Integration was performed with Isodat 3.0 software. TN in non-decalcified bulk sediments was measured with an NA 1500 NCS analyser (Fison Instruments, United Kingdom), at Utrecht University (Utrecht, the Netherlands). Fine-fraction (≤63 µm) soils and riverbed sediments were analysed with an NC2500 elemental analyser coupled to a Thermo Finnigan DeltaPlus isotope ratio mass spectrometer (ThermoQuest, Germany), at Vrije Universiteit Amsterdam (the Netherlands). The results were normalised to certified standards (acetanilide, benzoic acid and urea at NIOZ and USGS40, USGS41 and IAEA-601 at Vrije Universiteit Amsterdam), with an analytical uncertainty <0.1 % for TOC, <0.2 ‰ for δ¹³C and <3 % for TN, based on replicate analysis of standards and samples. The δ¹³C values were reported in the standard delta notation, relative to the international Vienna Pee Dee Belemnite (VPDB) standard for δ¹³C. $\mathrm{C}/\mathrm{N}$ ratios were reported as mass ratios.

2.4 Precipitation and regression analysis

The mean annual precipitation (MAP) in the Godavari basin was used to evaluate the control of drought stress on plant δ¹³C values, as prior studies found evidence for a relationship between MAP and δ¹³C values of C₃ plants around the world (Stewart et al., 1995; Diefendorf et al., 2010; Kohn, 2010). Nonetheless, field surveys and data compilations of C₃ vegetation in drought-stressed regions reported high inter- and intraspecies variation in C₃ plant δ¹³C values in response to MAP (Ma et al., 2012; Liu et al., 2013, 2014; Basu et al., 2021; Luo et al., 2021). The range of ∼ 500 to 1500 mm yr⁻¹ MAP in the Godavari basin was markedly lower than in tropical forests, where MAP is typically >2000 mm yr⁻¹ and where the majority of global C₃ biomass occurs (Kohn, 2010). Here, we focused on MAP in 2014 (Harris et al., 2014; Water Resources Information System, 2021), the growing season preceding the sampling campaign in the dry season in early 2015, considering that the majority (>80 %) of rainfall falls in the wet season and that dry to moist deciduous vegetation is prevalent, which grows new leaves over the wet season and sheds them at the end of the dry period (Kushwaha and Singh, 2005; Elliott et al., 2006). Long-term rainfall deficiencies have resulted in pronounced drought conditions in the upper basin (Kirkels et al., 2020b) (Fig. S2).

In order to deal with the uneven distribution of analysed C₃ plants over the MAP range in the Godavari basin, we used a binning approach for the C₃ plant δ¹³C values. The data were binned by calculating the average and standard error in C₃ plant δ¹³C values per MAP range of 100 mm yr⁻¹. These binned Godavari C₃ plant data were subsequently plotted against the average MAP of each bin and utilised for regression analysis to assess the relation between C₃ plant δ¹³C values and MAP. The correlation established by regression analysis, and the (sub)basin-specific MAP was subsequently used to correct the plant stable carbon isotope end-member values that are used in the C₃ / C₄ mixing model for drought effects. A cut-off value for MAP of 1750 mm yr⁻¹ was used to determine plant end-members without drought effects, as above this MAP the C₃ plant δ¹³C value can be considered constant, as there is no water limitation (Kohn, 2010).

2.5 Mixing model and Suess correction

The relative abundance of C₃ and C₄ plants was estimated based on the isotope mixing model by Phillips and Gregg (2001) using linear mass-balance equations and accounting for the variation in the C₃ and C₄ plants (i.e. sources) as well as in the soils or sediment (i.e. mixture):

$$\begin{array}{}\text{(1)}& {\displaystyle }\mathit{\%}{\mathrm{C}}_{\mathrm{3}}=\left[\right({\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{\mathrm{S}}-{\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{{\mathrm{C}}_{\mathrm{4}}})/({\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{{\mathrm{C}}_{\mathrm{3}}}-{\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{{\mathrm{C}}_{\mathrm{4}}}\left)\right]\cdot \mathrm{100}\mathit{\%},\text{(2)}& {\displaystyle }\mathit{\%}{\mathrm{C}}_{\mathrm{4}}=\mathrm{100}-\mathit{\%}{\mathrm{C}}_{\mathrm{3}}.\end{array}$$

In Eq. (1) ${\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{{\mathrm{C}}_{\mathrm{3}}}$ and ${\mathit{\delta }}^{\mathrm{13}}{\mathrm{C}}_{{\mathrm{C}}_{\mathrm{4}}}$ represent the bulk δ¹³C plant end-member values (‰), and δ¹³C_S is the (sub)basin-specific bulk δ¹³C value (‰) of soil or riverbed sediments. The δ¹³C_S values were concentration-weighted using the TOC content (wt %) of the individual samples in the (sub)basin, and error propagation was accounted for. This mixing model provided an estimation of the proportion of C₃ and C₄ plants, including the standard error of variance on these estimates. Alternative mixing approaches including the C₃ fraction woody cover, which accounts for vegetation structure and shading effects (e.g. Wynn and Bird, 2008; Cerling et al., 2011; Garcin et al., 2014), may be complicated by the fact that agricultural use (∼ 60 % of the basin) and deforestation since the 19th century have resulted in a more open landscape and have drastically reduced the area covered by native, closed-canopy forests, which is now limited to the region of the east tributaries.

C₃ and C₄ plant end-members to resolve the mixing model can be based on measurement of regionally occurring, modern vegetation in the Godavari basin (referred to as Godavari-based or regional end-members), which are representative of the prevailing habitat conditions. Alternatively, global end-members can be used based on C₃ and C₄ plants collected worldwide and reported in literature compilations. Commonly quoted global averages are −27 ‰ for C₃ plants (Cerling et al., 1997; Koch, 1998; Dawson et al., 2002) and −12 ‰ for C₄ plants (Koch, 1998; Dawson et al., 2002). However, the atmospheric δ¹³C of CO₂ has rapidly declined over the past decades due to fossil fuel burning, so δ¹³C data based on analyses of plants in the past require a correction for this so-called Suess effect. Unfortunately, the exact sampling year for these global averages was unknown, but considering that similar values have been reported since the late 1970s we estimated a maximum decrease in global plant δ¹³C values of ∼ 0.9 ‰ (atmospheric CO₂: ∼ −7.5 ‰ in 1978 to ∼ −8.4 ‰ in 2015; Keeling et al., 2001, 2017; Graven et al., 2017; NOAA, 2022). This correction for the Suess effect translated in to estimates of ∼ −27.9 ‰ and −12.9 ‰ for the global averages for C₃ and C₄ plants, respectively (Table 1). Alternatively, Kohn (2010) determined a value of ∼ −28.5 ‰ (corrected for Suess effect for the year 2000) representing global C₃ vegetation including equatorial and mid-latitude biomass, which was updated to a modern value of ∼ −28.9 ‰. To compare δ¹³C values measured in Godavari plants collected in early 2015 with those reported in earlier studies, we updated the latter to account for the Suess effect (Table 1).

Table 1 Suess correction of plant δ¹³C values.

^a Based on Keeling et al. (2001, 2017) and Graven et al. (2017). ^b Measured plants in the Godavari basin are modern.

The plant δ¹³C signal is subsequently transferred to soils or sedimentary deposits, where the δ¹³C_org signal is assumed to integrate long-term and/or spatial areas and thus incorporate/average the plant δ¹³C signal for a range of precipitation within this period/region. In order to compare the δ¹³C_org values of pre-aged soils and sediments with those of modern vegetation, a correction for the Suess effect is warranted. Analysis of Δ¹⁴C of OC in a selection of Godavari soils and sediments by Usman et al. (2018) revealed no distinct differences between the upper and lower basin nor between bulk soils and riverbed sediments. However, the large variation in Δ¹⁴C_OC values, potentially related to small contributions of very old OC from wind-blown coal dust from the open-pit mines in the north of the basin, made it difficult to determine the average age of OC in the Godavari basin based on these data. Assuming OC turnover rates of ∼ 10 years in tropical forest soils to ∼ 25–40 years in savanna ecosystems (Martin et al., 1990; Bird et al., 1996), we estimated an average age of ∼ 30 years for OC in Godavari soils and riverbed sediments. This estimate is at the upper end of recently determined biome-specific OC turnover rates for tropical forests and savannas, where precipitation was shown to have a major effect on soil OC turnover rates (e.g. Carvalhais et al., 2014; Hein et al., 2020). To enable direct comparison of δ¹³C values from plants with δ¹³C_org values from soils and sediments and employ these in the mixing model, we corrected the measured δ¹³C in modern vegetation for the Suess effect to the average age of soil/sediment OC (i.e. 30 years preceding the plant collection in 2015: 1985) (Table 1).

To evaluate the impact of C₄ and C₃ plant end-members on the reconstructed vegetation distribution, we compared four scenarios that included average global and regional vegetation end-members and C₃ plant end-members with and without a correction for a drought-induced ¹³C enrichment (i.e. drought correction). This included (1) Godavari C₄ and C₃ plant end-members, with drought correction of C₃ plant δ¹³C; (2) Godavari C₄ and C₃ end-members, with no drought correction of C₃ plant δ¹³C; (3) C₄ and C₃ end-members based on global data compilations, with the C₃ end-member according to Kohn (2010) that includes equatorial and low-latitude biomass; and (4) global C₄ and C₃ end-members based on global data compilations and commonly quoted in literature. All plant end-members were Suess corrected to the equivalent age of soils and sediments in the Godavari basin (Table 1).

2.6 Statistics

Spatial and seasonal differences were evaluated with (Welch's) one- and two-way ANOVA, (paired) t tests, and non-parametric Mann–Whitney and Kruskal–Wallis tests with the R software package for statistical computing (R4.0.4; RStudio, v. 1.2.5033) and SPSS (IBM, v. 27.0.1.0). The level of significance was p≤0.05. The reported values are the mean ± standard error (SE). Linear regression analysis (Pearson's R) was performed to obtain the correlation between δ¹³C and MAP. Spatial patterns were further investigated with ArcGIS software (Esri, v. 10.8.1).

3 Results and discussion

3.1 Modern C₃ and C₄ plants in the Godavari basin and control by MAP

The Godavari plants (n=96) showed two distinct groups, with bulk δ¹³C values that ranged from −12.7 ‰ to −15.1 ‰ for C₄ plants (n=16, 9 different species) and from −24.3 ‰ to −33.2 ‰ for C₃ plants (n=77, 38 different species) (Kirkels et al., 2022a) (Figs. 2, 3a). The sampled Godavari plants fell within the typical ranges for global C₄ (∼ −10.5 ‰ to −14.5 ‰; Cerling et al., 1997) and C₃ vegetation (∼ −20.5 ‰ to −37.5 ‰; Kohn, 2010) (downward corrected by ∼ −0.5 ‰ for fossil fuel burning). As in the Godavari basin only a very small area was covered by wet evergreen forest with a MAP of ∼ 1500–2000 mm yr⁻¹, we observed less negative δ¹³C values for C₃ plants compared to tropical rain forests with MAP typically exceeding 2000 mm yr⁻¹ (Kohn, 2010) (Figs. 1b, 2).

Figure 2 Map showing the spatial distribution of C₃ and C₄ plant δ¹³C values in the Godavari basin. The red dashed line indicates the boundary of the upper/lower basin. The points refer to the measured δ¹³C values and the 30-year average rainfall distribution (MAP; 0.25^∘, APHRODITE dataset, Asian Precipitation–Highly Resolved Observational Data Integration Towards Evaluation of Water Resources; Yatagai et al., 2012) is shown in the background.

Figure 3 Box-and-whisker plot of δ¹³C_(org) values measured in C₃ and C₄ plants, soils, SPM, and riverbed sediments collected in the dry and wet season in the upper and lower Godavari basin. The box represents the first (Q₁) and third (Q₃) quartiles, and the line in the box represents the median value; the whiskers extend to $\mathrm{1.5}\cdot (Q_{\mathrm{3}}-Q_{\mathrm{1}})$ values, and outliers are shown as points. Outcomes of the two-way ANOVA are indicated for plants, soils and SPM, and those of non-parametric tests (Mann–Whitney and Kruskal–Wallis) are for the sediments. The level of significance is as follows: NS for not significant, ^* for p≤0.05, ${}^{**}$ for p≤0.01 and ${}^{***}$ for p≤0.001.

The C₄ crops we collected in the Godavari basin, including Zea mays, Sorghum vulgare and Saccharum officinarum, had similar δ¹³C values as earlier reported for these species in the Godavari basin (Pradhan et al., 2014; Krishna et al., 2015). The Godavari C₄ plants we sampled had on average more negative δ¹³C values than those collected in the only other extensive field survey of Indian plants on the Gangetic Plain ($-\mathrm{14.0}\pm \mathrm{0.2}$ ‰ (± standard error: SE), n=16, (Fig. 3a) vs. $-\mathrm{12.7}\pm \mathrm{0.2}$ ‰, n=45; p≤0.001), where they found most negative δ¹³C values in areas with MAP < 1000 mm yr⁻¹ and observed an effect of MAP on C₄ plant δ¹³C values (Basu et al., 2015). In contrast, the Godavari C₄ plants showed no significant correlation with MAP (Eq. 3; Pearson's $R=-\mathrm{0.10}$; p=0.70) (Fig. 4). The absence of a correlation for the Godavari C₄ plants may be influenced by the relatively small sample size, and their main occurrence is only a limited part of the MAP range covered (i.e. ∼ 500–900 mm yr⁻¹ in 2014). Nevertheless, earlier studies also predominantly found no trends in C₄ plant δ¹³C values in response to MAP in dry ecosystems around the globe (<800 mm yr⁻¹; Schulze et al., 1996; Swap et al., 2004). This finding was attributed to the CO₂-concentrating mechanism in C₄ plants. This adaption to water loss due to evaporation in warm and dry climates may be influenced by the leaking of CO₂ from bundle sheath cells during extreme drought but functions relatively robustly for a wide range of environmental conditions, including drought stress (Murphy and Bowman, 2009). Taken together, we interpret that there was no basis for correction of the δ¹³C C₄ end-member for drought conditions.

Figure 4 Regression analysis of δ¹³C values against MAP (2014, previous growing season) for C₃ and C₄ plants in the Godavari basin and estimates based on global C₃ vegetation models by Diefendorf et al. (2010) and Kohn (2010). For Godavari C₃ plants, δ¹³C values are also binned per MAP of 100 mm yr⁻¹, and the mean ± standard error (SE; whiskers) is presented. The bin of 500–600 mm yr⁻¹ includes two samples with a MAP of 489 mm. The solid lines denote the linear fit for the Diefendorf et al. (2010), Kohn (2010) and C₃ binned correlation. Pearson's R, equation and R² are given for each correlation. Compilation-based δ¹³C values of vegetation sampled around the world are given for C₄ (−12.9 ‰) (yellow dashed) and C₃ plants (see Table 1). The latter is represented by a commonly quoted δ¹³C value of −27.9 ‰ (dark green, dotted) that is biased toward dry ecosystems according to Kohn (2010) and a δ¹³C value established by Kohn (2010) of −28.9 ‰ (dark green, dashed) that includes more equatorial and mid-latitude C₃ biomass. The level of significance is as follows: NS for not significant, ^* for p≤0.05, ${}^{**}$ for p≤0.01 and ${}^{***}$ for p≤0.001.

The average δ¹³C value of Godavari C₃ plants was significantly different in the upper and lower basin, where values are slightly less negative in the upper basin ($-\mathrm{28.0}\pm \mathrm{0.3}$ ‰, n=32) than in the lower basin ($-\mathrm{28.8}\pm \mathrm{0.2}$ ‰, n=45; p≤0.05) (Fig. 3a). This finding corresponds to the observed spatial gradient in MAP in the Godavari basin and thus suggests an effect of MAP on C₃ plant δ¹³C values. Indeed, long-term MAP (1901–2015) was markedly lower in the upper than the lower Godavari basin (p≤0.001), and this contrast became more extreme for the 5-year average and 2014 MAP, which resulted in drought conditions in the upper basin (Kirkels et al., 2020b) (Fig. S2). The least negative subbasin-averaged δ¹³C value for C₃ plants was found in the upper Godavari subbasin that received the least precipitation (δ¹³C: $-\mathrm{28.0}\pm \mathrm{0.3}$ ‰, n=30; MAP: 593±18 mm yr⁻¹), compared to the most negative value in the east tributaries that received significantly more precipitation (δ¹³C: $-\mathrm{30.1}\pm \mathrm{0.5}$ ‰, n=5; p≤0.05; MAP: 1530±142 mm yr⁻¹; p≤0.01) (Figs. 2, S2). The east tributaries are the only part of the Godavari basin that is covered by native, wet to moist forests where a denser canopy caused ample shading. This likely resulted in lower soil temperatures and higher moisture and humidity levels in the understorey, which generally favours C₃ vegetation (Cerling et al., 2011), which was indeed exclusively found in this Godavari subbasin (Figs. 1b, 2). Likewise, Garcin et al. (2014) reported very negative stable carbon isotopic values for dense tropical forests in Cameroon (>80 % tree cover), mainly controlled by water availability, although a “canopy effect” (van der Merwe and Medina, 1991), which involves recycling of ¹³C-depleted CO₂ in the understorey of closed-canopy forests and stable carbon isotope fractionation due to photosynthesis under low-light conditions, may have resulted in additional ¹³C-depletion in the C₃ leaves.

The individual Godavari C₃ plants revealed a small but significant effect by MAP on their δ¹³C values (Eq. 4; Pearson's $R=-\mathrm{0.34}$; p≤0.01) (Fig. 4). This finding supports earlier studies by Diefendorf et al. (2010) and Kohn (2010) that found a strong control by MAP on C₃ plant δ¹³C on a global scale and established quantified relationships between the fractionation of carbon isotopes and environmental conditions, including MAP. Application of these established relations for the Godavari basin revealed very similar trends as the Godavari C₃ plants (Eqs. 5 and 6, respectively; Fig. 4), although the Kohn (2010) relation was most similar in terms of slope and had a higher R² than the Diefendorf et al. (2010) correlation. For the individual Godavari C₃ plants, we noted considerable variation in δ¹³C values for any certain amount of precipitation, in line with earlier studies that found high inter- and intraspecies variation in C₃ plant δ¹³C values in response to MAP (Ma et al., 2012; Liu et al., 2013, 2014; Basu et al., 2021; Luo et al., 2021). Moreover, the individual Godavari C₃ plants were not evenly distributed over the entire precipitation range, making it more difficult to establish a correlation. Together, this resulted in a relatively weak linear correlation between MAP and individually measured C₃ plants (Eq. 4; R²=0.12), where MAP explained only ∼ 12 % of the variation in C₃ plant δ¹³C values.

Subsequent binning of C₃ plant δ¹³C values to overcome their uneven distribution over the range of MAP revealed a strong and significant correlation with MAP (Eq. 7; Pearson's $R=-\mathrm{0.90}$; p≤0.01) (Fig. 4). The slope of this binned C₃ plant correlation (i.e. −0.18 ‰ per 100 mm MAP) could be used to estimate the offset of measured plant δ¹³C values to those expected as a function of MAP. For the binned C₃ plants, MAP explained ∼ 82 % of the variation in δ¹³C values. This linear relation applies to the interval of precipitation of ∼ 500 to 1750 mm yr⁻¹, above which the C₃ plant δ¹³C value is assumed to be constant, as there is no water limitation (Kohn, 2010). For very dry ecosystems with MAP < 500 mm yr⁻¹, non-linear effects on C₃ plant δ¹³C values need to be considered due to extreme drought stress (Kohn, 2010, 2011; Freeman et al., 2011). Similar to the Godavari C₃ plants, a meta-analysis for low-latitude regions (11–30^∘ N) showed that the average C₃ plant δ¹³C value would change by ∼ −0.2 ‰ for every 100 mm increase in MAP for the interval of 500–1500 mm yr⁻¹ (Basu et al., 2019b). We note that although the slope of the binned Godavari C₃ plants was similar to that of the Kohn (2010) correlation, the ∼ 0.6 ‰ offset in intercept suggests consistently more negative δ¹³C values for Indian C₃ vegetation in relation to MAP than established based on global data compilations by Kohn (2010) and Diefendorf et al. (2010). Within the MAP interval of ∼ 1000–1500 mm yr⁻¹, Godavari C₃ plants had similar δ¹³C values ($-\mathrm{29.0}\pm \mathrm{0.4}$ ‰, n=11; p>0.20) to those collected on the Gangetic Plain situated in the Himalayan foreland ($-\mathrm{29.6}\pm \mathrm{0.2}$ ‰, n=76; Basu et al., 2015). The strong control by MAP on C₃ plant δ¹³C values we established for the Godavari basin suggests that this factor needs to be considered in mixing-model approaches to reconstruct vegetation in regions influenced by the Indian monsoon.

3.2 Tracing the plant δ¹³C signal along the plant–soil–river continuum

The bulk δ¹³C_org signal preserved in soil or sedimentary archives depends on the input and integration of the C₃- and C₄-plant-derived δ¹³C signal, where a temporal shift of ∼ 0.8 ‰ due to the Suess effect is considered regarding the ∼ 30-year turnover rate of OC in soils and sediments (Martin et al., 1990; Bird et al., 1996; Usman et al., 2018) (Table 1). Regardless, δ¹³C_org may also be influenced by hydroclimatic controls on OC degradation and stabilisation mechanisms in soils and in the river (e.g. Carvalhais et al., 2014; Ward et al., 2017; Hein et al., 2020; Eglinton et al., 2021). The river-transported OC, in the form of suspended (SPM) or riverbed sediments, may be a complex mixture depending on soil- and plant-/litter-derived OC sourcing from particular parts of the basin following rainfall distributions, soil mobilisation, aquatic primary production and hydrodynamic sorting processes within the river. This complexity warrants further exploration of the evolution and provenance of the δ¹³C_org signal along the plant–soil–river continuum in the Godavari basin.

3.2.1 Soils

The Godavari (bulk) soils had on average less negative δ¹³C_org values in the upper than in the lower basin ($-\mathrm{21.4}\pm \mathrm{0.5}$ ‰, n=22 vs. $-\mathrm{23.5}\pm \mathrm{0.5}$ ‰, n=25; p≤0.01) (Kirkels et al., 2022a) (Figs. 3b, 5b). This isotopic contrast corresponds with the vegetation distribution in the basin, with mixed C₃ and C₄ vegetation in the upper basin and more C₃ plants in the lower basin (Figs. 1b, 2). The least negative δ¹³C_org values were found in soils in the upper Godavari and north tributaries ($-\mathrm{21.3}\pm \mathrm{0.5}$ ‰, n=20 and $-\mathrm{22.0}\pm \mathrm{0.7}$ ‰, n=12, respectively) covered by thorny shrublands, dry deciduous forest and predominantly C₄ crops, followed by the middle Godavari ($-\mathrm{23.5}\pm \mathrm{0.5}$ ‰, n=4) in a transition zone, and most negative δ¹³C_org values were found in soils in the east tributaries and lower Godavari that were covered by moist/evergreen forests and C₃ crops ($-\mathrm{24.7}\pm \mathrm{0.6}$ ‰, n=6 and $-\mathrm{25.1}\pm \mathrm{0.6}$ ‰, n=5, respectively) (Figs. 1b, 5a, b, S3). These findings correspond with the general observation that the majority of the soil organic carbon derives from microbially processed plant residues, while the microbial biomass itself has been estimated to contribute only 1 %–5 % (Kögel-Knabner, 2002; Simpson et al., 2007). We note that it is challenging to determine the actual size of the microbial biomass, which is highly dependent on prevailing moisture levels and the availability of easy degradable carbon as an energy source and has a high spatial heterogeneity (Birge et al., 2015; Wiesmeier et al., 2019). Given our sampling strategy where soils were collected in the dry season, the low moisture levels likely limited the microbial biomass size and activity in the Godavari soils. Hence, the δ¹³C_org values in Godavari soils can be interpreted as a time-averaged plant signal on a decadal scale that reflects the long-term hydrological conditions that underlie this vegetation distribution.

Figure 5 Maps showing the spatial distribution of (bulk) δ¹³C_(org) values in the Godavari basin for (a) C₃ and C₄ plants, (b) (bulk) soils, SPM collected in the (c) dry and (d) wet season, and (bulk) riverbed sediments collected in the (e) dry and (f) wet season. The points refer to the measured δ¹³C_(org) values, and the 30-year average rainfall distribution is shown in the background.

A potential degradation-related ¹³C enrichment of the plant-derived OC in soils may be suggested by the ∼ 4 ‰ difference between Suess-corrected C₃ plant δ¹³C and soil δ¹³C_org values found in the C₃-dominated lower basin (Fig. 3b, Table 1). This offset is relatively large compared to typical 1 ‰–3 ‰ higher δ¹³C_org values due to soil OC degradation (e.g. Krull et al., 2005; Sreemany and Bera, 2020), which may result from preferential mineralisation of chemical compounds containing less ¹³C (i.e. with a lower δ¹³C value) and/or stable carbon isotope fractionation during microbial processing preferentially utilising ¹²C over ¹³C (Schmidt and Gleixner, 1998). Instead of a degradation-induced shift, deforestation since the late 19th century and agricultural expansion with predominantly drought-adapted C₄ crops (Ponton et al., 2012; CWC, 2014; Pradhan et al., 2014) may have contributed to relatively more C₄ input and thus enrichment of topsoil δ¹³C_org values. Regardless, soil OC degradation and stabilisation processes in tropical to subtropical biomes may have an opposite effect on soil δ¹³C_org signals. In mixed C₃ / C₄ ecosystems, C₄-plant-derived OC has been shown to contain more labile compounds and thus degrade more rapidly than C₃-plant-derived OC that contains compounds that are more difficult to degrade (Wynn, 2007; Wynn and Bird, 2007). However, C₄-derived OC has also been shown to be preferentially incorporated into fine fractions, as fine particles are presumed to have a higher ability to stabilise the labile, C₄-derived compounds onto mineral surfaces where they are better protected against degradation, whereas C₃-derived OC is preferentially added to the coarse fraction, thus leaving it less protected (Bird and Pousai, 1997; Wynn, 2007; Wynn and Bird, 2007). In the Godavari basin, Usman et al. (2018) reported similar Δ¹⁴C values for soil OC in the upper and lower basin that have different C₃ and C₄ plant covers, suggesting that the net effect of preferential degradation (more young OC) and stabilisation (more old OC) is minor. Indeed, extensive degradation of C₄-plant-derived OC is unlikely, given that the upper basin with most C₄ plants contains clay-rich, fine particles from weathering of the Deccan basalts (Giosan et al., 2017; Usman et al., 2018; Kirkels et al., 2022b), which would contribute to stabilising the C₄-derived OC. Fine soils (≤63 µm) sampled in the upper basin had on average slightly less negative δ¹³C_org values than the bulk soils ($-\mathrm{19.7}\pm \mathrm{1.0}$ ‰, n=8 vs. $-\mathrm{21.1}\pm \mathrm{0.7}$ ‰, n=8; p≤0.09), pointing towards preferential stabilisation of C₄-derived OC in the fine fraction (Fig. S4). This result highlights that differences between bulk and fine size fractions could be important for the δ¹³C_org signal preserved in soils.

3.2.2 Suspended particulate matter

For SPM collected in the Godavari River, bulk δ¹³C_org values were consistently more negative in the dry than in the wet season and more negative in the upper than the lower basin, but there was no significant interaction between the effects of the season and upper/lower basin (Fig. 3c). In the dry season, SPM δ¹³C_org values were significantly more negative than in Godavari soils ($-\mathrm{27.8}\pm \mathrm{0.3}$ ‰, n=40 vs. $-\mathrm{22.5}\pm \mathrm{0.4}$ ‰, n=47; p≤0.001), making a predominant contribution of soil material to the river in this season unlikely (Figs. 3b, c, 5b, c). Instead, the quiescent waters behind dams and very low river discharge favour aquatic primary production (Pradhan et al., 2014). Freshwater phytoplanktonic matter usually has relatively low δ¹³C_org values. The isotopic fractionation between phytoplankton and dissolved inorganic carbon (DIC) has been estimated at ∼ −23 ‰, resulting in typical phytoplankton-derived δ¹³C_org values between −31 ‰ and −35 ‰ in the Ganges–Brahmaputra as well as at the start of the dry season in the Godavari River (e.g. Aucour et al., 2006; Galy et al., 2008b; Krishna et al., 2015). The observation of slightly less negative δ¹³C_org values in our SPM collected at the end of the dry season may be explained by eutrophic conditions due to agricultural/wastewater inputs which fuelled intense aquatic production in the Godavari River (Balakrishna and Probst, 2005; Pradhan et al., 2014). During periods of high aquatic productivity, the carbon isotope fractionation becomes smaller (up to 0 ‰; Torres et al., 2012), and more ¹³C gets incorporated into the phytoplanktonic biomass. Aquatic primary production is also supported by the strong increase in %OC from soils to dry-season SPM (0.8 ± 0.1 %, n=47 vs. 11.4±1.1 %, n=39; p≤0.001), since phytoplankton-derived SPM is typically high in %OC (Aucour et al., 2006; Galy et al., 2008b). Notably, δ¹³C_org values of dry-season SPM became less negative near the Godavari's outflows into the Bay of Bengal (Figs. 5c, S3), suggesting mixing of freshwater and estuarine/marine phytoplankton in the delta, where the latter has typically less negative δ¹³C_org values (i.e. −22.8 ‰ to −24.4 ‰; Dehairs et al., 2000; Krishna et al., 2015; Gawade et al., 2018). This observation is consistent with changes in electrical conductivity and water isotopic values (δ¹⁸O) that showed seawater intrusion in the delta in the dry season (Kirkels et al., 2020b). Regardless, mixing of riverine and marine OC with different stable carbon isotopic values at the outflow complicates the tracing of the Godavari-derived OC signal from the river mouth to marine sedimentary deposits.

In the wet season, there was a strong isotopic contrast between SPM collected in the upper and lower basin ($-\mathrm{26.4}\pm \mathrm{0.8}$ ‰, n=14 vs. $-\mathrm{25.0}\pm \mathrm{0.3}$ ‰, n=26; p≤0.001) (Fig. 3c). The negative δ¹³C_org values in the upper basin suggest continuous aquatic production, allowed by the limited rainfall and abundant dams in this region that created standing waters and facilitated year-round aquatic productivity (Pradhan et al., 2014; Kirkels et al., 2020b). A few sites in the upper basin had remarkably less negative δ¹³C_org values, suggesting that some local soil or C₄ plant input occurred at locations where agricultural fields with exposed topsoils were situated next to the river (Fig. 5d). In the lower basin, wet-season SPM δ¹³C_org values varied from −21.3 ‰ to −27.1 ‰ and fell within the range of bulk soils in this region (−18.7 ‰ to −27.1 ‰). This resemblance suggests contribution of soil-derived OC to wet-season SPM in the lower basin. Water isotopic values and rainfall distributions confirmed substantial discharge in the wet season and identified the Weiganga and Pranhita rivers in the north tributaries and the Indravati River in the subbasin of the east tributaries as major source areas (Kirkels et al., 2020b). It is generally assumed that the high flow velocity and turbidity in the wet season would limit aquatic production as well as OC degradation during fluvial transport (Balakrishna and Probst, 2005; Acharyya et al., 2012). In contrast, Galy et al. (2008b) showed for the Ganges–Brahmaputra that ∼ 50 % of the wet-season SPM derived from the upper reaches of the basin was degraded during river transit and replaced by local input from the plains. Quantifying the extent of OC degradation during fluvial transport is thus not straightforward, as it depends on the stability of OC in the river, e.g. protected by mineral associations or not and on its residence time in a specific basin (e.g. Ward et al., 2017; Eglinton et al., 2021).

3.2.3 Riverbed sediments

For riverbed sediments collected in the Godavari River, bulk δ¹³C_org values were consistently less negative in the upper than in the lower basin ($-\mathrm{24.4}\pm \mathrm{0.6}$ ‰, n=21 vs. $-\mathrm{26.6}\pm \mathrm{0.2}$ ‰, n=53; p≤0.001), but there was no significant seasonal effect (Figs. 3d; 5e, f). The latter suggests that riverbed sediments represented a season-integrated δ¹³C_org signal vs. SPM that showed more seasonal variation. A Kruskal–Wallis test revealed that sediments collected in the upper basin in the dry season and those collected in the lower basin in the wet season were significantly different ($-\mathrm{24.4}\pm \mathrm{0.7}$ ‰, n=13 vs. $-\mathrm{26.6}\pm \mathrm{0.3}$ ‰, n=29; (Bonferroni-corrected) p≤0.05) (Fig. 3d). However, riverbed sediments showed a high spatial variation in δ¹³C_org values, which ranged from −20.0 ‰ to −29.0 ‰ in the upper and −23.3 ‰ to −28.7 ‰ in the lower basin (Fig. 5e, f). The strong influence of location (upper or lower basin) on sediment δ¹³C_org values corresponds to the vegetation distribution with more C₄ plants and drought-stressed C₃ plants in the upper basin, leading to less negative plant δ¹³C and soil δ¹³C_org values being transferred to the riverbed sediments. Riverbed sediment δ¹³C_org values were consistently more negative than soil δ¹³C_org values in the upper and lower basin (Fig. 3b, d), suggesting input from an additional, ¹³C-depleted source, likely phytoplankton- or C₃-plant-derived OC.

Possible contributions of additional sources can be further explored using a source diagram, where the relation between δ¹³C_org values and $\mathrm{C}/\mathrm{N}$ ratios reveals source-specific distributions which help to identify the provenance of OC sources in aquatic ecosystems (e.g. Lamb et al., 2006). For the Godavari basin, riverbed sediments plotted close to Godavari C₃ plants, albeit at lower $\mathrm{C}/\mathrm{N}$ ratios (Fig. 6). This may suggest that slightly degraded C₃-plant-derived OC is selectively transported by/stored in the sediments of the lower basin. The lower $\mathrm{C}/\mathrm{N}$ ratios suggest slight degradation of this OC, although there could also be a small contribution of phytoplankton-derived OC that settled onto the riverbed. Notably, the sediments of the upper basin plotted generally closer to the soils of the upper basin, suggesting an inherited soil δ¹³C_org signal in these sediments, as soil input was diminished by limited rainfall in this region during our sampling campaigns.

Figure 6 Source plot with relation between δ¹³C_org values and $\mathrm{C}/\mathrm{N}$ mass ratios for modern C₃ and C₄ plants and (bulk) soils and riverbed sediments collected in the dry and wet season in the Godavari basin. $\mathrm{C}/\mathrm{N}$ ratios are retrieved from Kirkels et al. (2021). Freshwater phytoplankton is defined for Indian monsoonal rivers in the dry season (solid line) ($\mathrm{C}/\mathrm{N}$: 2.5–8; data from Balakrishna and Probst, 2005; δ¹³C_org: −31 ‰ to −35 ‰; data from Aucour et al., 2006; Galy et al., 2008b; Krishna et al., 2015) and based on SPM collected at the end of the dry season in this study (dashed pink line). For C₃ and C₄ plants, Suess-corrected δ¹³C values (equivalent to the age of soils and sediments; see Table 1) are indicated by the dotted lines.

3.2.4 Hydrodynamic sorting

Fluvial transport of fine particles is of particular interest, as they are generally enriched in OC and their size facilitates effective offshore transport to sedimentary deposits (e.g. Bianchi et al., 2018). At the same time, coarse- and fine-grained particles have been shown to be sensitive to hydrodynamic sorting effects within the river, resulting in depth-specific OC distributions and/or preferential transport of certain OC components (Galy et al., 2008b; Bouchez et al., 2014; Feng et al., 2016; Repasch et al., 2022). Similar to the Godavari soils, riverbed sediments collected in the wet season revealed distinctly less negative δ¹³C_org values for the fine fractions than the corresponding bulk sediments ($-\mathrm{24.4}\pm \mathrm{0.3}$ ‰, n=25 vs. $-\mathrm{25.7}\pm \mathrm{0.4}$ ‰, n=25; p≤0.001) (Figs. 5f, S4). This difference corresponds to previous findings that C₄-plant-derived OC is preferentially associated with finer fractions, whereas C₃-plant-derived OC is contained in the coarser fraction of soils (Bird and Pousai, 1997; Wynn and Bird, 2007) and river(-dominated) sediments (e.g. France-Lanord and Derry, 1994; Bianchi et al., 2002). This finding suggests differential transport of C₃- and C₄-derived OC in the fine and bulk Godavari riverbed sediments.

Furthermore, wet-season SPM collected along river depth profiles in the middle and lower Godavari showed in general more negative δ¹³C_org values with depth (Fig. S5). Also, the riverbed sediments dredged at each location had more negative δ¹³C_org values than the SPM in the water column above. This trend in δ¹³C_org values suggests that C₃-derived OC with typically more negative δ¹³C_org values was transported in coarse-grained sediments and SPM near the riverbed, similar to findings at peak discharge for the Himalayan-derived Ganges–Brahmaputra (Galy et al., 2008b) and the Bermejo River draining the central Andes (Repasch et al., 2022). Alternatively, a ¹³C-depleted phytoplankton source was unlikely given that the high flow velocity and turbidity in the monsoon season prevent light penetration and generally limit algae production (Balakrishna and Probst, 2005; Acharyya et al., 2012). Indeed, analysis of lignin, a macromolecule that is exclusively produced by plants, revealed a C₃-derived “woody undercurrent” in Godavari riverbed sediments collected in the wet season in the lower basin (Pradhan et al., 2014), similar to earlier findings in the Madre de Dios and Mississippi rivers with strong seasonal variability in their hydrology (Bianchi et al., 2002; Feng et al., 2016). This information corroborates differential transport of fine- and coarse-grained particles and their associated OC sources in the Godavari basin. This insight has important implications for the δ¹³C_org signal that is finally exported to marine sedimentary deposits in front of the Godavari's mouth, where C₄-derived OC with less negative δ¹³C_org values may be overrepresented, as fine fractions are transported farther offshore and in greater quantities than coarse-grained particles (Goñi et al., 1997, 1998; Bianchi et al., 2002). Indeed, Holocene marine sediments collected in front of the Godavari's mouth contained no woody particles (Ponton et al., 2012; Giosan et al., 2017; Usman et al., 2018), in contrast to the Bengal Fan fed by the Himalayan-derived Ganges–Brahmaputra, which covers a steep altitudinal gradient and carries coarse sediments far offshore at high-flow conditions, where wood particles were found in sediments spanning the last 19 Myr (Lee et al., 2019).

3.3 C₃ / C₄ vegetation reconstruction and end-member analysis

Mixing-model estimates based on the δ¹³C_org signal in Godavari soils revealed a high %C₄ plants in the upper Godavari subbasin (43±3 %; scenario 1) and in general in the upper basin (43±3 %), which decreased toward the lower basin (31±3 %), from the north tributaries (40±5 %), middle Godavari (32±4%) and east tributaries (29±4 %) to the lowest %C₄ plants in the lower Godavari subbasin (19±4 %) (Fig. 7a). This result is consistent with the vegetation structure in the Godavari basin, with a mixed C₃ / C₄ vegetation in the dry upper basin and dominant C₃ vegetation in the lower basin (Figs. 1b, 2) (Giosan et al., 2017; Usman et al., 2018). Pollen assemblages of the modern vegetation corroborate these findings, as they found ∼ 50 % Poaceae (here C₄ grasses) in the upper basin based on surface sediments of Lonar Lake (Prasad et al., 2014; Riedel et al., 2015) and ∼ 30 % Poaceae for the whole Godavari basin based on marine surface sediments in front of the Godavari mouth in the Bay of Bengal (Zorzi et al., 2015), which was close to our estimate of 37±2 % C₄ plants for the whole basin. Notably, the estimated %C₄ plants in the C₃-dominated lower basin was relatively high (31±3 %) and may have been influenced by degradation-related enrichment of the soil δ¹³C_org signal leading to a potential overestimation of the %C₄ plants. The C₃ / C₄ abundances based on the δ¹³C_org signal of riverbed sediments collected in the wet season also reflected the vegetation gradient (Figs. 7b, 1b). The estimated %C₄ plants was consistently lower for riverbed sediments than for soils (−1 % to −21 %; scenario 1), particularly in the north tributaries (−19 %) and the middle Godavari (−21 %) subbasins. Hydrodynamic sorting and preferential transport of C₃-plant-derived OC in the coarse-grained riverbed sediments may have influenced the C₃ / C₄ plant estimates, resulting in a potential underestimation of the %C₄ plants in these basins.

Figure 7 The estimated proportion (%) of C₄ plants determined by a linear mixing model (by Phillips and Gregg, 2001) vs. MAP for different (sub)basins, based the OC-weighted δ¹³C_org signal in (a) soils and (b) riverbed sediments collected in the wet season. The whiskers represent ± standard error (SE). Four scenarios are included: (1) Godavari C₄ and C₃ plant end-members, with drought correction of C₃ plant δ¹³C; (2) Godavari C₄ and C₃ end-members, with no drought correction of C₃ plant δ¹³C; (3) C₄ and C₃ (based on Kohn, 2010) plant end-members based on compilations of global vegetation; and (4) C₄ and C₃ plant end-members as commonly quoted in literature (e.g. Cerling et al., 1997; Koch, 1998; Dawson et al., 2002) based on global vegetation (see Table 1). For scenario (2) the mean δ¹³C value for Godavari C₃ plants with no drought stress (i.e. −29.3 ‰) was determined at a MAP of 1750 mm yr⁻¹ using Eq. 7 (see Fig. 4).

In order to assess the influence of drought stress and the use of region-specific plant δ¹³C end-members on C₃ / C₄ vegetation mixing-model estimates, we compared four scenarios (Fig. 7). Comparing end-members with and without drought correction of the Godavari C₃ plant δ¹³C end-member (scenario 1 and 2, respectively) revealed a differential effect of MAP on the estimated %C₄ plants. For MAP above 1500 mm yr⁻¹ when water was probably no limiting factor for C₃ plants (Kohn, 2010), the estimated %C₄ plants was ∼ 2 % lower without correction. But for MAP below 1500 mm yr⁻¹, the estimated %C₄ plants was 2 %–6 % higher without this correction for MAP, suggesting a potential overestimation of the %C₄ plants. Notably, drought correction of the C₃ plant δ¹³C end-member is increasingly important with decreasing MAP.

Comparing mixing-model results for Godavari-specific end-members with drought correction and global average end-members for C₄ plants (−12.1 ‰) and C₃ plants (−28.1 ‰; based on Kohn, 2010) (scenario 1 and 3, respectively) revealed that the estimated %C₄ plants was 1 %–9 % lower for the latter, suggesting a potential, minor underestimation. The estimated %C₄ plants was 5 %–9 % lower for the wetter north and east tributaries, while the difference was minor (≤ 5 % lower) for the other subbasins as well as for the whole, upper and lower basin. Notably, comparing Godavari-specific end-members with drought correction and commonly quoted global end-members for C₄ (−12.1 ‰) and C₃ plants (−27.1 ‰) (scenario 1 and 4, respectively) showed that the estimated %C₄ plants was 7 %–15 % lower for the latter, suggesting a substantial underestimation when using these global end-members. This would result in a %C₄ plants of only 36±3 % for the upper basin, where C₄ grasses and crops are most abundant (Pradhan et al., 2014; Prasad et al., 2014; Riedel et al., 2015; Giosan et al., 2017) and only 21±3 % C₄ plants in the C₃-dominated lower basin. Finally, the underestimation of %C₄ plants was even larger (11 %–14 %) for all (sub)basins when comparing the latter global end-members with Godavari-specific end-members without drought correction (scenario 4 and 2, respectively). Our results from the Godavari basin highlight how important it is to make informed choices about the plant end-members used for vegetation reconstructions, considering characteristics of the regional vegetation and environmental controls such as MAP in regions where moisture limiting conditions may prevail (on a seasonal scale).

4 Conclusion – implications for vegetation reconstructions

Our analysis of contemporary C₃ and C₄ plants, soils, and sediments from the Godavari basin has resulted in three important considerations for the reconstruction of C₃ / C₄ vegetation distributions using bulk δ¹³C analysis. Firstly, our C₃ and C₄ plant δ¹³C data for the Godavari revealed the importance of making informed choices about the plant end-members used for vegetation reconstructions, considering the regional conditions including vegetation species, structure, land use and environmental controls. In particular in tropical to subtropical regions with generally high vegetation biomass and more negative δ¹³C C₃ plant values, the use of end-members representing these negative C₃ plant δ¹³C values would result in a higher estimation of %C₄ plants when applied to areas with less precipitation (MAP < 2000 mm), such as the Godavari basin. Given that extensive datasets of the stable carbon isotopic values of commonly occurring plants in certain areas are not generally available, an alternative approach is to derive carbon stable isotope plant end-member values from data compilations of global vegetation, with consideration for the conditions prevailing in the area of interest. For the Godavari basin, using plant end-members based on the measurement of Godavari plants resulted in a ∼ 1 %–9 % higher estimated %C₄ plants compared to using a compilation-based end-member value from Kohn (2010) that was considered representative of the vegetation in this region (Table 2). Correction for the Suess effect to account for fossil fuel burning over the last decades had a relatively minor impact (generally <2 %, depending on the plant end-member chosen) on the estimated %C₄ plants (Table 2). The offset we observed for the Godavari basin can thus be attributed to the choice of plant end-members used in the mixing model.

Secondly, it is well-known that C₃ plant δ¹³C values become less negative with decreasing amounts of precipitation, which may result in an underestimation of the %C₄ plants. Incorporating a drought correction for the C₃ plant end-member resulted in maximally 6 % lower %C₄ plants for the different Godavari subbasins, but the extent of this effect will be most pronounced in C₃-dominated ecosystems with low MAP. For the Godavari basin, we established that the impact of drought correction on the estimated %C₄ plants differed depending on the amount of MAP (Table 2). The effect was relatively minor (∼ 2 % lower) for C₃-dominated areas that received >1500 mm yr⁻¹ and where water was no limiting factor but became positive (∼ 2 % higher) for areas with 1000–1500 mm yr⁻¹ rainfall and even larger (4 %–6 %) for areas with a MAP of ∼ 500–1000 mm yr⁻¹ and mixed C₃ and C₄ vegetation.

Table 2 Summary of processes/factors influencing the C₃ / C₄ mixing-model estimates.

^a With/without drought correction of the regional end-member. ^b Equivalent to the OC turnover in soils and sediments (i.e. 1985). ^c C₃ plant (global) end-member of −27.1 ‰; Godavari-based or −28.1 ‰ (global) end-member.

Thirdly, C₃ / C₄ vegetation reconstructions based on the δ¹³C_org signal preserved in soils, riverine and/or marine sediments should take potential alterations in this signal along the plant–soil–river continuum into account. Differential degradation and stabilisation of C₃- and C₄-derived OC in soils and sediments as well as phytoplankton contributions in the river and hydrodynamic sorting effects that separate coarse and fine fractions and their associated OC may result in a complex interplay affecting the mixing-model estimates (Table 2). For the Godavari basin, our results suggest that the δ¹³C values of C₃ and C₄ plants are relatively well-maintained throughout the transport chain, so δ¹³C_org values of soils and riverbed sediments can be used as proxy for C₃ / C₄ vegetation distributions.

Finally, our results in the Godavari basin highlight that precipitation plays an impactful role in the reconstruction of vegetation by influencing plant δ¹³C end-members and OC decomposition rates and driving OC transport by rivers. Our results thereby support the assumption that in monsoon-influenced regions, wetter conditions/periods require no drought correction of the C₃ plant end-member and generate more erosion and rapid transport downstream, which generally limits degradation during river transit. This implies that the δ¹³C_org signal exported to sedimentary deposits would reflect the C₃ / C₄ vegetation distribution in the basin.