Package 'HMR'

Trace gas flux estimation with static chamber data

Description

HMR reads static chamber data from a semicolon/comma separated ASCII text file and analyzes the selected data series by either user selected (with decision support) models or by user configured automatically selected models. Results are exported to a semicolon/comma separated ASCII text file.

Usage

HMR(filename, series = NA, dec = '.', sep = ';', SatPct = NA, SatTimeMin = NA,
pfvar = NA, pfalpha = 0.05, LR.always = FALSE, FollowHMR = FALSE,
IfNoValidHMR = 'No flux', IfNoFlux = 'No flux', IfNoSignal = 'No flux')

Arguments

filename	The name of an HMR data file. It is assumed that the data file folder has previously been set by the setwd command, eg. setwd('C:/My HMR applications'), or through the 'File|Change dir...' menu. The required fixed format of HMR data files is described in the details section below.
series	Data series in the data file to be analyzed by HMR. The default, series=NA, requests analysis of all available data series, whereas eg. series=c('A','B') requests analysis of series named 'A' and 'B' in the data file.
dec	The decimal separator used in the data file. Options are dec='.' and dec=',' and the default is dec='.'.
sep	The column separator used in the data file. Options are sep=';' and sep=',' and the default is sep=';'.
SatPct	The chamber saturation percentage assumed to be met no earlier than SatTimeMin after deployment; see details on flux limiting below. The percentage must be >0 and <100 or NA (no saturation assumption). The default is SatPct=NA.
SatTimeMin	The earliest time for chamber saturation SatPct assumed for for flux limiting; see details below. The earliest saturation time must be >0 or NA (no saturation assumption). The default is SatTimeMin=NA.
pfvar	The assumed variance of replicate measurements of the ambient trace gas concentration at the chamber site in case of no gas emission used for prefiltering of data series before analysis; see details below. The variance must be >0 or NA (no prefiltering). The default is pfvar=NA.
pfalpha	The statistical significance level used for the prefiltering test; see details below. Data series with prefiltering p-value less than pfalpha are classified as 'Noise', and pfalpha must be >0 and <1. The default is pfalpha=0.05.
LR.always	A logical indicating whether data series should be analysed by linear regression in addition to the selected analysis. The default is LR.always=FALSE.
FollowHMR	A logical indicating whether or not automatic model selection should be applied. This cancels all plotting and interactive user selection of analysis. The automatic selection can be configured by options IfNoValidHMR, IfNoSignal and IfNoFlux, and interacts with flux limiting and prefiltering as described below. The default is FollowHMR=FALSE.
IfNoValidHMR	With automatic model selection, ie. FollowHMR=TRUE, the analysis specified by IfNoValidHMR is performed when no valid non-linear HMR model can be fitted; see further details below. Options are 'LR' and 'No flux' and the default is IfNoValidHMR='No flux'.
IfNoFlux	With automatic model selection, ie. FollowHMR=TRUE, the analysis specified by IfNoFlux is performed when the HMR function recommends 'No flux'; see further details below. Options are 'LR' and 'No flux' and the default is IfNoFlux='No flux'.
IfNoSignal	With automatic model selection and prefiltering, ie. FollowHMR=TRUE and pfvar are specified, the analysis specified by IfNoSignal is performed when the prefiltering test classifies the data series as 'Noise'; see further details below. Options are 'LR' and 'No flux' and the default is IfNoSignal='No flux'.

Details

HMR basic methodology

The HMR function implements the methods introduced in Pedersen et al. (2010) supplemented with methods for flux limiting, noise/signal prefiltering and options for automatic model selection investigated by Pullens et al. (2023).

$\; \; \; \; \;$ The HMR approach is based on the exponential model ('HMR'):

$C_t=\varphi+f_0 \frac{\exp(-\kappa t)}{-\kappa h}$

Here $C_t$ denotes the chamber concentration at time $t$, $\varphi$ denotes the chamber equilibrium concentration, $f_0$ denotes the initial flux to be estimated, $\kappa>0$ is a non-linearity parameter, and $h=V/A$ denotes the effective chamber height, ie. the chamber volume divided by the cross-sectional area. The HMR model can be derived from the diffusion model suggested by Hutchinson and Mosier (1981) (augmented to allow for horizontal gas transport by Pedersen et al. 2010) or from a standard two compartment model (Seber and Wild 1989, ch. 8) with first order trace gas transport in and out of the chamber compartment (in: emission; out: leakage through or horizontal transport below chamber walls). It can also simply be seen as a non-linear model with exponential curvature used for fitting exponentially curved concentration data and with:

$f_0=h\frac{dC_t}{dt}|_{t=0}$

The HMR approach does, however, insist on parameter values that correspond to positive fitted values of the initial and the equilibrium concentration, ie. $C_0$ and $\varphi$. Further details can be found in Pedersen et al. (2010). For small values of $\kappa$, the HMR model is close to the linear model ('LR'):

$C_t=\varphi+f_0\frac{t}{h}$

For large values of $\kappa$, it is close to the constant model ('No flux'):

$C_t=C_0$

The HMR function estimates the model parameters by minimzing the mean squared error (MSE) criterion. If the optimal value of $\kappa$ is indeterminate by the MSE criterion, ie. the MSE continues to decrease for $\kappa$ approaching either zero or infinity, the HMR function instead recommends the corresponding limiting model for flux estimation, ie. 'LR' and 'No flux', respectively. In these limiting cases, the estimation uncertainty in the estimate of $\kappa$ (zero or infinity) is not reflected in the standard error of the estimated flux (nor the 95%-confidence interval or the p-value) implying a discontinuity when going from very small/large estimated values of $\kappa$ in the non-linear model ('HMR') to the limiting linear ('LR') or constant model ('No flux').

HMR flux limiting

The family of non-linear HMR models is clearly too large as all values of $\kappa>0$ are permissible and large values of $\kappa$ may correspond to huge and unrealistic fluxes. With small data series, such excessive fluxes may appear optimal by the MSE criterion when random patterns in data, erroneously, exhibit strong exponential curvature, so it may be important to limit the permissible values of $\kappa$. A priori, $\kappa$ can not be bounded by any fixed number, as it depends on the time unit, the chamber design, local characteristica and more, but it may be limited by making assumptions about the chamber saturation. The non-linear model can be rewritten like:

$C_t=C_0+(\varphi-C_0)(1-\exp(-\kappa t))$

Hence, $1-\exp(-\kappa t)$ denotes the saturation fraction at time $t>0$ after deployment. According to the model, complete chamber is not obtained until $t=\infty$, but for a given saturation fraction, $0<p<1$, and assuming that this saturation is not obtained before time-point $T$ after deployment, an upper limit for $\kappa$ is implied by the equation:

$\kappa \leq \frac{1}{T} \log(\frac{1}{1-p})$

Even very crude assumptions about the chamber saturation may limit the permissible values of $\kappa$ in an important way. Example: As complete chamber saturation does not happen within finite time in the theoretical model, one may use eg. 99% saturation in the model as a proxy for what is considered complete chamber saturation in real measurement situations. Assuming that this saturation will not occur until after $T=1$ or $T=2$ limits the permissible values of $\kappa$ to $\kappa\leq 4.6$ and $\kappa\leq 2.3$, respectively.

$\; \; \; \; \;$ Flux limiting with HMR is controlled by options SatPct and SatTimeMin that specifies the values of the saturation fraction (in %) and the assumed earliest time-point for the saturation percentage, respectively. If flux limiting is active, ie. the globally optimal value of $\kappa$ according to the MSE criterion is larger than the user specfified maximal value of $\kappa$ (whereby the MSE criterion recommends 'No flux'), then 'LR' should be considered for estimating the flux as the mathematical justification for 'No flux' is violated.

HMR data prefiltering

The MSE criterion is unable to detect if the variation in chamber concentrations is low relative to the variation in data that can be expected merely from variation in ambient trace gas concentration measurements at the site in case of no trace gas emission. Hence, and particularly with small samples, the MSE criterion may be fooled by random patterns in data and estimate a, potentially large, flux even if there is none. This is particularly problematic with the non-linear HMR model. To avoid this, external information is needed about the expected variation in measured chamber concentrations at the site in case of no flux. This variation depends on the analytical laboratory precision and on the natural variation of true trace gas concentrations in replicate samples from the site in case of no flux. Although this abstract quantity is presumably unknown, assessment of its value may prove useful for prefiltering of data series. For instance, the requested variability could be estimated from replicated samples of ambient air before chamber deployment, although this may lead to an overestimation due to the potential presence of trace gas emission at the site. Alternatively, one may consider the natural variation of true trace gas concentrations in replicate samples from the site to be negligible, in which case the requested variability is simply the analytical laboratory precision. This may of course underestimate the variability. In this way the abstract variability measure may be framed. If the frame is narrow, an intermediate value may be a good assessment of the required abstract variability, whereas the two framing values can be used for, respectively, conservative (overestimated variability: classifies too few data series as 'Signal') or liberal (underestimated variability: classifies too many data series as 'Signal') prefiltering.

$\; \; \; \; \;$ Given that a reasonable measure of the variation of replicate concentrations measurements at the site without flux can be provided, this can be used for prefiltering of data series, ie. the classification of series into 'Noise' or 'Signal'. The idea is then to avoid erroneous flux estimation for data series identified as 'Noise'. The HMR prefiltering test assumes that replicate concentration measurements at the site without trace gas emission follow a normal distribution with variance $\sigma^2_0$. The prefiltering classification is performed by a one-sided statistical test of the null hypothesis of no variation in data in excess of what can be expected from concentration measurements from the site in case of no flux against the alternative hypothesis of larger variation in data than can be expected at the site in case of no flux: The prefiltering p-value is given by:

$p=1-F_{\chi^2(n-1)}((n-1)\frac{s^2}{\sigma^2_0})$

Here $n$ denotes the sample size, $s^2$ denotes the sample variance of measured chamber concentrations, and $F_{\chi^2(n-1)}$ denotes the cumulative distribution function of the $\chi^2$-distribution with $n-1$ degrees of freedom.

$\; \; \; \; \;$ With prefiltering at statistical significance level $\alpha$, the HMR function classifies data series as 'Signal' if $p<\alpha$ and otherwise as 'Noise'. The level of prefiltering can be controlled by the selected values of $\sigma^2_0$ and $\alpha$. With statistical significance level $\alpha$, $100\alpha\%$ of data series will, erroneously, be classified as 'Signal' due to statistical type I errors, so increasing $\alpha$ will reduce the noise filtering rate for a given value of $\sigma^2_0$. Increasing the assumed variation of concentration measurements at the site without flux, ie. $\sigma^2_0$, increases the noise filtering rate as more data series will be classified as 'Noise'.

$\; \; \; \; \;$ Prefiltering with HMR is controlled by options pfvar and pfalpha that specify $\sigma^2_0$ and $\alpha$, respectively.

HMR user selection of analysis

To assist the user in selecting the appropriate method for flux estimation for a given data series, the HMR function displays – organised in a 2x2 matrix – plots of the criterion function, the various models fits, decision support from the flux limiting and prefiltering methods, and a panel with buttons for user selection of method for flux estimation. The upper left plot displays the criterion function over the range of numerically feasible values of $\kappa$ – possibly further limited by the user through flux limiting, cf. sections above. Green parts of the curve represent valid values of $\kappa$, red parts represent invalid values (cf. above), and the optimal value according to the MSE criterion is indicated by a blue square. The upper right plot is a zoom into the upper left plot, and the lower left plot shows the fit of the possible models with the model selected by the MSE criterion displayed in the headline. Moreover, decision support from the prefiltering and flux limiting methods is written in red text in this plot. The lower right panel contains buttons for user selection of the method for flux estimation (select by left mouse button click). Pressing the cancel button interrupts the HMR function.

HMR automatic selection of analysis

Although users are encouraged to do sequentially user selection of model for flux estimation, HMR does facilitate automatic selection of analysis (FollowHMR=TRUE). The automatic selection is based on exactly the same decision support, however, with consistent unsubjective choices configured by the user through options IfNoValidHMR, IfNoFlux and IfNoSignal (cf. above). Note that IfNoSignal has precedence over IfNoValidHMR and IfNoFlux when relevant, and that 'LR' is selected irrespective of the user selection for IfNoFlux if flux limiting is active, ie. the globally optimal value of $\kappa$ according to the MSE criterion is larger than the user specfified maximal value of $\kappa$ (whereby the MSE criterion recommends 'No flux'). The automatic HMR decision tree with option FollowHMR=TRUE is illustrated in the figure below (only visible with HTML help, help(HMR,help_type='html'), and in the PDF manual).

HMR data processing order

Firstly, HMR analyzes if the data series can be fitted by valid HMR models. If not, user selection between 'LR' and 'No flux' is supported by results from the prefiltering test if selected, and automatic model selection is controlled by the user selected values of IfNoValidHMR and IfNoSignal. If valid HMR models can be fitted, the range of values of $\kappa$ is potentielly further limited by user selected flux limiting assumptions. User selection of models is then supported by the prefiltering and flux limiting methods, if selected, whereas automatic model selection is controlled by the user specified values of IfNoFlux and IfNoSignal.

HMR data files

HMR data files are semicolon or comma separated files organised in five columns with, respectively, the data series names, chamber volumes, chamber cross-sectional areas, observation time-points, and the observed chamber concentrations. Semicolon/comma separated files can for instance be created and edited by ASCII text editors or spreadsheet software. The following 36 lines shows an HMR data file containing four data series of static chamber data:

Series;V;A;Time;Concentration
k0d;140.625;0.5625;0;15.60
k0d;140.625;0.5625;10;15.62
k0d;140.625;0.5625;20;16.53
k0d;140.625;0.5625;30;16.90
k0d;140.625;0.5625;40;17.40
k0d;140.625;0.5625;50;17.69
k0d;140.625;0.5625;60;18.64
k0d;140.625;0.5625;70;18.36
k0d;140.625;0.5625;80;19.14
k0d;140.625;0.5625;110;18.83
k0d;140.625;0.5625;120;19.27
k10d;140.625;0.5625;0;0.3517
k10d;140.625;0.5625;10;0.3523
k10d;140.625;0.5625;20;0.3660
k10d;140.625;0.5625;30;0.3673
k10d;140.625;0.5625;40;0.3603
k10d;140.625;0.5625;50;0.3623
k10d;140.625;0.5625;60;0.3580
k10d;140.625;0.5625;70;0.3650
k10d;140.625;0.5625;80;0.3700
k10d;140.625;0.5625;90;0.3673
k10d;140.625;0.5625;110;0.3647
k10d;140.625;0.5625;120;0.3693
F2T2;2.0101;0.0201;0;10.87
F2T2;2.0101;0.0201;20;19.49
F2T2;2.0101;0.0201;54;24.99
F2T2;2.0101;0.0201;85;27.24
F2T2;2.0101;0.0201;119;33.13
F2T2;2.0101;0.0201;155;30.14
F2V2;2.0101;0.0201;0;9.940
F2V2;2.0101;0.0201;28;31.64
F2V2;2.0101;0.0201;60;48.88
F2V2;2.0101;0.0201;91;58.08
F2V2;2.0101;0.0201;123;76.16
F2V2;2.0101;0.0201;162;106.8

Apart from the required header (with optional ASCII character content), the five columns contain:

Column 1

Text labels that identify data series. Hence, labels must be identical within and different between data series. In the sample data above, the first column identifies four data series named k0d, k10d, F2T2, and F2V2.

Column 2

The chamber volumes, $V$. Chamber volumes must be identical within data series. In the sample data above, $V=140.625$ [$L$] for data series k0d and k10d, and $V=2.0101$ [$L$] for data series F2T2 and F2V2.

Column 3

The chamber cross-sectional areas, $A$. Chamber cross-sectional areas must be identical within data series. In the sample data above, $A=0.5625$ [$m^2$] for data series k0d and k10d, and $A=0.0201$ [$m^2$] for data series F2T2 and F2V2.

Column 4

The measurement time-points within data series in strictly increasing order and with the first time-point equal to zero. At least three observation time-points is required per data series. In the sample data above, time-points are in minutes and cover, approximately, two-four hour periods per data series.

Column 5

The measured chamber concentrations corresponding to the time-points in the fourth column. Chamber concentations must be strictly positive. In the sample data above, the fifth column contains measured nitrous oxide concentrations [$\mu g/L$].

Missing values (NA's) or empty lines or columns are not allowed in HMR data files.

HMR data and flux physical units

For maximal flexibility, HMR has no requirements for the physical units of input data. The chosen units do, however, determine the unit of the estimated flux, which has the physical unit of $(VC)/ (At)$, where $t$ and $C$ denote, respectively, time and concentration. Some examples:

$V$ $[L]$, $A$ $[m^2]$, $t$ $[h]$, $C$ $[\mu g/L]$	$\Rightarrow$	$f_0$ $[\mu g/m^2/h]$
$V$ $[L]$, $A$ $[m^2]$, $t$ $[min]$, $C$ $[\mu L/L]$	$\Rightarrow$	$f_0$ $[\mu L/m^2/min]$
$V$ $[m^2]$, $A$ $[km^2]$, $t$ $[s]$, $C$ $[kg/m^3]$	$\Rightarrow$	$f_0$ $[kg/km^2/s]$

Value

HMR first examines the specified function arguments. If problems are detected, HMR halts and issues the comment 'Error in input parameters'. Rules for the HMR function arguments are provided above in the Arguments section. If specified function arguments are valid, HMR continues to examine the provided data file. If fatal problems are detected, HMR halts and issues the comment 'Data file could not be read'. Typical problems with data files are:

(i)	Inconsistency between the decimal and column separators used in the data file and the standards for these separators on the computer. The easiest way to resolve this issue is to specify the separators used in the data file as function arguments to HMR, ie. dec and sep, respectively.
(ii)	Empty rows below or columns beside valid HMR data are not allowed in HMR data files.

If specified function arguments are valid, and the data file contains at least one valid data series, HMR continues to analyze data. The output is then a data frame with one row per analysed data series and variables:

Series	Name of the data series.
f0	The estimated flux.
f0.se	The standard error of the estimated flux.
f0.p	The p-value for the null hypothesis of zero flux.
f0.lo95	The lower end-point of the 95%-confidence interval for the flux.
f0.up95	The upper end-point of the 95%-confidence interval for the flux.
Method	The method used for estimating the flux ('HMR', 'LR', 'No flux' or 'None').
Warning	A character string with a warning message in case of estimation problems.
Prefilter	The prefiltering classification ('Signal', 'Noise' or 'None').
Prefilter.p	The prefiltering p-value.
SatCrit.Warning	A character string with a warning message if flux limiting is active.
LR.f0	The flux estimated by linear regression. (Only present if LR.always=TRUE.)
LR.f0.se	The standard error of the flux estimated by linear regression. (Only present if LR.always=TRUE.)
LR.f0.p	The p-value for the null hypothesis of zero flux calculated by linear regression. (Only present if LR.always=TRUE.)
LR.f0.lo95	The lower end-point of the 95%-confidence interval for the flux calculated by linear regression. (Only present if LR.always=TRUE.)
LR.f0.up95	The upper end-point of the 95%-confidence interval for the flux calculated by linear regression. (Only present if LR.always=TRUE.)
LR.Warning	A character string with a warning message if linear regression estimated a negative predeployment concentration. (Only present if LR.always=TRUE.)

The data frame is also exported to a semicolon/comma separated file with the name of the data file preceded by 'HMR - ' and located in the data file folder. The exported file uses the same column and decimal separators as the data file.

Author(s)

Asger R. Pedersen, Ph.D. in statistics, SEGES Innovation, Aarhus, Denmark

References

Hutchinson, G.L. and Mosier, A.R. (1981). Improved soil cover method for field measurement of nitrous oxide fluxes. Soil Science Society of America Journal, 45, pp. 311-316

Seber, G.A.F. and Wild, C.J. (1989). Nonlinear regression. Wiley, New York

Pedersen, A.R., Petersen, S.O. and Schelde, K. (2010): A comprehensive approach to soil-atmosphere trace-gas flux estimation with static chambers. European Journal of Soil Science, 61, pp. 888-902

Pullens, J.W.M., Abalos, D., Petersen, S.O. and Pedersen, A.R. (2023). Identifying criteria for greenhouse gas flux estimation with automatic and manual chambers: A case study for N2O. European Journal of Soil Science, 74, e13340. https://doi.org/10.1111/ejss.13340

Examples

## Not run: 
# Suppose the sample data above are located on a Windows machine in the file
# 'C:\My HMR applications\N2O.csv'

# Start by setting the data file folder:
setwd('C:/My HMR applications')

# Notice that R uses '/' in folder declarations, whereas Windows uses '\'.

# Analyse all data series with default settings:
HMR(filename='N2O.csv')

# Produces the following (slightly edited) output when HMR recommendations are followed:
  Series        f0     f0.se      f0.p    f0.lo95   f0.up95
1    k0d 1.872e+01 4.188e+00 2.084e-03  9.062e+00 2.838e+01
2   k10d 1.797e-01 1.385e-01 2.269e-01 -1.337e-01 4.930e-01
3   F2T2 4.509e+01 1.470e+01 5.469e-02 -1.697e+00 9.188e+01
4   F2V2 5.584e+01 3.586e+00 9.926e-05  4.588e+01 6.579e+01
       Method Warning Prefilter Prefilter.p SatCrit.Warning
     1    HMR    None      None          NA              NA
     2    HMR    None      None          NA              NA
     3    HMR    None      None          NA              NA
     4     LR    None      None          NA              NA
# The non-linear 'HMR' analysis was recommended by the MSE criterion for three data series,
# whereas the 'LR' analysis was recommended for the fourth.

# The output was also exported to the semicolon-separated file:
# 'C:\My HMR applications\HMR - N2O.csv'

# Analyse all data series with flux limiting assuming 90% chamber saturation not before 60 minutes:
HMR(filename='N2O.csv',SatPct=90,SatTimeMin=60)

# Produces the following (slightly edited) output when HMR recommendations are followed:
  Series        f0     f0.se      f0.p    f0.lo95   f0.up95
1    k0d 1.872e+01 4.188e+00 2.084e-03  9.061e+00 2.838e+01
2   k10d 2.689e-02 9.257e-03 1.571e-02  6.261e-03 4.751e-02
3   F2T2 4.509e+01 1.470e+01 5.469e-02 -1.697e+00 9.188e+01
4   F2V2 5.584e+01 3.586e+00 9.926e-05  4.588e+01 6.579e+01
       Method Warning Prefilter Prefilter.p
     1    HMR    None      None          NA
     2     LR    None      None          NA
     3    HMR    None      None          NA
     4     LR    None      None          NA
                                                SatCrit.Warning
          1                                                None
          2 Flux limited by saturation assumption - consider LR
          3                                                None
          4                                                None
# The chamber saturation assumption excluded the MSE optimal value of 'kappa', and HMR
# therefore recommends 'LR'.

# The output was in both HMR analyses above exported to the semicolon-separated file
# named 'C:\My HMR applications\HMR - N2O.csv'. Hence, several analyses of the same
# data file overwrites the output file, so to save a particular output file it has to
# be renamed before the next HMR analysis of the same data file.

## End(Not run)

---