Study area and field sampling<\/p>\n
Spanning about 30% of the entire QTP, this 780,000\u2009km2 study area is situated between the boundary of the Loess Plateau (Gansu Province, 103\u00b0 E) and the border of China (Tibet Autonomous Region, 78\u00b0 E), and includes the headwater of eight major river basins: the Yellow, Yangtze, Lancang\u2013Mekong, Nu\u2013Salween, Derung (Irrawaddy), Yarlung Tsangpo\u2013Brahmaputra, Marja Tsangpo (Ganges) and Indus Rivers (Fig. 1a<\/a> and Supplementary Table 1<\/a>). This region spans gradients of climate, topography, vegetation cover, permafrost extent and lithology, with elevations between 1,600 and 8,000\u2009m above sea level. The rivers drain the major tectonic terranes of the QTP that contain a wide range of sedimentary, igneous and metamorphic rocks51<\/a>,52<\/a> (Extended Data Fig. 1b<\/a>). In the northeast of the study area, the Yellow and Yangtze Rivers\u00a0drain marine carbonate and siliciclastic sediments of the eastern Kunlun\u2013Qaidam Terrane and the Triassic deep-marine sediments of the Songpan\u2013Ganzi flysch complex51<\/a>. The Lancang, Nu and the southwestern Yangtze catchments drain mostly terrestrial siliciclastic and shallow marine carbonate rocks of the eastern Qiangtang and northern Lhasa Terrane51<\/a>. In particular, the sediments of the Qiangtang Terrane contain abundant evaporites45<\/a>. By contrast, the catchments of the Yarlung Tsangpo and Indus Rivers are underlain by a series of plutonic and volcanic rocks of the southern Lhasa Terrane apart from the marine Tethyan sedimentary series of the Himalayas (Extended Data Fig. 1b<\/a>). Permafrost extent across the study area varies from continuous to isolated with the most extensive permafrost extent concentrated in the Yangtze and eastern Yellow Rivers, and high-elevation areas of all other river basins53<\/a>. The climate is characterized by a cold and dry winter season and an ice-free season between April and October with heavy monsoon rains during the study period54<\/a>. Potential glacial influence is minimized by sampling streams at least 20\u2009km downstream from glaciers. Flow distance and tributaries collectively reduce potential biogeochemical signatures of glacial meltwaters.<\/p>\n We collected 175 individual samples from 50 rivers across this area during the daytime in spring (May\u2013June), summer (July\u2013August) and autumn (September\u2013October) between 2016 and 2018, and in spring ice-out (early April) in 2023. The Yellow River basin was visited eight times; the Yangtze River\u00a0basin five times; and the Lancang, Nu and Yarlung River basins were all sampled on four dates, the Derung and Indus River basins twice, and the Marja River basin once. We measured, compiled, interpolated and\/or calculated all variables first for each of the 175 samples (Supplementary Table 2<\/a>) and then calculated mean values and their standard deviations for sites with multiple samples to obtain one estimate for each of the 50 sampling sites (Supplementary Table 5<\/a>).<\/p>\n Topography, hydrology, lithology and climate parameters<\/p>\n All topographic analyses were performed on a 90-m resolution digital elevation model (DEM) from the Shuttle Radar Topography Mission (SRTM)55<\/a> using version pre2.5 of the TopoToolbox56<\/a>. We extracted a river network and estimated, for each sample location, drainage-basin-averaged slopes, Sbsn, and normalized steepness index, ksn, with a reference concavity of 0.45, a segment length of 1,000\u2009m and a drainage-area threshold of 100\u2009km2. For all rivers in the network, we also extracted the stream order (with drainage-area threshold of 5\u2009km2) and calculated the total length of rivers, \\({L}_{{\\mathrm{riv}}_{i}}\\) with a given order, that is, upstream of each sampling point. All these analyses were performed on a DEM projected in ArcGIS-Pro using the inbuilt Asia North Lambert Conformal Conic projection. To provide an estimate of the total surface area of rivers upstream of a sampling point, Ariv, we multiplied the river length in each stream order, \\({L}_{{\\mathrm{riv}}_{i}}\\), with an average width for rivers in that stream order, \\({\\bar{W}}_{i}\\):<\/p>\n $${A}_{\\mathrm{riv}}=\\mathop{\\sum }\\limits_{1}^{N}{L}_{{\\mathrm{riv}}_{i}}{\\bar{W}}_{i}$$<\/p>\n \n (1)\n <\/p>\n Where N is the highest stream order in each catchment. For stream orders 3\u20138, \\({\\bar{W}}_{i}\\) and its uncertainty were based on the mean and standard deviation of 100 individual width measurements (only 30 measurements for stream order 8). Measurements from hydrologic stations at our sampling sites were combined with widths of randomly selected rivers in the study area mapped on Google Earth (Supplementary Table 6<\/a>). Average width for the first- and second-order streams were estimated by extrapolating an exponential fit to rivers of orders 3\u20136, with an uncertainty based on the prediction band at these stream orders (Supplementary Fig. 1<\/a>). We used only stream orders 3\u20136 for the extrapolation, because many higher-order rivers in the study area are inset into narrow bedrock gorges and break the exponential increase in stream-widths18<\/a> (Supplementary Fig. 1<\/a>). Our interpolation predicted the first-order stream widths of 4.3 (1.3\u201313.8)\u2009m, which covers the range reported in ref.\u200957<\/a> (1.9\u2009\u00b1\u20091.1\u2009m), but is wider than the value reported more recently in ref.\u200958<\/a> (0.32\u2009\u00b1\u20090.07\u2009m), leading to a potential overestimate in CO2 emissions.<\/p>\n The same DEM was analysed in TopoToolbox after projecting it with the Asia North Albers Equal Area Conic projection in ArcGIS-Pro to extract drainage areas, Absn, and to estimate an areal fraction of permafrost extent in the landscape upstream of each sampling point. We based the latter analysis on a published map of permafrost extent53<\/a> that classifies surfaces as being underlain by continuous (\u226590%), discontinuous (\u226550% and <90%), sporadic (\u226510% and <90%) or isolated (<10%) permafrost. Because this classification is based on ranges, we assigned values of 95%, 70%, 30% and 5%, respectively, to each permafrost class (0% for no permafrost), and then calculated the average pixel value upstream of each sampling point in TopoToolbox. Finally, we re-classified all catchment-averaged permafrost-extent values into the same permafrost categories. We generated a TIFF image from the global lithologic map59<\/a> to extract the fraction of carbonate lithologies upstream of each sampling point in TopoToolbox. In the fraction of carbonate rocks, we included the lithologic categories \u2018sc\u2014carbonate sedimentary rocks\u2019 and \u2018sm\u2014mixed sedimentary rocks\u2019. We did not consider \u2018mt\u2014metamorphics\u2019 here, but we found that including these would not make a substantial difference to our analysis. We estimated the fraction of wetlands, the mean annual air temperature and the mean annual precipitation values for the watersheds draining each of our sites by matching the coordinates of each site with the corresponding watershed in HydroATLAS60<\/a>, in which we obtained the relevant attributes.<\/p>\n For sampling sites within the cross-section of the hydrometric station, we obtained a measurement of discharge, Qspl (m3\u2009s\u22121), and suspended sediment concentrations (Supplementary Table 2<\/a>) at the time of sampling. We were also able to obtain average annual discharge, Qannual (m3\u2009yr\u22121), for the most downstream sampling sites of all sampled catchments (basin IDs 9, 16, 18, 19, 22, 27, 31, 37, 45, 46, 47, 48, 49 and 50) and annual suspended sediment for the four largest outlet basins (basin IDs 7, 16, 22 and 27; 74% of the total study area). Note that suspended sediment concentrations were not available for the sampling site with ID9, but for the more upstream sampling site ID7.<\/p>\n Sample collection and analysesCO2 concentrations and fluxes<\/p>\n Triplicate river-water samples for \\({p}_{{\\mathrm{CO}}_{2}}\\) were measured in situ from different locations at each site using a Qubit Dissolved CO2 System equipped with an internal pump and a \\({p}_{{\\mathrm{CO}}_{2}}\\) probe (S157-P, Qubit Biology). The \\({p}_{{\\mathrm{CO}}_{2}}\\) probe was coated with a polytetrafluoroethylene membrane sleeve that is permeable to CO2 gas, but not to water. The system was calibrated with standard CO2 gas ranging from 0\u2009\u03bcatm to 10,000\u2009\u03bcatm. Recalibration of the system was done before each round of fieldwork. Water\u2013air CO2 emission rates were measured with \\({p}_{{\\mathrm{CO}}_{2}}\\) collection simultaneously using four floating chambers deployed for 1\u2009h. At each transect, chambers were held in place above a range of water depths from shallow water within 2\u2009m of the river bank to deep water in the mid-channel. After mixing the contents in the chamber three times, 50\u2009ml of gas was retrieved from the chambers after 0\u2009min, 2\u2009min, 5\u2009min, 10\u2009min, 30\u2009min and 60\u2009min of measurement and transferred to air-tight gas sampling bags for concentration determination in the laboratory by gas chromatograph equipped with a flame ionization detector (Agilent 7890B GC-FID). Characteristics and construction of the chambers followed previously tested setups that showed minimally biased results61<\/a>. Local ambient air samples were also taken and used to back-calculate CO2 concentrations in water in equilibrium with the atmosphere.<\/p>\n Riverine CO2 emission rates, \\({F}_{{{\\rm{CO}}}_{2}}\\) (mmol\u2009m\u22122\u2009day\u22121) were derived as follows:<\/p>\n $${F}_{{{\\rm{CO}}}_{2}}=k\\times ({C}_{{\\rm{water}}}-{C}_{{\\rm{eq}}})$$<\/p>\n \n (2)\n <\/p>\n where k is the gas transfer velocity (m\u2009day\u22121), Cwater is the water gas concentration (mol\u2009m\u22123) and Ceq is the gas concentration (mol\u2009m\u22123) in water in equilibrium with the local atmosphere corrected for temperature-induced changes in solubility according to Henry\u2019s law. We assumed that \\({F}_{{{\\rm{CO}}}_{2}}\\) is solely diffusive (that is, CO2 ebullition is negligible). We obtained k from a linear regression of \\({p}_{{\\mathrm{CO}}_{2}}\\) against time to eliminate possible bias due to gas accumulation in the chamber headspace. Fluxes from all chambers were averaged to yield one measurement per sample.<\/p>\n Physicochemical analyses<\/p>\n Water samples were collected simultaneously with gas samples for physicochemical analyses at each site. Water temperature, pH and dissolved oxygen were measured in situ with portable field probes (Hach HQ40d). Air temperature, air pressure and wind speed were measured in situ with a portable anemometer (Testo 480).<\/p>\n Water samples were filtered through 0.22\u2009\u00b5m polyethersulfone syringe filters to measure cation and anion concentrations. We rinsed each polyethylene bottle with filtered water before filling, and acidified each bottle for cation analysis with ultrapure HNO3. We measured major cations (K+, Na+, Mg2+ and Ca2+) and anions (Cl\u2212 and SO42\u2212) on an ion chromatograph (Dionex ICS-1100). NH4+ and NO3\u2212 were determined colorimetrically with an Autoanalyzer-3 (Bran+Luebbe) following the methods in ref.\u200962<\/a>. We estimated HCO3\u2212 by charge balance. Dissolved silica (Si) was analysed on a Varian 720 ICP-OES. DOC concentrations were measured on duplicate filtered water samples acidified with H2SO4 following ref.\u200918<\/a>. DIC concentrations were calculated from the sum of HCO3\u2212 and \\({p}_{{\\mathrm{CO}}_{2}}\\) concentrations. To determine the POC content of the suspended sediments, water samples were filtered in duplicates through pre-weighed combusted 0.45\u2009\u00b5m filters. The filters were then dried at 105\u2009\u00b0C for at least 12\u2009h. Before POC analysis, the samples were acidified with 2\u2009M HCl to remove inorganic carbon. No rinsing step was applied after acidification to avoid washing fine particulate material or acid-labile organic matter off the filters. The filters were subsequently dried again at 105\u2009\u00b0C for over 5 days to ensure complete removal of residual acid and moisture. Afterwards, the filters were reweighed, and POC mass on the filters was quantified using a Vario Micro Cube elemental analyser (Elementar) with the furnace temperature set at 980\u2009\u00b0C.<\/p>\n Sulfur and oxygen isotope analyses<\/p>\n Stable isotopic composition of sulfur, oxygen and carbon is reported in the text using delta notations:<\/p>\n $${{\\rm{\\delta }}}^{34}{\\rm{S}},\\,{{\\rm{\\delta }}}^{18}{\\rm{O}}\\,\\mathrm{or}\\,{{\\rm{\\delta }}}^{13}{\\rm{C}}=\\frac{{R}_{\\mathrm{sample}}}{{R}_{\\mathrm{standard}}}-1$$<\/p>\n \n (3)\n <\/p>\n where Rsample and Rstandard are ratios of heavy to light isotopes (34S:32S, 18O:16O or 13C:12C) in samples (Rsample) and in standards (Rstandard) of standard mean ocean water or Vienna PeeDee Belemnite, respectively.<\/p>\n Stable isotopes of \u03b434S-SO42\u2212 and \u03b418O-SO42\u2212 of dissolved SO42\u2212 were measured following ref.\u200963<\/a> for 34 river samples collected during the fall 2018 sampling campaign. In brief, sulfate was eluted off a Dowex resin column with 20\u2009ml 0.8\u2009M HCl and mixed with BaCl2 to precipitate barite (BaSO4). The precipitate was subsequently cleaned with 6\u2009M HCl and rinsed three times with ultrapure water. After removing impurities with 10\u2009ml 0.05\u2009M DTPA, the barite was then cleaned three times with ultrapure water and dried at 70\u2009\u00b0C. For \u03b434S-SO42\u2212 analysis, 400\u2009\u03bcg of barite was combusted in excess oxygen with vanadium pentoxide in a Flash EA coupled by continuous flow and analysed by an isotope ratio mass spectrometry (Thermo Scientific 253 Plus). For \u03b418O-SO42\u2212 analysis, 180\u2009\u03bcg of barite was pyrolysed in a thermal conversion element analyser and passed through continuous He flow into an isotope ratio mass spectrometer (Thermo Scientific 253 Plus) using a ConFlo 3. We determined \u03b418O-H2O using a wave scanning cavity ring-down spectrometer (Picarro L2130i) as well. Additional triplet \u03b434S-SO42\u2212, \u03b418O-SO42\u2212 and \u03b418O-H2O data for rivers within the Yellow, Yangtze, Lancang, Nu, Yarlung Tsangpo and Indus catchments were obtained from the literature (Supplementary Table 7<\/a>). We had measurements from only 34 of the 50 sampling sites. For the remaining 16 sites, we substituted data from the same or closest sites reported in the literature. We then calculated the average for sites with multiple measurements to obtain a single estimate for each of the 50 sampling sites.<\/p>\n Carbon isotope analysis<\/p>\n We collected samples for \u03b413C\u2013\u039414C analyses of DOC, DIC and CO2 both near the river bank and in the mid-channel from the sampled rivers from 2017 to 2023 (except 2019 because of COVID restrictions). For \u039414C and \u03b413C measurements of DOC and DIC, the river-water samples were filtered using 0.45\u2009\u03bcm polysulfone filters, then the near-bank and mid-channel samples were mixed and collected in acid-washed (10% HCl v\/v, 24\u2009h) 1\u2009l HDPE Nalgene bottles that were rinsed three times with filtered river water before collection. Four submerged bubble traps were deployed from the near bank to the mid-channel to collect gas samples for \u039414C-CO2 and \u03b413C-CO2 analyses. The traps remained in place for 2\u20134\u2009h before sampling. Accumulated gas was transferred into 10-l air-tight Tedlar PVF gas sampling bags to prevent contamination by contemporary atmospheric CO2. \u03b413C was analysed on an isotope ratio mass spectrometer (Thermo Scientific 253 Plus), and \u039414C was measured at the accelerator mass spectrometry facility at Peking University, Guangzhou Institute of Geochemistry and Beta Analytic Laboratory.<\/p>\n We also compiled 92 concurrent \u03b413C\u2013\u039414C measurements each for DOC and DIC within the QTP from the literature. Finally, we were able to obtain 378 total observations of fluvial DOC (n\u2009=\u2009163), DIC (n\u2009=\u2009147) and CO2 (n\u2009=\u200968) \u039414C data, including our own measurements outlined above.<\/p>\n Endmembers of carbon isotopes<\/p>\n For Fig. 2a<\/a>, we used the following \u03b413C ranges for endmember carbon sources: terrestrial OC derived from globally dominant C3 plant photosynthesis (\u221232\u2030 to \u221224\u2030)64<\/a>; autochthonous, in-stream OC production (\u221223\u2030 to \u221217\u2030)13<\/a>; and rock sources, including kerogen (\u221232\u2030 to \u221217\u2030), carbonate rocks (\u22122.5\u2030 to 0\u2030) and solid Earth degassed CO2 (\u22126.5\u2030 to 10\u2030 dependent on source rocks)64<\/a>.<\/p>\n The Stable Isotope Mixing Models in R (simmr) package was used to estimate possible contributions of three endmember sources for all \u039414C samples as detailed described in ref.\u200946<\/a>. Endmember sources for the QTP were further constrained by using measured \u039414C values obtained from QTP soils (Supplementary Table 8<\/a>). A mass balance fitting the potential contributions of the endmember sources was solved by using a Markov Chain Monte Carlo method:<\/p>\n $${\\rm{C}}{14}_{\\mathrm{sample}}={\\rm{C}}{14}_{\\mathrm{topsoil}}+{\\rm{C}}{14}_{\\mathrm{deep}\\mathrm{soil}}+{\\rm{C}}{14}_{\\mathrm{rock}}$$<\/p>\n \n (4)\n <\/p>\n where C14sample represents the 14C content of the different waterborne carbon components of each sample (in pmC), which is the sum of the mass balance among the proportional contributions of three endmember sources. To solve this mass balance, the relative contributions of the three endmember sources were assumed to sum to 1. We divided the soil column into top (<50\u2009cm) and deep (\u226550\u2009cm) pools to where the topsoil source (C14topsoil) is a simplified representation of the uppermost soil carbon accumulation in the active layer, reflecting the seasonally thawed, often younger soils. The deep soil source (C14deep soil) represents older carbon from the deep active layer to the permafrost, including the Holocene and late Pleistocene aged sources from decaying over 5,000\u2009years bp carbon, and decaying more than 10,000\u2009years bp carbon, respectively. The rock source (C14rock) represents \u039414C-dead carbon sources from kerogen and carbonates.<\/p>\n Unmixing major elements<\/p>\n For all 175 samples, we used the inverse model mixing elements and dissolved isotopes in rivers (MEANDIR)65<\/a> to unmix contributions of solute sources to our measured riverine chemistry.<\/p>\n Considered solutes and endmembers<\/p>\n In the inversion, we considered measured concentrations of sodium [Na+], calcium [Ca2+], magnesium [Mg2+], chloride [Cl\u2212], sulfate [SO42\u2212] and silicate [Si], as well as the isotopic composition of sulfur in sulfate (\u03b434S-SO42\u2212). All elemental concentrations were normalized to the sum of the major ions [Na+], [Ca2+], [Mg2+] and [SO42\u2212] (refs.\u200945<\/a>,65<\/a>). We considered contributions from atmospheric sources and from the weathering of evaporites, sulfides, carbonates and silicates. Anthropogenic inputs were considered negligible because the area is far from major industrial activity or human settlements, and we did not consider the inputs from hot springs (Supplementary Text\u00a03<\/a>). The size of the study area and lithologic heterogeneity is large, so we used endmember compositions of silicate and carbonate weathering that reflect global averages23<\/a>,65<\/a>,66<\/a> (Supplementary Text\u00a04<\/a> and Supplementary Table 9<\/a>). The range of precipitation endmembers was derived from the 5th and 95th percentiles of average precipitation compositions in each of the drainage basins (Supplementary Table 10<\/a>). A range of different evaporite minerals contribute to the solute budget of the QTP rivers, which complicates the distinction between evaporite and non-evaporite solute sources45<\/a>. To address this challenge, we defined a very broad evaporite endmember that allowed for a wide range of compositions (Supplementary Table 9<\/a>). Moreover, to distinguish evaporite and silicate weathering sources, we included silica in the inversion and defined a range of allowed silica concentrations for the silicate endmember using the 5th and 95th percentiles of a compilation of unmixed silicate weathering concentrations from small mountain streams for which evaporite contributions are negligible35<\/a>. Finally, we used published \u03b434S compositions of sulfides and sulfur-bearing evaporites from our study area, to unmix sulfur sources (Supplementary Table 11<\/a>).<\/p>\n MEANDIR model parameters<\/p>\n MEANDIR uses a Monte Carlo approach to estimate the range of fractional endmember contributions that fit the measured data. At each iteration, the model picks a set of endmember compositions randomly from user-defined endmember ranges and then fits the fractional contributions of endmembers that best fit the observed data. The final result is commonly obtained as the median of a given fraction or number of model runs with the lowest misfit between the model and the observations. Whereas this approach produces a solution to all samples, it also allows very poor fits to the data. Hence, we iterated over each sample until we had 100 solutions that fit all of the ion data within \u00b115% and the \u03b434S-SO42\u2212 data within \u00b12\u2030. We allowed a maximum of 100,000 iterations to run before the inversion was stopped. Sixteen out of the 175 samples did not have all of the major element analyses data and 10 samples did not produce the necessary 100 solutions within the defined error bounds during the inversion. Hence, only 149 samples were considered in the endmember unmixing and carbon balance. All MEANDIR input parameters are reported in the\u00a0Supplementary Code<\/a>.<\/p>\n Treatment of potassium<\/p>\n Inclusion of potassium [K+] in the inversion doubles the number of samples that cannot be well-modelled. Because [K+] and [Na+] are well-correlated (Supplementary Table 2<\/a>), we made the simplifying assumption that fractional contributions to [K+] scale linearly with the fractional contributions to [Na+]. Because the concentrations of potassium are generally low, the uncertainty introduced by this approximation is small.<\/p>\n Areal CO2 fluxes from weathering<\/p>\n From the cation charge contributed by carbonate and silicate weathering and the anion charges contributed by sulfate from sulfide oxidation, we estimated the balance of CO2 sinks and sources from weathering in the landscapes drained by each sampled river. Any carbonate that was already precipitated in streams or soils upstream of the sampling point did not contribute to our CO2 budget and is not considered here. We followed previous approaches28<\/a>,67<\/a> to estimate the following fluxes of CO2, \\({F}_{{{\\rm{CO}}}_{2}}\\), in mass of carbon per area per time (tC\u2009km\u22122\u2009yr\u22121) with positive numbers being CO2 production and negative numbers being CO2 sequestration:<\/p>\n 1.<\/p>\n CO2 consumed by silicate weathering, \\({F}_{{(\\text{C}{\\text{O}}_{2})}_{\\text{sil}}}\\):<\/p>\n $${F}_{{({\\mathrm{CO}}_{2})}_{\\mathrm{sil}}}=-(1-{f}_{\\mathrm{sulf}})\\times {M}_{{\\rm{C}}}\\times {[\\mathrm{Cat}]}_{\\mathrm{sil}}^{\\mathrm{eq}}\\times \\frac{{Q}_{{\\rm{w}}}}{{A}_{\\mathrm{bsn}}}$$<\/p>\n \n (5)\n <\/p>\n 2.<\/p>\n CO2 consumed by carbonate weathering, \\({F}_{{(\\text{C}{\\text{O}}_{2})}_{\\text{carb}}}\\):<\/p>\n $${F}_{{({\\mathrm{CO}}_{2})}_{\\mathrm{carb}}}=-0.5(1-{f}_{\\mathrm{sulf}})\\times {M}_{{\\rm{C}}}\\times {[\\mathrm{Cat}]}_{\\mathrm{carb}}^{\\mathrm{eq}}\\times \\frac{{Q}_{{\\rm{w}}}}{{A}_{\\mathrm{bsn}}}$$<\/p>\n \n (6)\n <\/p>\n 3.<\/p>\n CO2 produced by sulfuric acid weathering of carbonate, \\({F}_{{(\\text{C}{\\text{O}}_{2})}_{\\text{sulf}}}\\):<\/p>\n $${F}_{{({\\mathrm{CO}}_{2})}_{\\mathrm{sulf}}}=0.5{f}_{\\mathrm{sulf}}\\times {M}_{{\\rm{C}}}\\times {[\\mathrm{Cat}]}_{\\mathrm{carb}}^{\\mathrm{eq}}\\times \\frac{{Q}_{{\\rm{w}}}}{{A}_{\\mathrm{bsn}}}$$<\/p>\n \n (7)\n <\/p>\n where fsulf () is the fraction of weathering by sulfuric acid, \\({f}_{\\text{sulf}}=\\frac{{[{\\mathrm{SO}}_{4}^{2-}]}_{\\text{sulf}}^{\\text{eq}}}{{[\\mathrm{Cat}]}_{\\text{sil}}^{\\text{eq}}+{[\\mathrm{Cat}]}_{\\text{carb}}^{\\text{eq}}}\\), \\({[\\mathrm{Cat}]}_{\\text{sil}}^{\\text{eq}}\\) and \\({[\\mathrm{Cat}]}_{\\text{carb}}^{\\text{eq}}\\) (charge equivalents in mol\u2009m\u22123) are the cation charges contributed by silicate and carbonate weathering, respectively, Mc (g\u2009mol\u22121) is the molar mass of carbon, Qw (m3\u2009yr\u22121) is the annual average discharge, and Absn (m2) is the drainage basin area. Equations (5<\/a>)\u2013(7<\/a>) assume that the sulfuric acid partitions between silicate and carbonate phases in the proportion of the fraction of cation charges contributed by each reaction. Depending on the conditions in the soil, this simplification may not apply68<\/a>. However, changing this assumption would not affect the net CO2 balance from weathering reactions, that is, the sum of equations (5<\/a>)\u2013(7<\/a>), but would instead trade off CO2 drawdown fluxes between silicates and carbonates. We note that in the calculation of the CO2 balance, we excluded one (out of 175) sample with fsulf\u2009>\u20091, which is impossible. This estimate must have arisen because of the high evaporite load that was not correctly estimated by the inversion model. We also excluded one sample that had an anomalously high sulfate concentration compared with all other 175 samples (almost 10 times higher than the 75th percentile of all samples and twice as high as the next highest concentration), and disproportionately dominated the upscaling described in the next sections. Finally, we note that we did not have annual discharge estimates for all sampling sites (see section \u2018Topography, hydrology, lithology and climate parameters\u2019), and instead upscaled Qw (m3\u2009yr\u22121) from the discharge estimate at the time of sampling, Qspl (m3\u2009s\u22121), for each of the 175 samples by simply multiplying the value by the number of seconds in a year. Therefore, Qspl for the 50 sampling sites is effectively a weighted average flux. Note that these fluxes were used only to compare relative changes across the permafrost gradient. To obtain the total fluxes from the study area, we used actual measured annual averages for the most downstream sites of the catchments (see sections \u2018Topography, hydrology, lithology and climate parameters\u2019 and \u2018Total carbon fluxes from the study area\u2019).<\/p>\n Areal CO2 fluxes from river emissions<\/p>\n We estimated areal CO2 emission fluxes from each of the sampled catchments (in mass of carbon per area per time) based on CO2 emission flux rates (\\({F}_{{\\mathrm{CO}}_{2}}\\)) measured with the floating chambers as<\/p>\n $${F}_{{({{\\rm{CO}}}_{2})}_{{\\rm{rem}}}}=365.25\\times {F}_{{{\\rm{CO}}}_{2}}\\times {M}_{{{\\rm{CO}}}_{2}}\\times \\frac{{A}_{{\\rm{riv}}}}{{A}_{{\\rm{bsn}}}}$$<\/p>\n \n (8)\n <\/p>\n where \\({F}_{{({{\\rm{CO}}}_{2})}_{{\\rm{rem}}}}\\) (mol\u2009m\u22122\u2009day\u22121) is the upscaled riverine emission (rem) flux.<\/p>\n Total carbon fluxes from the study area<\/p>\n We estimated the total CO2 fluxes from weathering and river emissions as well as downstream fluxes of DOC and DIC across the study area in TgC\u2009yr\u22121. To this end, we summed the contributions of CO2 fluxes from the 14 most downstream sampling sites of each sampled basin. We calculated the summed weathering, DOC, and DIC fluxes based on the annual discharge estimates, Qannual multiplied by the average concentrations in the basin, and we used the standard deviations of the concentrations to estimate the uncertainty of the summed value (Supplementary Table 5<\/a>). Here, we were able to use actual measured annual discharge measurements and did not use the Qw values upscaled from individual discharge measurements at the time of sampling. Note that we used the solute and carbon concentrations of basin 7 instead of the further downstream basin 9 because we did not have any major element chemistry for basin 9. To account for uncertainty in the stream-area calculations, we estimated the summed CO2 emission fluxes, \\({F}_{{(\\text{C}{\\text{O}}_{2})}_{\\text{rem}}}\\), with 100,000 Monte Carlo runs. First, we randomly picked values of river widths from a normal distribution for stream orders 3\u20138 (Supplementary Table 5<\/a>) and from a log-normal distribution for the interpolated stream orders 1\u20132. Multiplied by stream lengths (equation (1<\/a>)), these values yield 100,000 estimates of stream areas upstream of each sampling site. We then randomly picked 100,000 \\({F}_{{({{\\rm{CO}}}_{2})}_{{\\rm{rem}}}}\\) values from a normal distribution based on the average and standard deviation for each sample site (Supplementary Table 5<\/a>) and solved equation (8<\/a>) for each of the 100,000 Monte Carlo runs. The final reported value is the median and standard deviation from all 100,000 runs. We also multiplied annual suspended sediment fluxes by POC content of the largest four of the outlet basins to obtain POC fluxes. Because our measurements were limited to the ice-free season between April and October, and because our first- and second-order stream widths are higher than expected, we probably overestimate CO2 emissions over the year.<\/p>\n Statistical analysis<\/p>\n Using R (v.4.4.2), we performed a PCA to explore the main environmental drivers and their effects on our variables of interest. The drivers considered were catchment size, average ksn, permafrost cover, average annual precipitation, average annual air temperature, wetland cover and fraction of carbonate rocks within the catchment. These variables were scaled and centred before the PCA. We then used the two main axes of the PCA (PC1 and PC2) to explore how they were related to several variables using linear regression models. Apart from the PCA, we used an MLR model, in which each response variable was regressed against all key environmental variables (Extended Data Table 1<\/a> and Supplementary Table 4<\/a>). All predictor and response variables were standardized before analysis to facilitate the comparison of effect sizes using standardized regression coefficients. We used a bidirectional stepwise variable selection procedure based on the Akaike information criterion to identify the most parsimonious set of drivers for each response variable. Model performance and the significance of individual predictors were assessed using adjusted R2 and P\u00a0values.<\/p>\n","protected":false},"excerpt":{"rendered":"Study area and field sampling Spanning about 30% of the entire QTP, this 780,000\u2009km2 study area is situated…\n","protected":false},"author":2,"featured_media":505462,"comment_status":"","ping_status":"","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[7],"tags":[9850,2026,13989,61,60,2027,82],"class_list":["post-505461","post","type-post","status-publish","format-standard","has-post-thumbnail","category-science","tag-carbon-cycle","tag-humanities-and-social-sciences","tag-hydrology","tag-ie","tag-ireland","tag-multidisciplinary","tag-science"],"_links":{"self":[{"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/posts\/505461","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/users\/2"}],"replies":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/comments?post=505461"}],"version-history":[{"count":0,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/posts\/505461\/revisions"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/media\/505462"}],"wp:attachment":[{"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/media?parent=505461"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/categories?post=505461"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/www.newsbeep.com\/ie\/wp-json\/wp\/v2\/tags?post=505461"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}