View source: R/fluxEstimates.R
calcClosedChamberFlux | R Documentation |
Calculate gas flux and its uncertainties for a non steady-state canopy chamber.
calcClosedChamberFlux(ds, colConc = "CO2_dry",
colTime = "TIMESTAMP", colTemp = "TA_Avg",
colPressure = "Pa", volume = 1, area = 1,
fRegress = c(lin = regressFluxLinear,
tanh = regressFluxTanh), fRegressSelect = regressSelectPref1,
concSensitivity = 1, maxLag = 50, minTLag = 0,
useFixedTLag = NA, debugInfo = list(useOscarsLagDectect = FALSE,
omitAutoCorrFit = FALSE, omitEstimateLeverage = FALSE,
tLagFixed = NA, isStopOnError = FALSE),
...)
ds |
data.frame with concentration and time column of a chamber measurement of one replicate |
colConc |
column name of CO2 concentration [ppm] |
colTime |
column name of time [s], must be of class POSIXct or numeric or integer |
colTemp |
column name of air temperature inside chamber [degC] |
colPressure |
column name of air pressure inside chamber [Pa] |
volume |
volume inside the chamber im [m3] |
area |
area of the exchange surface [m2] |
fRegress |
list of functions to yield a single flux estimate, see details |
fRegressSelect |
function to select the regression
function based on fitting results. Signature and return must correspond
to |
concSensitivity |
measurement sensitivity of concentration. With concentration change below this sensitivity, only a linear model is fit |
maxLag |
number of initial records to be screened for a breakpoint, i.e. the lag (higher for water vapour than for CO2) |
minTLag |
possibility to specify a minimum lag-time in seconds |
useFixedTLag |
scalar numeric value in seconds, if not-NA, use the specified lag-time instead of estimating the lag-time from the concentration data |
debugInfo |
rarely used controls, mostly for debugging
|
... |
further arguments to |
The function fRegress
must conform to regressFluxSquare
, i.e.
return a vector of length 2: the flux estimate and its standard deviation.
Optionally, it may return the model fit object in attribute "model"
If several functions are given, then the best fit is selected according
to function with argument fRegressSelect
, by default to the AIC criterion.
Fit an expoenential curve by using function regressFluxExp
.
There are two kinds of uncertainty associated with the flux. The first comes from the uncertainty of the slope of concentration increase. The second comes from the leverage of starting and end points of the regression (estimated by a bootstrap) return value sdFlux is the maximum of those two components
data.frame (tibble) of one row with entries
flux |
the estimate of the CO2 flux into the chamber [mumol / m2 / s] |
fluxMedian |
the median of the flux bootsrap estimates [mumol / m2 / s] |
tLag |
time of lag phase in seconds |
lagIndex |
index of the row at the end of lag-time |
autoCorr |
autocorrelation coefficient, NA if model with autocorrelation could not be fitted or had higher AIC than a model without autocorrelation |
AIC |
AIC goodness of fit for the model |
sdFluxRegression |
the standard deviation of the flux by a single regression of CO2 flux |
sdFluxLeverage |
the standard deviation of the flux by leverage of starting or end values of the time series |
iFRegress |
index of the best (lowest AIC) regression function |
times |
integer vector of predictor times for the model fit (excluding and possibly first record) |
Thomas Wutzler, Oscar Perez Priego
RespChamberProc
#data(chamberLoggerEx1s)
library(dplyr)
ds <- chamberLoggerEx1s
ds$Pa <- chamberLoggerEx1s$Pa * 1000 # convert kPa to Pa
conc <- ds$CO2_dry <- corrConcDilution(ds)
#trace(calcClosedChamberFlux, recover) #untrace(calcClosedChamberFlux)
resFit <- calcClosedChamberFlux(ds)
select(resFit, flux, sdFlux) %>% unlist()
#plotResp(ds, resFit)
.tmp.compareFittingFunctions <- function(){
resLin <- calcClosedChamberFlux(ds, fRegress = list(regressFluxLinear))
resPoly <- calcClosedChamberFlux(ds, fRegress = list(regressFluxSquare))
#trace(regressFluxExp, recover) #untrace(regressFluxExp)
resExp <- calcClosedChamberFlux(ds, fRegress = list(regressFluxExp) )
#trace(regressFluxTanh, recover) #untrace(regressFluxTanh)
resTanh <- calcClosedChamberFlux(ds, fRegress = list(regressFluxTanh ))
times <- ds$TIMESTAMP
times0 <- as.numeric(times) - as.numeric(times[1])
times0Fit <- times0[times0 >= resLin$stat["tLag"] ]
plot( resid(resTanh$model, type = "normalized") ~ times0Fit ) # residual plots
qqnorm(resid(resTanh$model, type = "normalized")); abline(0,1)
#length(times0Fit)
plot( ds$CO2_dry ~ times0, xlab = "time (s)", ylab = "" )
mtext("CO2_dry (ppm)",2,2,las = 0)
abline(v = resLin$stat["tLag"], col = "grey", lty = "dotted")
lines( fitted(resExp$model) ~ times0Fit , col = "red" )
lines( fitted(resLin$model) ~ times0Fit , col = "grey" )
lines( fitted(resTanh$model) ~ times0Fit , col = "purple" )
lines( fitted(resPoly$model) ~ times0Fit , col = "blue" )
legend("topright", inset = c(0.02,0.02), legend = c("exp","lin","tanh","poly")
, col = c("red","grey","purple","blue"), lty = "solid")
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.