We obtained data for meta-analysis from published and unpublished sources (Supplementary Table 1). We conducted literature searches using ISI Web of Science and by personal communication. Data selection continued until September 2016. The studies included in our assessment all included a direct comparison of the enriched water and bubble methods with common samples. In all, 13 studies met the requirements for inclusion (n = 368 observations) along with two unpublished data sets (Supplementary Table 1, Supplementary Figure 1). ¹⁵N₂ fixation rates were either provided personally by authors, retrieved from tables published in supplemental material or manually digitized from figures in the published studies using the software WebPlotDigitizer (version 3.10) (Ankit Rohatgi, http://arohatgi.info/WebPlotDigitizer/).

Before performing the meta-analysis, we determined the statistical dispersion and variability within the bubble and enriched water data sets. We calculated the mean absolute deviation (MAD) of each data set as follows:

where n is the number of observations and x- is the mean of the individual observations x_i.

Meta-analysis and meta-regression were conducted with R 3.1.2 using the “metafor” package (Viechtbauer, 2010), metafor package of R [http://www.metafor-project.org/], R Development Core Team 2013). To assess the two methods, we used the logarithmically transformed response ratio (Hedges et al., 1999), lnRR, where RR is the ratio of rates (R) measured with the enriched water and bubble methods (RR = R_{enrichedwater}:R_bubble), or the effect size of individual experiments), as well as the corresponding pooled standard deviation. All lnRR values were weighted by the reciprocal of their sampling variance, followed by a random effects model to compute the overall mean effect size, which is equivalent to Cohen's d parameter (Cohen, 1977). The random effects model, based on the DerSimonian-Laird estimator (DerSimonian and Laird, 1986) calculates the between-study variance (σ²) and weights each study by the inverse sum of the individual study variance (ν_i) and the between- study variance. Mean effect sizes (i.e., Cohen's d) were statistically significant different if their 95% confidence intervals did not overlap zero.

Apart from calculating the mean effect sizes for all 13 studies and 368 observations, we also divided observations into two groups: short incubation time (0–12 h) and long incubation times (24 h). For short incubation times three data sets were available (Mohr et al., 2010; Klawonn et al., 2015; Benavides and Wannicke et al., unpubl.,). Unfortunately, the experiments of Mohr et al. (2010) and Klawonn et al. (2015) lacked replication from single time points. To be able to include this data set in our meta-analysis, we pooled values collected at t = 0 to 12 h. Results from these publications for t = 24 h had to be omitted from the meta-analysis due to lack of replication. Furthermore, we added a subgroup analysis for each meta-analysis to determine mean effect sizes excluding unpublished data sets.

To explore heterogeneity in the meta-analysis, we calculated Cochrane Q-tests for heterogeneity (Cochran, 1954). A significant Q-statistic (p-value of < 0.1) indicates that there is heterogeneity within the mean effect size and a greater variance among individual effect sizes than expected by sampling error. We furthermore, examined the sensitivity of the meta-analysis by examining publication bias (the probability that statistically significant (p < 0.05) results are more likely to be published than non-statistically significant results) using a contour-enhanced funnel plot (Peters et al., 2008; Supplementary Figure 2). The funnel plot represents a scatter plot of the effect estimates from individual studies against the standard error of the effect estimate. Specifically, a contour-enhanced funnel plots display the area of statistical significance on a funnel plot (Peters et al., 2008) to improve the correct identification of the presence or absence of publication bias. Publication bias would be expected when the usual funnel plot is asymmetrical.

Random-effects meta-regression analysis using a linear mixed-effects model was used to evaluate the association between incubation time and the mean effect size of ¹⁵N₂ fixation.

For interpretation of the meta-analysis we converted Cohen's d value to Cohen's U₃ parameter (Cohen, 1977) to give a measure of the degree of separation (i.e., % non-overlap) of data produced by the bubble and enriched water methods according to:

where Φ is the cumulative distribution function of the standard normal distribution, and δ is the population value of Cohen's d.

In addition, we calculated the overlapping coefficient (OVL) of data from the two methods by converting Cohen's d using the following formula (Reiser and Faraggi, 1999)

where Φ is the cumulative distribution function of the standard normal distribution, and δ the population Cohen's d.

The probability of superiority (CL) (Ruscio and Mullen, 2012) i.e., probability that a sample picked at random from the treatment group will have a higher score than a sample picked at random from the control group, was calculated the following (Ruscio, 2008)

where Φ is the cumulative distribution function of the standard normal distribution, and δ the population Cohen's d.

Dissolution of ¹⁵N₂ gas according to Montoya et al. (1996) resulted in a maximum ¹⁵N atom% enrichment of 9.1 % (mean value for 24 h of incubation, Figure 1). Fitting of data resulted in an atom% enrichment of 7.8% (solid line, Figure 1). The time to reach 50% of the maximum ¹⁵N atom% enrichment was 1.7 h and the shift from an exponential rise in dissolved ¹⁵N enrichment to a plateau with only minimal changes in ¹⁵N atom% enrichment (>60% of maximum) occurred after 4 h. After 8 h, 90% of ¹⁵N atom% equilibration was reached. Injection of water pre-enriched with ¹⁵N₂ gas according to Mohr et al. (2010) resulted in a stable ¹⁵N atom% enrichment over 24 h with a mean value of 7.5 ± 0. 9% (Figure 1).

When using the bubble method, the rate of ¹⁵N₂ fixation will be systematically underestimated during the equilibration phase of ¹⁵N₂ gas with the dissolved pool of N₂ (Figure 2A). The overall magnitude of the underestimate during an incubation will depend on the length of the incubation and the timing and duration of ¹⁵N₂ fixation activity during the incubation period (Figure 2B), which in turn will reflect the nature (light depended or light independent) of the diazotrophic community present. The maximum error of −72% occurs when the time lag between gas injection and start of the ¹⁵N₂ fixation is zero (T_i = 0 h) and the duration of ¹⁵N₂ fixation is 1 h (T_f = 1 h, red dot in Figure 2A). In contrast, a 12 h period of active fixation beginning 6 h after the injection of ¹⁵N₂ gas will result in an error of −0.2% (green square in Figure 2A, green solid line in Figure 2B, T_i = 6 h and T_f = 12 h).

When organisms are able to fix N₂ continuously, the error would be −12% (blue diamond in Figure 2A, T_i = 0 h and T_f = 12 h) over a 12 h incubation and −6% over a 24 h incubation (T_i = 0 h and T_f = 24 h).

Variability among the replicates of individual studies was quite large as indicated by the standard deviation for the single studies and MAD of the two methods applied (Figure 3). Specifically, the MAD for the enriched water method is 7 and the MAD for bubble method is 5.

To determine a mean effect size (i.e., Cohen's d), we firstly resolved the response ratio (lnRR) of N₂ fixation for each study (Figure 4A). Subsequently, we calculated the overall mean effect size across all studies which resulted in a significant Cohen's d of 0.406 ± 0.110, (p < 0.001; Figure 4B). The mean effect size of ~0.4, detected in our analysis suggests that for a method comparison, 65 % of measurements from the enriched water group will be above the mean of the bubble group (Cohen's U₃, for subgroup analysis 73%) and that 84% of the measurements in the two groups will overlap (subgroup analysis 76%). In addition, there is a 61% chance that a sample picked at random from the enriched water group will have a higher value than one picked at random from the bubble group (subgroup analysis 66%).

A subgroup analysis excluding the unpublished data sets revealed a significant mean effect size of 0.631 ± 0.125 (p < 0.001).

We also performed a meta-regression to evaluate the influence of time on the overall mean effect size. No significant impact was detected (data not shown). Moreover, we checked for publication bias, which is expected when scatter of data in the funnel plot is asymmetric. In our analysis, assessment of the contour-enhanced funnel plot indicates an asymmetrical scatter of data and potential publication bias introduced by a lack of non-significant studies (Supplementary Figure 2).

We furthermore, did separate meta- analysis for observations with short incubation times (0–12 h) and long incubation times (24 h only, Supplementary Figures 3A,B). The overall mean effect size of the three studies with short incubation time was not significant with a Cohen's d of 0.057 ± 0.188, (p > 0.1) (Supplementary Figure 3A). Excluding the unpublished data set of Benavides and Wannicke et al., revealed a non-significant mean effect size of 0.557 ± 0.382 (p > 0.1). Furthermore, the meta- analysis for long incubation times of 24 h resulted in a Cohen's d of 0.406 ± 0.101 (p < 0.001) (Supplementary Figure 3B). Excluding unpublished studies resulted in a mean effect size of 0.68 ± 0.136 (p < 0.001).

We found no significant correlations with ocean province or temperature in exploring which factors might influence the mean effect size (data not shown).

Our error estimation of the bubble method during a 24 h experiment reveals a negligible error of −0.2% assuming a diazotroph community that fixes only during 12 h daytime and starting of nitrogen fixation 6 h after the injection of ¹⁵N₂ gas i.e., the addition was done 6 h before sunrise. Considering that the error introduced by using a gas-tight syringe to inject ¹⁵N₂ gas of ±1% (according to the manufacturer, Hamilton USA), the error introduced by using bubble injection is insignificant. The error introduced by using the bubble method will increase to −6% when diazotrophs fix continuously over 24 h (assuming there is no time lag between bubble injection and start of active fixation). Overall, it is important to adjust incubation times relative to the onset of active N₂ fixation (which is in turn is depended on the dominating diazotrophs present), as has been indicated before (Mohr et al., 2010; Wilson et al., 2012).

Statistical dispersion, as represented by the MAD from the mean value, appears to be higher in the data set based on the enriched water method (mean MAD of all studies 7) than in the data set of studies that used the bubble method (mean MAD of all studies 5). That is, measured N₂ fixation rates appear to be more consistent when determined using the bubble method. A larger dispersion of data in experiments using the enriched water method might be introduced in by the process of preparing the enriched water for later usage, i.e., degassing water of different volumes and varying accuracy of degassing. In addition, Wilson et al. (2012) presented evidence for a dispersion of data in experiments using the enriched water method, which was connected to abiotic factors influencing the incubation itself. They noted an interesting contrast between the two methods when comparing samples incubated at sea using either an in situ array or shipboard incubators. Specifically, they found that the enriched water method produced depth-integrated N₂-fixation estimates that were 30% greater when samples were incubated aboard ship in deck incubators than when incubated on an in situ array. In contrast, depth-integrated rates measured using the bubble method were not significantly different between incubations carried out on deck and on an in situ array. The 30 % difference among replicates of the enriched water method also overlaps with the divergence of the enriched water and bubble methods referred to by Mohr et al. (2010) and Wilson et al. (2012) in their methods assessments.

Our meta- analysis revealed a large congruence in the estimates of ¹⁵N₂ fixation rate produced using the two experimental methods. The 84% overlap of rate estimates make it very difficult at this stage to estimate any sort of a global factor to quantify the degree of underestimation of ¹⁵N₂ fixation rates when using the bubble method. Our literature review moreover revealed that a thorough comparison over a 24-h cycle is needed with only three studies on sort incubation times of 0–12 h. A larger comparative analysis is clearly necessary especially in view of the elevated dispersion (i.e., larger MAD) of ¹⁵N₂ fixation rates measured using the enriched water method.

As Großkopf et al. (2012) have pointed out, there is an indication for a species-specific potential for underestimation of N₂ fixation rates using the bubble method, especially when buoyant diazotrophs are presented. Thus, in habitats dominated by filamentous species like Trichodesmium or Nodularia, underestimation seems to be less severe compared to habitats dominated by unicellular species (UCYN, e.g., Zehr et al., 1998), as well as non-diazotrophic Bacteria and Archaea (e.g., Riemann et al., 2010). In the latter, incubation times have to be adjusted to the equilibration time of ¹⁵N₂ when using the bubble method. Alternatively, as proposed in the sub-chapter below (see “Experimental Recommendations”) the determination of the final ¹⁵N₂ (i.e., substrate) enrichment in the incubation bottle enables a concerted calculation of N₂ fixation rates.

In our analysis, we have included two data sets that are currently unpublished (Benavides and Wannicke et al., Fabian et al., Supplementary Table 1). The literature dealing with the inclusion of unpublished data in meta-analysis clearly recommends inclusion of this sort of “gray literature” (e.g., Cook et al., 1993; McAuley et al., 2000). Exclusion from meta-analysis can lead to exaggerated overestimation of treatment effects (e.g., Cook et al., 1993; McAuley et al., 2000). This is also reflected in our analysis of funnel plot symmetry, where studies displaying no significant effect (Cohen's d of ~ 0) are under-represented and/or missing in our analysis. Overall, publication bias is the greatest threat to the validity of meta-analysis, because combining only the identified published studies uncritically, may lead to an incorrect, unusually one sided, conclusion.

A number of experimental factors have a strong influence on the precision and accuracy of the determination of N₂-fixation rates. Firstly, the sensitivity of any experiment using ¹⁵N₂ depends on the amount of tracer added to the dissolved pool of N₂. For example, addition of 1 mL of ¹⁵N₂ per liter of sample will produce an equilibrium enrichment of ~5–10 atom% ¹⁵N, depending on the size of the ambient pool of N₂. In contrast, the procedure proposed by Großkopf et al. (2012) results in an enrichment of only 2 atom%. Since N₂ availability does not limit N₂-fixation activity, greater additions may easily be realized to increase the substrate labeling, thereby increasing the sensitivity of the rate measurement. This is especially important in systems where rates are expected to be low, for example in aphotic deep waters. In general, we recommend adding sufficient ¹⁵N₂ to raise the ¹⁵N content of the dissolved N₂ pool to 9–10 at% as noted by Montoya et al. (1996).

Secondly, the natural variability of δ¹⁵N of the particulate nitrogen (PN) pool sets a lower limit to rate measurements. If the variability in δ¹⁵N of PN is high at the start of the incubation (t₀) and the final increase in δ¹⁵N values of the PN in the incubation bottles is low due to low N₂ fixation rates, then N₂-fixation activity may not be detectable. For example, Wasmund et al. (2015) have nicely explored this issue, discussing ¹⁵N₂ fixation rate measurements in the Benguela upwelling region where they compared initial (t₀) and final δ¹⁵N measurements of samples incubated with ¹⁵N₂. The mean values of the two batches of filters differed only by 0.9‰, leading Wasmund et al. (2015) to conclude that ¹⁵N₂ fixation rates were too low to resolve with the tracer method. Nowadays, mass spectrometers clearly perform analytical precisions of 0.2‰ and better. Therefore, the detection limit of enriched PN is well below 4‰, as original proposed in the paper by Montoya et al. (1996).

Finally, the two experimental approaches differ fundamentally in the degree and nature of experimental manipulation of the sample. The bubble method involves minimal handling (a thourough mixing of water and gas bubble after injection has to be guaranteed, by using e.g., a continuously rotating) and perturbation of the system, but can lead to a systematic underestimate of N₂-fixation rate if a significant fraction of the overall activity during the experiment occurs during the isotopic equilibration phase. In contrast, the enriched water method requires extensive processing in advance to prepare the ¹⁵N₂-labeled water used to inject tracer into the experimental bottle as described by Klawonn et al. (2015) and Mulholland et al. (2012). Ideally, the enriched water should be obtained from the same location and depth as the experimental sample, which imposes a significant cost in time and handling while setting up the experiment, with the added risk of contamination with ammonium, DON, or other dissolved species during handling. If the enriched water is prepared in advance from artificial seawater or water obtained at a different station, it represents an addition of alien water and dissolved species and should be minimized in volume to less than a few percent of the bottle volume. The potential for contamination and the degree to which the addition of enriched water may affect estimates of N₂ fixation rates are very difficult to constrain but can clearly compromise the reliability of the final rate estimates. Another approach using the bubble addition followed by the removal of the bubble after only few hours and retrieval of a subsample for determination of ¹⁵N₂ atom% enrichment requires the tedious determination of the ¹⁵N₂ atom% enrichment for each incubation bottle (Jayakumar et al., 2017). Although this approach is easy to handle, it requires efficient handling and subsequent stable isotope analysis.

Our dissolution experiments investigating the isotopic equilibration in seawater along with the theoretical error calculation both suggest that incubation times longer than about 6 h are minimally affected by the equilibration of the added ¹⁵N₂ gas and the dissolved pool of N₂ in the experimental bottle (Figures 3A, 4A). In practical terms, 24 h incubations are frequently recommended and used because these experiments integrate over a full day/night cycle.

A final recommendation to improve the accuracy of N₂ fixation measurements and potentially help resolve the source of variability among replicates is collection and preservation of a water sample from each experimental bottle for determination of the final ¹⁵N₂ (i.e., substrate) enrichment. This would improve the accuracy of the enriched water, as well as the bubble method, both of which typically rely on solubility calculations to estimate the size of the ambient pool of N₂, which in turn determines the actual ¹⁵N₂ enrichment of the dissolved pool.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.