Spatially Resolved Measurements in Tropical Reservoirs Reveal Elevated Methane Ebullition at River Inflows and at High Productivity

An increasing number of rivers is being dammed, particularly in the tropics, and reservoir water surfaces can be a substantial anthropogenic source of greenhouse gases. On average, 80% of the CO2‐equivalent emission of reservoirs globally has been attributed to CH4, which is predominantly emitted via ebullition. Since ebullition is highly variable across space and time, both measuring and upscaling to an entire reservoir is challenging, and estimates of reservoir CH4 emission are therefore not well constrained. We measured CH4 ebullition at high spatial resolution with an echosounder and bubble traps in two reservoirs of different use (water storage and hydropower), size and productivity in the tropical Brazilian Atlantic Rainforest biome. Based on the spatially most well‐resolved whole‐reservoir ebullition measurements in the tropics so far, we found that mean CH4 ebullition was twice as high in river inflow areas as in other parts of the reservoirs, and more than 4 times higher in the eutrophic compared to the oligotrophic reservoir. Using different upscaling approaches rendered similar whole‐reservoir CH4 ebullition estimates, suggesting that highly spatially‐resolved measurements may be more important for constraining reservoir‐wide CH4 estimates than choice of upscaling approach. The minimum sampling effort was high (>250 and >1,700 thirty‐meter segments of hydroacoustic survey to reach within 50% or 80% accuracy, respectively). This suggests that traditional manual bubble‐trap measurements should be abandoned in favor of highly resolved measurements in order to get spatially representative estimates of CH4 ebullition, which accounted for 60% and 99% of total C emission in the two studied reservoirs.

reservoirs can be substantial. Reservoir water surfaces constitute about 1.3% of the global anthropogenic emission of carbon (Deemer et al., 2016). Some hydropower reservoirs have even been found to be comparable in terms of gross greenhouse gas emission (i.e., not excluding emissions that occurred even prior to dam construction) per produced energy unit to coal or natural gas combined cycle plants (Abril et al., 2005;de Faria et al., 2015;Demarty & Bastien, 2011;Hertwich, 2013), while other hydropower reservoirs seem to emit much less (Scherer & Pfister, 2016;Teodoru et al., 2012).
Most of the greenhouse gas emission from reservoir surfaces occurs as methane (CH 4 ), carbon dioxide (CO 2 ), and nitrous oxide (N 2 O). The most recent compilation of greenhouse gas emission from reservoirs worldwide estimated that reservoirs account for an emission of 0.8 Pg CO 2 equivalents per year, and suggested that on average 79% of the reservoir water surface emission occurs as CH 4 , 17% as CO 2 , and 4% as N 2 O (Deemer et al., 2016). These calculations assume a 100-year timespan to express all reservoir emission as CO 2 -equivalent using the radiative forcing of CH 4 and N 2 O, which is 34 and 282 times higher than that of CO 2 , respectively (Myhre et al., 2013). Whereas CO 2 is highly soluble in water and emitted primarily via diffusion, the low solubility of CH 4 means that it is frequently emitted via ebullition (bubbles). On average, 65% of reservoir CH 4 emission occurs via ebullition and only 35% via diffusion (Deemer et al., 2016). It is therefore important to understand the dynamics of CH 4 ebullition as a pathway of greenhouse gas emission from reservoir surfaces.
Unfortunately, CH 4 ebullition is difficult to accurately quantify because it is highly variable across space and time, even within a single system Linkhorst et al., 2020;Maeck et al., 2013). Obtaining measurements representative for an entire water body is not trivial, especially when even massive sampling efforts usually only cover a fraction of the reservoir during a limited time period. Ultimately, these constraints severely inhibit upscaling CH 4 ebullition to a whole reservoir and beyond the time period of measurement. Current techniques for measurement of CH 4 ebullition can cover space but lack coverage in time (hydroacoustics: DelSontro et al., 2010DelSontro et al., , 2011Tušer et al., 2017), or cover time but not space (bubble traps: Descloux et al., 2017;dos Santos et al., 2017). While eddy covariance can produce accurate point measurements of methane emission at a very high temporal resolution and integrated over space (Deshmukh et al., 2014;Eugster et al., 2011;Podgrajsek et al., 2016), the footprint from a single eddy covariance tower in a fixed location is relatively small, which cannot adequately represent large or dendritic systems. Current procedures for upscaling in space either apply an average of point measurements to the entire system (dos Santos et al., 2017), or apply randomization techniques to the sampling strategy to reduce biases during upscaling , or use a categorization of sub-areas (e.g., inflow and non-inflow areas; DelSontro et al., 2010DelSontro et al., , 2011 or bathymetry (Descloux et al., 2017). These procedures have been applied to reservoirs that are tens to thousands of km 2 in size. In how far these different upscaling approaches affect the system-wide estimate of CH 4 ebullition, however, is at present unknown.
In this study, we measured CH 4 ebullition at high spatial resolution using hydroacoustics in two reservoirs in the Brazilian Atlantic Rainforest biome. The reservoirs differ in size (9.5 and 35 km 2 ), purpose (water storage and generation of electricity) and productivity (oligotrophic and eutrophic). We investigated how various upscaling approaches from highly spatially resolved measurements affect whole-reservoir ebullition estimates and also calculated the minimal sampling effort. We used bubble trap deployments during dry and wet season for an estimation of seasonal differences in CH 4 ebullition. Finally, we combined our whole-reservoir CH 4 ebullition estimates with whole-reservoir CH 4 and CO 2 diffusion estimates (Paranaíba et al., 2018(Paranaíba et al., , 2021 to calculate the total C emission for these two contrasting tropical reservoirs.

Sampling Sites and Strategy
We measured CH 4 ebullition in two contrasting reservoirs in Brazil: the oligotrophic flood regulation and water supply reservoir Chapéu d'Uvas (CDU) in the state of Minas Gerais, and the eutrophic storage hydropower reservoir Funil (FUN) in the state of Rio de Janeiro. Both reservoirs are located in the Atlantic Rainforest biome and are exposed to a climate at the threshold of subtropical and tropical, with hot wet summers and colder drier winters. CDU has a mean surface area of 9.5 km 2 , and had a maximum water column depth of 34 m during our sampling campaign ( Figure S1). The total regulation amplitude in CDU between 2008 and 2019 was 12 m, and the mean seasonal regulation amplitude (i.e., change over 1 year with two major hydrological peak seasons) was 7 m. The catchment area of CDU is 310 km 2 , and land cover is primarily grassland and Atlantic Rainforest ( Figure S2; Table S1). FUN has a mean surface area of 35 km 2 and had a maximum depth of 60 m during our sampling campaign ( Figure S1). The total regulation amplitude in FUN between 2008 and 2019 was 21 m, and the mean seasonal regulation amplitude was 14 m. The catchment area of FUN is 13,518 km 2 , primarily consisting of cultivated terrestrial areas and Atlantic Rainforest ( Figure S3; Table S1). Both reservoirs can experience low oxygen conditions in their bottom waters (Figure S4). During our 2-week sampling campaigns, we measured CH 4 ebullition hydroacoustically with an echosounder and with bubble traps in May 2016 in CDU (high stagnant water level) and in October 2016 in FUN (falling water level), and with only bubble traps in March 2015 in FUN (wet season, rising water level) and in April 2015 in CDU (dry season, low stagnant water level). Both CDU and FUN are dendritic in shape, and several smaller rivers flow into the reservoirs in addition to their respective main inflow (Figures S2 and S3). The overall sampling strategy was to survey different parts of the reservoirs in order to gain data for calculation of a system-wide CH 4 ebullition estimate. Thus, we surveyed large parts of the main stem of each reservoir as well as various types of side bays. Because sediment deposition in river inflow areas has been linked to high CH 4 ebullition Maeck et al., 2013;Quadra et al., 2020;Sobek et al., 2012), we also surveyed similar numbers of bays with and without a river inflow.

CH 4 Ebullition Measurements Using Hydroacoustics and Bubble Traps
We used a portable scientific echosounder (Simrad EY60, Kongsberg Maritime AC) with a downward-facing 7° split-beam transducer operating at 120 kHz. The transducer was mounted to a boat, which we navigated in zigzag lines across the reservoirs (Figure 1). We divided the echosounder track into segments of ∼30 m hydroacoustic track, and used the mean ebullitive CH 4 flux per segment for further calculations. Depending on water quality, the transducer was operated at a ping rate of 5-10 Hz. We used only the 3 m closest to the sediment for data analysis, as the highest quality data is found at depth. In shallow regions, data within 1.5 m of the transducer face were excluded due to near-field disturbance. We calculated bubble density per water volume that was covered by the echosounder beam, and then converted bubble density into bubble gas flux, following DelSontro et al. (2011DelSontro et al. ( , 2015. We analyzed the resultant hydroacoustic data with Sonar5 Pro (Lindem Acquisition), which allows for distinction between bubbles and fish or plankton based on their travel path, velocity and signal intensity (DelSontro et al., 2015;Ostrovsky et al., 2008). We divided each hydroacoustic file into segments of 300 ± 5 pings, which equates to distances of ∼30 m, and used the mean flux per segment for further calculations. We accounted for bubble dissolution during rise using the model LINKHORST ET AL.  described in McGinnis et al. (2006), assuming an average CH 4 content of 70% in the bubbles (according to bubble trap measurements as described below). Hydroacoustic surveys were completed in various regions over 6 days between 7 and 20 May 2016 in CDU, and over 8 days between 11 and 21 October 2016 in FUN. The hydroacoustic data were later used for all upscaling approaches and to calculate whole-system CH 4 ebullition.
As bubble traps we used inverted funnels (custom-built conical nylon frames covered with impermeable tarp) that were 50 cm in diameter and 43 cm tall, and collected the gas in 220-mL glass bottles that were screwed into threads at their narrow end at the top. Traps, including the 220-mL bottles, were completely submerged in water during deployment. Two 500-mL polyethylene terephthalate bottles were attached to each trap to keep them floating just underneath the water surface. Each trap was attached to a buoy ∼1.5 m away that was anchored to the sediment; thus, the trap could move freely around the buoy and was not impacted by anchor weight. Each sampling day, we deployed 5-10 bubble traps at different sites (e.g., one inflow and one non-inflow bay) in one area of the reservoir for 3-5 h. At the end of each deployment, we sampled the gas from the traps' glass bottles using 60-mL syringes. We measured the gas volume from the volume of displaced water in the glass bottles, and estimated the CH 4 concentration of the captured gas on an ultraportable greenhouse gas analyzer (UGGA, Los Gatos Research; see Paranaíba et al. (2018) for details on manual gas injection of discrete samples into the UGGA). In CDU, a total number of 75 bubble traps were deployed in various regions of the reservoir on 8 days between 7 April and 4 May 2015, and 33 bubble traps on 5 days between 7 and 18 May 2016. In FUN, 48 bubble traps were deployed on 9 days between 15 and 25 March in 2015, and 89 bubble traps on 7 days between 1 and 21 October 2016 (see Figures S5 and S6 for their distribution). The bubble trap measurements were solely used for quantification of the variability in CH 4 ebullition between two sampling campaigns during 2 consecutive years, and for calculation of mean CH 4 concentration in the captured bubbles for the hydroacoustic approach. Because of their limited spatial coverage, bubble trap data have not been integrated in the upscaling approaches or in the calculation of whole-system CH 4 ebullition.

Bathymetry
We combined all bathymetric data from reservoir-wide surveys during 2014-2016 that were acquired with two different measurement devices: the Simrad EY60 echosounder as used in this study for determination of ebullitive CH 4 flux, and an Innomar SES-2000 compact sub-bottom profiler, both with a resolution of ∼0.1 m. All measurements were interpolated with inverse distance weighting (IDW, ArcGIS version 10.6.1). We adjusted the depth measurements according to the water level during the respective surveys, such that the bathymetric maps represent May 2016 for CDU and October 2016 for FUN.

Downsampling
To calculate the minimum sampling effort needed for high-accuracy CH 4 ebullition measurements (according to Wik et al., 2015), we simulated hypothetical sampling scenarios. We calculated how many sampling points (i.e., 30-m hydroacoustic segments) were needed to reach within ±20% of the mean CH 4 ebullition calculated from all the 30-m segments of the hydroacoustic surveys (CDU: n = 2,696; FUN: n = 2,982) with high confidence (95% probability). We randomly chose n number of ebullition values (CDU: 10-2,600; FUN: 10-2,900) to compare against the true mean CH 4 ebullition for each reservoir 100 times.

Upscaling to Whole-Reservoir CH 4 Ebullition
We upscaled the hydroacoustically measured CH 4 ebullition in space using four different approaches (I-IV) to calculate whole-reservoir CH 4 ebullition during the sampling period.
For Approach I (IDW interpolation), we interpolated all fluxes as retrieved from the hydroacoustic measurements (in segments of ∼30 m, see Figure 1) using the IDW algorithm across the whole area of the reservoir. We then applied a grid with cells of 0.000700 × 0.000748° (decimal degrees, dd), which is equivalent to 77.5 × 77.5 m (6006.25 m 2 ) at 21° latitude (see Text S1, supporting information, for details on optimal cell size determination). For each reservoir, we then extracted the mean ebullitive CH 4 flux from the pixels of the interpolation within each grid cell. We multiplied each grid cell mean with its individual grid cell area and summed up all area-resolved means for each reservoir to retrieve the total ebullitive CH 4 emission from each reservoir water surface.
For Approach II (binned by bathymetry), we related, for each reservoir separately, the hydroacoustic CH 4 ebullition fluxes (in ∼30 m segments) to the mean water column depth within each segment. Such a relationship between ebullition and water column depth is frequently mentioned in the literature on low-latitude reservoir ebullition (Abril et al., 2005;Deshmukh et al., 2014;Keller & Stallard, 1994). We made depth bins of 1.5-3 m, 3-5 m, and from then on in 5 m bins, and determined the mean ebullitive CH 4 flux for each of these depth bins. By applying the above-mentioned 77.5 × 77.5 m grid on the bathymetric maps, we extracted the mean water depth within each grid cell. We then applied the mean ebullitive CH 4 flux per depth bin to the respective grid cells of the same depth range across the whole reservoir surface area.
For Approach III (grouped into inflow and non-inflow areas and binned by bathymetry), river inflow areas were identified and delimited based on the rationale that a river inflow area is a part of the reservoir within which the unidirectional water flow of the incoming river successively slows down, until it ceases and is replaced by the lake-like hydrodynamics of the standing water in the reservoir; the slowing-down of water flow allows the suspended sediment load of the inflowing river to settle, which can promote the build-up of methanogenic sediments (e.g., Maeck et al., 2013;Sobek et al., 2012). We assumed that the outer limit of a river inflow area (where unidirectional water flow ceases) is proportional to the ratio between the volume of water entering the bay with river inflow and the volume of the bay that receives the water. We then used visual field observations in one well-studied bay of CDU, and a bathymetric feature (shallow sill) in FUN to determine this ratio, and applied it to all other river inflow bays of each respective reservoir in order to determine the outer limit of each river inflow bay. Thereby, we assumed that local and temporally variable hydrographic and hydrodynamic conditions have a minor impact on the delimitation of the river inflow area compared to the volumes of inflowing water and the volume of the receiving bay. In the absence of robust criteria for river inflow delimitation from the literature and given the unfeasibility of within this study establishing validated hydrodynamic models for these complex dendritic reservoirs, we suggest that the combination of basic physical reasoning and field observations make our approach a defendable approximation. In detail, we determined the incoming rivers from Google Earth (version 7.3.2.5776, image from March 5, 2019), and calculated the catchment area of each of these rivers in ArcGIS using an Aster Global Digital Elevation Model (DEM) layer as provided by USGS Earth Explorer (May 23, 2019; data acquisition date: October 17, 2011, last update November 3, 2016). From a well-studied inflow-bay in CDU and from the main inflow in FUN, we calculated the water volumes of the river-inflow areas (i.e., the part of the bay which we assumed to be under the influence of the incoming river of this bay) using the surface volume tool in ArcGIS. In CDU, we delimited the river-inflow area from the non-inflow (i.e., lake-like) reservoir area based on field observations of water flow during 15 field campaigns over one year in a northern bay with river inflow (Linkhorst et al., 2020) ( Figure S7: river inflow 6); downstream of the placed delimitation, we could not observe a visible unidirectional flow of surface water, and assumed that this corresponds to the end of the river inflow area. In FUN, we delimited the river inflow area behind a large lateral land dam, that is, a levee, at a shallow sill in the main stem in the northernmost large bend of the reservoir ( Figure S8: river inflow 1). It is likely that the particle load of the large main river inflow in FUN is transported this far but does not pass the shallow sill ( Figure S1).
We then calculated the ratio V R /A R between the catchment area A R (m 2 ) of the river coming into the respective reference bays in CDU and FUN, and the water volume V R (m 3 ) that we assumed to be under influence of these respective rivers (Tables S2 and S3). For each of the other inflow-bays, we then multiplied this V R / A R ratio with the respective catchment area A X of each inflowing river x to calculate the volume of water V X that is influenced by the inflowing river: For each of the inflow-bays for which we had hydroacoustic measurements, we extracted all ebullitive CH 4 fluxes per segment of ∼30 m that fell into the surface water area of the estimated V X . We grouped the mean ebullitive CH 4 flux values from all inflow bays that we covered during our hydroacoustic survey into depth categories of 1.5-3 m, 3-5 m, and from then on 5-m bins and calculated the mean ebullitive CH 4 flux for each of the depth bins. We then applied the mean flux per depth from measurements in the reference inflow bays to all inflow bays, by applying the mean of each depth bin to the respective depth of the gridded bathymetric map ( Figure S1). Similarly, we calculated the mean ebullitive CH 4 flux per depth category for the rest of each reservoir, which we considered as non-inflow areas, and applied them to the respective depth of the gridded bathymetric map of all non-inflow areas. We then calculated the whole-reservoir mean flux considering both inflow and non-inflow areas together.
For Approach IV (grouped into hydromorphological sections and binned by bathymetry), we divided each reservoir into five different hydromorphological sections, according to the major morphometric features visible on a map: dam area (A), secondary tributaries (B and C), main stem (D), and main tributary (E) ( Figures S7 and S8). The delimitations of these sections were purposely arbitrary (i.e., did not follow the distinction of inflow and non-inflow areas as used in Approach III) to simulate the potential effect of using a rough but fast approach to divide the reservoir into sub-sections. We then established, for each of the sections separately, a relationship of CH 4 ebullition over binned depth for the measured hydroacoustic data. We calculated the mean CH 4 ebullition per depth bin and applied the means per depth bin on the gridded bathymetric map of each section separately. Finally, we calculated the whole-reservoir mean flux considering all sections together.

Calculation of Whole-Reservoir C Emission
We multiplied the mean CH 4 ebullition and the water surface area within each individual grid cell and summed them up to estimate total reservoir CH 4 ebullition for each upscaling approach. Total C emission from each reservoir was calculated as CH 4 ebullition plus CH 4 and CO 2 diffusion, which was taken from Paranaíba et al. (2021). CO 2 ebullition is negligible in these reservoirs (unpublished data). We applied a global warming potential of 34 to CH 4 , according to Myhre et al. (2013), to convert CH 4 to CO 2 equivalents.
Calculations were done in R version 3.5.2 and Matlab version R2018a, and Wilcoxon tests (one-sided, unpaired wilcox.test) in R version 3.5.1. All geographical analyses and the maps were made in ArcGIS version 10.6.1, and all other graphics were created in R version 3.5.2 and assembled in Inkscape 0.92.

Hydroacoustics
Mean CH 4 ebullition, as measured from the hydroacoustic surveys per ∼30 m segment (Figure 1), was 87.2 mg C m −2 d −1 in CDU in May 2016 and 210 mg C m −2 d −1 in FUN in October 2016. Even though the mean flux between reservoirs differed by a factor of 3, the range of CH 4 ebullition rates from individual segments was similar, ranging over 6 orders of magnitude (CDU: 0-22,149 mg C m −2 d −1 ; FUN: 0-21,747 mg C m −2 d −1 ). In both reservoirs, mean CH 4 ebullition was higher in areas under the influence of a river inflow than in the rest of the reservoir. Mean CH 4 ebullition was significantly higher in inflow areas (152.7 and 475.6 mg C m −2 d −1 ) than in non-inflow areas (63.1 and 217.5 mg C m −2 d −1 ) in CDU and FUN, respectively ( Figure 2; Wilcoxon for CDU: p = 3.6e-08, W = 337,037 and for FUN: p < 2.2e-16, W = 1,090,955). In CDU, we found the highest CH 4 ebullition in the shallowest (1.5-3 m) regions and a decreasing trend toward the deepest regions of the reservoir-both in inflow and non-inflow areas ( Figures S9 and S10). In FUN, we found no consistent pattern of CH 4 ebullition over depth ( Figure S9). Instead, we found contrasting depth patterns for inflow and non-inflow areas. In inflow areas of FUN, we found highest CH 4 ebullition in the deepest (25-31 m) regions, with a trend of decreasing CH 4 ebullition toward the shallowest (1.5-3 m) regions ( Figure S10). In non-inflow areas of FUN, there was no evident trend over depth, but the shallowest depth bin had the highest CH 4 ebullition ( Figure S10). The V R /A R ratio that we estimated for the model bays was similar for the two reservoirs, 0.0110 for CDU and 0.0167 for FUN (see Tables S2 and S3).

Upscaling Approaches
Despite the clearly different spatial patterns of ebullition that the upscaling approaches produced (Figures 3  and 4), there was little variability between the upscaled mean areal and mean whole-reservoir ebullitive fluxes that these different approaches yielded (mean areal flux in CDU: 48.9-49.2 mg C m −2 d −1 , and in FUN: 189-214 mg C m −2 d −1 ; whole-reservoir flux in CDU: 436-463 kg C d −1 , and in FUN: 6,531-7,569 kg C d −1 ; Table 1a). Variability of mean fluxes in this section are given as the standard error of the mean (SEM).
Approach I (IDW interpolation) returned a mean ebullitive CH 4 emission interpolation across the whole reservoir of 49.2 ± 2.4 mg C m −2 d −1 (range of grid cells: 0-29,475 mg C m −2 d −1 ) for CDU, and 202 ± 3 mg C m −2 d −1 (range of grid cells: 0-27,821 mg C m −2 d −1 ) for FUN (Figures 3 and 4; Table 1a). In CDU, we found maximum CH 4 ebullition in proximity of incoming rivers, and in the main stem of the middle of the reservoir (Figure 3). CH 4 ebullition from the main river inflow (122 mg C m −2 d −1 ), however, was similar to ebullition from all inflow bays (130 mg C m −2 d −1 ). In FUN, we found maximum CH 4 ebullition in proximity of incoming rivers, especially along the main stem near the main river inflow (Figure 4). CH 4 ebullition was lowest in most regions close to the dam and along the main reservoir stem.
Approach II (binned by bathymetry; depth groups: Figure FUN). In FUN, the main stem in proximity of the main river inflow is comparatively shallow (<10 m; Figure S1) and consequently showed high CH 4 ebullition (range: 293-405 mg C m −2 d −1 ; Figure 4).
Approach III (grouped into inflow and non-inflow areas and binned by bathymetry; depth groups: Figure S10) resulted in a spatial ebullition pattern similar to Approach II in CDU (Figure 3), but not in FUN where there are many shallow bays without river inflows (Figure 4). Approach III consequently resulted in a similar mean areal flux for CDU as that from Approach II (49.1 ± 0.6 mg C m −2 d −1 , Table 1a; range of grid cells: 3-150 mg C m −2 d −1 , Figure 3) but a slightly higher mean areal flux for FUN (214 ± 2 mg C m −2 d −1 , Table 1a; range of grid cells: 15-1,950 mg C m −2 d −1 , Figure 4).
Approach IV (grouped into hydromorphological sections and then binned by bathymetry; depth groups: Figures S11 and S12) resulted in a decreasing trend in CH 4 ebullition with increasing depth in all sections of CDU ( Figure S11), while depth patterns differed largely between sections in FUN ( Figure S12). For both reservoirs, Approach IV produced a similar spatial pattern of CH 4 ebullition as that of Approach III (Figures 3 and 4). We found a mean CH 4 ebullition of 48.9 ± 0.6 mg C m −2 d −1 (Table 1a;  1-485 mg C m −2 d −1 , Figure 3) in CDU and 213 ± 2 mg C m −2 d −1 (Table 1a; range of grid cells: 6-1,950 mg C m −2 d −1 , Figure 4) in FUN.

Temporal Variability of CH 4 Ebullition From Bubble Trap Measurements
The ebullitive flux as obtained from all bubble traps ( Figures S5 and S6) was significantly higher in 2016 than in 2015 for both CDU (mean ebullition: 297 and 54 mg C m −2 d −1 , respectively; Wilcoxon: p = 1.248e-05, W = 605) and FUN (mean ebullition: 340 and 7 mg C m −2 d −1 , respectively; Wilcoxon: p = 6.064e-12, W = 586.5; Figure S13). Note that while the bubble traps were placed in the same spots for both campaigns in FUN, their distribution in CDU differed ( Figures S5 and S6) such that only the categories inflow and non-inflow can be compared between sampling occasions in CDU. In CDU, CH 4 ebullition was an order of magnitude higher in May 2016 than in April 2015, both in the inflow bays (403 and 62 mg C m −2 d −1 , respectively) and non-inflow bays (363 and 49 mg C m −2 d −1 , respectively; Figure S13). In FUN, CH 4 ebullition was 2 orders of magnitude higher in  we found higher CH 4 ebullition with decreasing depth in the bubble trap data of CDU but no consistent depth pattern for FUN ( Figure S14), similar to the hydroacoustic data ( Figures S9-S12).

Spatial Pattern and Magnitude of Upscaling Approaches
The four upscaling approaches produced different spatial patterns of CH 4 ebullition (Figures 3 and 4), but the resulting estimates of reservoir-wide CH 4 ebullition were remarkably similar between approaches (Table 1a). This indicates that the choice of upscaling method is of lesser importance if spatially well-resolved measurements are available for estimating reservoir-wide CH 4 ebullition, and emphasizes that measurement campaigns with a large number of spatially distributed samples are indispensable for assessing reservoir CH 4 emission.
Approach I (IDW interpolation) is similar to the original data since it interpolates between actual measurements and does not use any observed patterns or relationships, such as depth, to interpolate CH 4 ebullition. Therefore, Approach I is not helpful in inferring values of regions not measured in this reservoir or any other, even if relevant driver data (e.g., depth) were available. Approaches II-IV use observed patterns of CH 4 ebullition from different categories of water depth, either for the entire reservoir (depth groups: Figure S9) or stratified by inflow versus non-inflow areas (depth groups: Figure S10) or hydromorphological section LINKHORST ET AL.
10.1029/2020GB006717 9 of 16 (depth groups: Figures S11 and S12), to derive values of non-measured areas. Because of the underlying relationships with water depth, these approaches resulted in relatively similar spatial CH 4 ebullition patterns (Figures 3 and 4). The similarity between Approaches II and III for both reservoirs, but particularly for CDU, may be attributed to the fact that river inflow areas are generally shallow ( Figure S1), so that water depth and proximity to river inflow are correlated. The similarity between Approaches III and IV, particularly in FUN, may be further related to the observation that CH 4 ebullition predominantly occurred in the main river inflow, which was considered a separate hydromorphological section in Approach IV. Interestingly, the relationship between water depth and ebullition was inverse in the inflow areas of FUN compared to CDU ( Figure S10). We attribute this to the relatively greater depth of the main river inflow to FUN and to the shallow sill that was used to delimit the river inflow area. We suggest that high sedimentation rates of presumably labile organic matter from this inflowing eutrophic river leads to enhanced ebullition in this relatively deep region because the availability of organic matter for extensive CH 4 production and accumulation in the sediments supersedes the reduction in potential ebullition due to increased hydrostatic pressure. Eventually, the basis of all upscaling approaches is the same set of measurements, and hence, the various approaches may lead to different spatial attributions of high-ebullition areas, but not to any noteworthy difference in the total CH 4 ebullition estimate.
In assessing the outcome of upscaling approaches, the magnitude of the resulting estimate is only one criterion; another criterion is how well the upscaling procedure reflects the spatial patterns observed in nature. In CDU, we found a consistent pattern of variation in CH 4 ebullition with bathymetry for both inflow and non-inflow areas ( Figure S10). Therefore, Approach II, III, and IV all led to satisfactory results for CDU, reflecting overall similar spatial patterns of CH 4 ebullition, and returned similar total emission estimates (Table 1a; Figure 3). In FUN, however, CH 4 ebullition varied differently with bathymetry for inflow and non-inflow areas ( Figure S10). We therefore argue that the best way to interpolate whole-reservoir CH 4 ebullition for the FUN data set is to differentiate between inflow and non-inflow areas and then upscale according to depth (Approach III). Dividing the reservoir into different hydromorphological sections (Approach IV) returned similar results as Approach III for FUN, both in terms of magnitude and spatial patterns. In contrast, Approach II, which involved upscaling according to bathymetry but without division of the reservoir, returned high CH 4 ebullition in all bays (Figure 4), which is clearly not what we observed in FUN (Figure 1). It appears that choosing Approach III or IV would have worked for upscaling of both studied reservoirs, with respect to both spatial patterns and the overall magnitude of emission.

Perspective: Application of Upscaling Approaches
With respect to upscaling ebullition measurements in other reservoirs worldwide, we recommend upscaling by bathymetry and river inflows (Approach III), which is a method that incorporates fundamental mechanistic understanding of ebullition (bathymetry: Bastviken et al., 2004;sedimentation patterns: Del-Sontro et al., 2011, Sobek et al., 2012 and is therefore applicable to other reservoirs as well. Approach III returned spatial patterns that were congruent with observed patterns in both CDU and FUN and uses simple relationships between bay volume and sub-catchment size, which can be applied to any reservoir where bathymetry data are available. However, the delimitation of river inflow areas was based on local considerations, and even if the similarity in V/A ratio (0.0110 for CDU and 0.0167 for FUN; Tables S2 and S3) was encouraging, further research is required to assess how far these ratios may be applicable to other systems. Knowledge of the actual sediment deposition could improve the accuracy of determined river inflow delimitation. Unfortunately, sub-bottom profiling did not work for detecting sediment thickness in many regions of the studied reservoirs because of the similar densities of pre-flooding soil and post-flooding sediment in CDU and because of gas bubbles in the sediment acting as strong reflectors that obscured sediment structures in both CDU and FUN. Alternatively, future studies could use hydrodynamic modeling for defining the extent of river inflow areas. Ultimately, the overall CH 4 ebullition estimates from the various upscaling approaches in our study were similar, indicating that the exact placement of the delimitation of inflow areas might not dramatically influence the overall emission estimate and that simpler upscaling approaches (e.g., Approach I, II, or IV) could be used depending on available data.

Reservoir Productivity and Ebullition
FUN had substantially higher fluxes in its main inflow stem than in the rest of the reservoir. While CDU also showed higher fluxes in the inflow bays, the difference to the rest of the reservoir was not as strong as that in FUN. The catchment of FUN is much larger (13,518 km 2 ; Figure S3) than that of CDU (310 km 2 ; Figure S2). Whereas the catchment of FUN contains many urban and agricultural areas, CDU has low population density in its grassland-and forest-dominated catchment, and the source of its main inflow is only ∼50 km away from its impoundment. The catchment characteristics of FUN translate into a high productivity in the reservoir (Pacheco et al., 2015), and FUN probably also experiences a larger sediment load from its catchment compared to CDU. Presumably, a large sediment load leads to a high sediment accumulation rate and ebullition-prone sediments in the main river inflow area of FUN (Berberich et al., 2020;Maeck et al., 2013;Sobek et al., 2012). This suggests that when upscaling CH 4 ebullition, inflows with a particularly large or distinctly different sub-catchment should be treated separately. In reservoirs with small and rather homogeneous sub-catchments, this may not be necessary.
The two studied reservoirs not only differed in terms of patterns related to CH 4 ebullition and thus suitable upscaling approaches, but also in the magnitude of CH 4 ebullition. Mean CH 4 ebullition was 4.4 times higher in FUN than in CDU (214.3 vs. 48.9 mg C m −2 d −1 ; Table 1b). We compare here CDU and FUN from our sampling campaigns in May 2016 and October 2016, but we argue that the different time of sampling is unlikely to affect our conclusion of higher ebullition in FUN compared to CDU. Drops in water level trigger the release of bubbles that have accumulated in the sediment (Harrison et al., 2017;Linkhorst et al., 2020;Varadharajan & Hemond, 2012), but in CDU, the water level was stagnant ( Figure S15), suggesting that the release of bubbles from the sediment was probably at steady state with the production of CH 4 in the sediment. In FUN, we measured during falling water level ( Figure S15), which has been shown to induce a rapid release of accumulated gas bubbles in reservoirs (e.g., Harrison et al., 2017). However, the water level had been falling at a relatively constant rate to a total of ∼10 m in the months preceding our measurements in FUN ( Figure S16). It is thus likely that the rapid release of any gas bubbles accumulated in the sediment were released prior to our measurement campaign, and we can therefore assume that bubble release was at steady state with CH 4 production in the sediment.
Temperature is another important driver of seasonal differences in CH 4 ebullition (e.g., Wilkinson et al., 2015), and air temperature was lower during the CDU measurements in May 2016 than during the FUN measurements in October 2016 (average of 18.2°C and 25.2°C, respectively). Also, bottom water temperatures tend to be lower in CDU than in FUN (19.6°C and 25.7°C in April 2016 in CDU and October 2016 in FUN, respectively; Figure S4), and this relatively small temperature difference could explain a two-fold difference in CH 4 ebullition (using a Q10 of methanogenesis of 4; Bastviken, 2009) between CDU and FUN, but not the observed four-fold difference.
We therefore attribute the difference in CH 4 ebullition to the different trophic states of the two reservoirs, in line with recent evidence suggesting that highly productive systems such as FUN tend to have higher rates of CH 4 emission than less productive systems such as CDU (Beaulieu et al., 2019;Deemer et al., 2016). Accordingly, experimental work has shown that organic matter from aquatic production, that is, phytoplankton or macrophyte debris, can be rapidly and efficiently transformed to CH 4 in anoxic sediments (Grasset et al., 2018(Grasset et al., , 2019. For comparison, we applied the empirical model of CH 4 ebullition from DelSontro et al. (2018), relating CH 4 ebullition to chlorophyll a, to our reservoirs. The model returned a CH 4 ebullition of 35 mg C m −2 d −1 for CDU and of 75 mg C m −2 d −1 for FUN (Table S5). While for CDU, the modeled estimate (35 mg C m −2 d −1 ) is within the same order of magnitude as our estimate (48.9 mg C m −2 d −1 ), it is nearly 3 times lower than our estimate for FUN (modeled: 75 mg C m −2 d −1 ; measured: 214.3 mg C m −2 d −1 ). This difference may be related to the fact that the DelSontro et al. (2018) model only explains ∼30% of the variability in CH 4 ebullition, or to a lack of reservoirs similar to FUN in the database from which the model was derived. The discrepancy may also be related to the extremely limited spatial and temporal resolution of the applied chlorophyll a measurements (see Table S5 for details) that are perhaps not accurately prepresenting productivity in the system. Nevertheless, this exercise corroborates our conclusion that higher productivity in FUN results in higher CH 4 ebullition in FUN than in CDU.

Temporal Variability of CH 4 Ebullition
In this study, we concentrated on capturing spatial variability of ebullition for upscaling purposes, and generally assumed negligible temporal variability in ebullition over the two-week sampling campaigns in each reservoir. However, each variation in space assessed via a mobile echosounder, is inherently also including a variation in time. Accordingly, previous studies in CDU have shown that diel, daily, and seasonal variability in CH 4 ebullition can be substantial (Linkhorst et al., 2020). However, several lines of evidence suggest that the temporal variability that undoubtedly existed during our sampling campaigns was overruled by a stronger spatial variability. In Figure S19, we present all CH 4 ebullition data acquired during the spatial hydroacoustic surveys aggregated by sampling hour and sampling date, and no apparent hourly or daily trends in CH 4 ebullition over the respective 2-week sampling campaigns were visible. With respect to the variability between different sampling days and differences between day and night, we monitored CH 4 ebullition with four bubble traps at fixed locations in one bay of FUN for 6 days ( Figure S20), and sampled the traps every morning (8:30-10:00 h) and every evening (18:00-19:00 h). While higher ebullition was observed on some of the sampling days during daytime (18,19,20,and 21 October), other days showed no difference between day and night emissions (15, 16 October). Thus, these measurements do not indicate a consistent signal of diel variability in the ebullition data. Further, we correlated mean daily ebullition (from hydroacoustic measurements) with daily water level and water level change ( Figure S15), as well as with daily mean atmospheric pressure and temperature ( Figure S21). We also correlated hourly ebullition with hourly mean atmospheric pressure and temperature and with the change in atmospheric pressure from the previous hour to the sampling hour ( Figure S21). The linear regression results are reported in Table S4. Again, we found no significant relationships that could explain any potential hourly and daily trend in CH 4 ebullition, suggesting that temporal variability during our sampling campaigns did not leave a discernible imprint in the ebullition data, and corroborating that the ebullition measurements largely captured only variability in space.
Our previous comparison of spatial and temporal variability in CDU (Linkhorst et al., 2020) suggested to sample each hydrological season (rising and falling water level) for annual upscaling of CH 4 ebullition. In the present study, we used CH 4 ebullition measurements from bubble trap deployments at two different sampling occasions to illustrate the magnitude of seasonal variability. We found that CH 4 ebullition was about 2 times higher during the sampling occasions in 2016 than in 2015, for both CDU and FUN. We speculate this to be connected to pressure drops during falling water level  in FUN as we found severalfold higher CH 4 ebullition during falling water level in October 2016 than during rising water level in March 2015 (Figures S13 and S16). In CDU, the water level was much higher during the campaign in 2016 than in 2015, but nevertheless, it was rather stable at both time points ( Figure S16). Despite the stability, CH 4 ebullition was higher in May 2016 than in March 2015 ( Figure S13), suggesting that factors other than hydrostatic pressure change must have been at play. Regardless of the cause, the difference in CH 4 ebullition between sampling occasions was large and, therefore, the whole-reservoir CH 4 ebullition estimates presented here (Table 1a) are valid for the time of sampling only. For upscaling of spatially resolved surveys to annual whole-reservoir ebullition, we recommend highly resolved spatial surveys during both hydrological seasons of rising and falling water level (Linkhorst et al., 2020).

Daily Reservoir C Emission
Using CH 4 ebullition data of this study and CH 4 diffusion data of subsequent measurements by Paranaíba et al. (2021), the mean areal CH 4 emission that we estimated for CDU (78 ± 1 mg C m −2 d −1 for CH 4 ebullition + diffusion; mean ± SEM) was below the global mean reservoir surface CH 4 emission rate (120 mg C m −2 d −1 for CH 4 ebullition + diffusion; Deemer et al., 2016). The mean areal CH 4 emission estimate for FUN (215 ± 2 mg C m −2 d −1 for CH 4 ebullition + diffusion; mean ± SEM) was at the higher end of reported reservoir emission estimates, as would be expected considering the trophic state of this reservoir (Deemer et al., 2016).
Based on CH 4 ebullition data of this study and CO 2 and CH 4 diffusive emission measurements from Paranaíba et al. (2021), we calculated a total daily reservoir surface C emission ± SEM from CDU (during our sampling in April-May 2016) and from FUN (during our sampling in October 2016) of 35 ± 0.3 and 345 ± 4 t CO 2 -eq d −1 , respectively (Table 1c). The above-mentioned emission values include solely diffusive and ebullitive emission from the reservoir water surface. Other emission pathways such as degassing from the turbines, emission from drawdown zones or through macrophytes may add further emission. For instance, drawdown zone CO 2 emission as reported for CDU by Almeida et al. (2019), may increase our daily C emission estimate of CDU (35 t CO 2 -eq d −1 ) by another 3 t CO 2 -equivalents per day.
According to a global study on reservoirs (Deemer et al., 2016), 82% of C emission expressed as CO 2 -equivalents from reservoirs occur as CH 4 , and 18% as CO 2 . While CH 4 diffusion in our study accounted for 38% and 0.4%, and CO 2 diffusion for 2% and less than 0.1% in CDU and FUN, respectively, CH 4 ebullition accounted for 60% and 99% of total C emission (in CO 2 equivalents) during our sampling campaigns in CDU and FUN, respectively. The fact that more than half of the greenhouse gas emission was emitted via CH 4 ebullition in one reservoir and almost 100% in the other, highlights the need for high-quality CH 4 ebullition estimates for assessing the overall role of reservoirs in the C cycle.
A previous study in FUN had shown total reservoir surface emission via diffusion and ebullition of CH 4 and CO 2 of 30.11 t CO 2 -eq d −1 (dos Santos et al., 2017). Our estimate of 345 t CO 2 -eq d −1 is an order of magnitude higher, which may be due to differing spatial and temporal coverage of measurements. Dos Santos et al. (2017) reported CH 4 emission rates of 0.18 t d −1 (diffusion + ebullition; we assume the number is given as tons of CH 4 ) from an aerial coverage of 14 m 2 (18 deployments of bubble traps with area of 0.69 m 2 ; 34 deployments of floating chambers with area of 0.047 m 2 ), but distributed over several sampling campaigns over two years. In comparison, we report here 10.09 t CH 4 d −1 from ebullition alone (CH 4 ebullition + diffusion: 10.14 t CH 4 d −1 ), based on an immediate spatial coverage by the echosounder beam of 77,067 m 2 but only during one sampling campaign. Our bubble trap data indicate that we probably would have observed lower CH 4 emissions if we had measured during a different season, since mean ebullition rates were 49 times lower during rising water level than during falling water level (7.1 vs. 347.8 mg C m −2 d −1 ; Figure S13), which is when our hydroacoustic survey was conducted (October 2016). Accordingly, reducing our hydroacoustically measured CH 4 ebullition by a factor of 49 returns a flux of 0.21 t CH 4 d −1 via ebullition alone (or ∼0.25 t CH 4 d −1 including diffusion) from FUN during rising water level, which is in the same order of magnitude as estimates from dos Santos et al. (2017). That temporal variability should be considered in upscaling CH 4 emission as much as spatial variability, is shown in a previous study in CDU, where we compared the importance of spatial and temporal CH 4 ebullition variability from monthly hydroacoustic surveys in one bay over one year (Linkhorst et al., 2020). However, given that the measurements by dos Santos et al. (2017) were gathered during four campaigns, it seems unlikely that all of these four campaigns were during low-ebullition periods, and more likely that their spatial coverage was too low to capture areas with high ebullition, in line with our analysis of minimum sampling effort ( Figures S17 and S18).

Importance of Adequate Spatial Coverage
This is an unprecedented study in terms of spatial coverage of CH 4 ebullition. From the total survey lengths of 34.6 and 63.2 km (and mean beam widths of 1.68 and 1.22 m), we calculated an area of 58,009 and 77,067 m 2 that was covered by the echosounder beam in CDU and FUN, respectively. In terms of total hours of echosounder surveying, this translates in 27 h for CDU, and 39 h for FUN. While the immediate coverage corresponds to only 0.6% and 0.2% of the total surface area of CDU and FUN, respectively, this coverage is much higher than in previous studies in tropical reservoirs. Further, our measurements covered a substantial share of the spatial heterogeneity of CH 4 ebullition as we focused on covering different hydromorphological regions of the reservoirs that have been described as particularly high or low in CH 4 ebullition (i.e., shallow and deep areas, river and non-river inflow areas, areas in which sedimentation rate was considered high or low; DelSontro et al., 2011;McGinnis et al., 2006;Sobek et al., 2012). Also, the similarly large range in ebullition measured on each day of measurement, spanning over 3-4 orders of magnitude ( Figure S19), indicates that we probably captured a significant share of the total range of ebullition in these reservoirs.
We found that the minimum sampling effort to reach an accurate CH 4 ebullition estimate is substantial. In order to reach within ±20% of our mean CH 4 ebullition as calculated from all 30-m segments, at least 1,700 locations need to be sampled in CDU and at least 1,900 locations in FUN, corresponding to roughly 65% of all the locations covered in our survey (Figures S17 and S18). Furthermore, we found that surveying fewer than 250 locations with an area of ∼45 m 2 (30 m hydroacoustic segments times average beam width of 1.5 m) will lead to less than a 50% chance of reaching accurate CH 4 ebullition estimates. Most importantly, this clearly shows that attaining good estimates of reservoir CH 4 ebullition require substantial sampling efforts and high spatial coverage. Our study found that while CH 4 ebullition fluxes covered the same range for inflow and non-inflow areas ( Figure S19), mean CH 4 ebullition was two-fold higher in inflow areas compared to non-inflow areas; therefore, the sampling effort should cover both types of regions to a similar extent. That CH 4 ebullition made up 60% and 99% of total C emissions (in CO 2 -eq) in our two reservoirs, shows that it is imperative to adequately resolve CH 4 ebullition when calculating total C emissions and estimating the climate impact of (tropical) reservoirs.

Conclusions
This study supports the demand for greenhouse gas emission studies of reservoirs to be performed at high-spatial resolution and during different seasons. Greenhouse gas emission from reservoir water surfaces is currently estimated to account for 1.3% of the global anthropogenic C emission. In order to make this number more reliable as reservoir construction for water storage is gaining importance in face of Climate Change, and as mankind is increasingly turning to hydropower as a renewable and supposedly sustainable energy source, our study clearly shows that more studies that adequately resolve reservoir greenhouse gas emissions are needed. This applies particularly to CH 4 ebullition, which accounted for 60% and 99% of total reservoir surface C emission from the two reservoirs in our study.
The density of ebullition measurements was much more important than the technique used to upscale to system-wide estimates, with at least 250 measurement locations (i.e., segments of 30 m hydroacoustic survey, for reservoirs with water surfaces of 9.5-35 km 2 ) required for a 50% probability to reach accurate CH 4 ebullition estimates. This is a nearly impossible number of locations to achieve with manual bubble traps, meaning that spatially resolved measurement techniques such as hydroacoustics are preferable; however, manual bubble traps remain useful for surveying temporal variability or in shallow regions (<2 m) that cannot be adequately covered by echosounders. Along with the necessary density in ebullition measurements required for accurate estimates, one should consider major differences in reservoir heterogeneity that may promote substantial differences in ebullition magnitude. It is vital to consider proximity to river inflows and dams, as well as depth relationships with ebullition, when developing surveys that are representative of the true CH 4 ebullition potential of the reservoir. In combination with a probabilistic sampling design such as applied by Beaulieu et al. (2016), in which different types of areas are surveyed proportionally to their coverage area in the reservoir, the accuracy of the outcome can be further improved.
This study focused on two reservoirs in the (sub)tropical region of Brazil, but our results can be relevant for reservoirs in other climates and altitudes. While the magnitude of fluxes will be different and largely regulated by temperature and supply of reactive particulate organic C (i.e., eutrophication), the methods and approaches we have applied here are also useful in other systems.

Data Availability Statement
Data that were used to generate the results of this article can be retrieved from DiVA, Uppsala University's platform for digital publishing (http://uu.diva-portal.org/smash/record.jsf?pid=diva2%3A1512421).