Perform parametric curve fitting or non-parametric estimation of mNIRS kinetics on vector data.
Usage
# S3 method for class 'monoexponential'
process_kinetics(
y,
x = NULL,
data = NULL,
method = c("monoexponential", "sigmoidal", "half_time", "peak_slope"),
end_fit_span = 20,
x0 = 0,
verbose = TRUE,
...
)
# S3 method for class 'sigmoidal'
process_kinetics(
y,
x = NULL,
data = NULL,
method = c("monoexponential", "sigmoidal", "half_time", "peak_slope"),
end_fit_span = 20,
x0 = 0,
verbose = TRUE,
...
)
# S3 method for class 'half_time'
process_kinetics(
y,
x = NULL,
data = NULL,
method = c("monoexponential", "sigmoidal", "half_time", "peak_slope"),
end_fit_span = 20,
x0 = 0,
verbose = TRUE,
...
)
# S3 method for class 'peak_slope'
process_kinetics(
y,
x = NULL,
data = NULL,
method = c("monoexponential", "sigmoidal", "half_time", "peak_slope"),
end_fit_span = 20,
x0 = 0,
verbose = TRUE,
...
)
Arguments
- y
A numeric vector of the response variable, or the name of the variable within a dataframe.
- x
An optional numeric vector of the predictor variable, or the name of the variable within a dataframe. Defaults to using the index of
x = seq_along(y)
.- data
An optional dataframe containing the predictor and response variables named in
x
andy
. Names forx
andy
must be in quotations.- method
Indicates which model to evaluate the kinetics event (see Details for method parametrisation).
method = "monoexponential"
A four-parameter monoexponential association function in the form:
ifelse(x <= TD, A, A + (B - A) * (1 - exp((TD - x) / tau)))
.method = "sigmoidal"
A four-parameter generalised logistic (sigmoidal) function in the form:
A + (B - A) / (1 + exp((xmid - x) / scal))
.method = "half_time"
A non-parametric estimate of the time to recover half of the total kinetics amplitude.
method = "peak_slope"
A non-parametric estimate of the peak dynamic rate of change in the mNIRS signal, calculated from rolling local linear regression slope within a window in units of
x
defined bypeak_slope_span
.
- end_fit_span
(Default
end_fit_span = 20
) A numeric value indicating the local window in units of the predictor variablex
after the kinetics extreme (peak or trough) value to look for subsequent greater extremes. The kinetics model will be fit to the data up to the first local extreme with no subsequent greater extremes within a trailing window in units ofx
defined byend_fit_span
, or the limits of the data.- x0
(Default = 0) A numeric value indicating the value of the predictor variable
x
representing the start of the kinetics event.- verbose
A logical.
TRUE
(the default) will return warnings and messages which can be used for data error checking.FALSE
will silence these messages. Errors will always be returned.- ...
Additional arguments.
peak_slope_span
A numeric value defining the window in units of the predictor variable
x
for rolling local calculations used formethod =
"peak_slope".align = c(
"center", "left", "right"
)
Specifies the window alignment of
peak_slope_span
as "center" (the default), "left", or "right". Where "left" is forward looking, and "right" is backward looking by the windowpeak_slope_span
from the current time. Only used formethod =
"peak_slope".fixed parameters
Parameters (coefficients) of the parametric models ("monoexponential" and "sigmoidal") can be defined a priori and fixed, to exclude them from the model fitting optimisation. e.g.,
A = 10
will define the functionSSmonoexp(x, A = 10, B, TD, tau)
.
Value
A list of class mnirs.kinetics
with components:
method
The kinetics method used.
model
The model object.
equation
The equation of the kinetics model used.
data
A dataframe of original and fitted model data.
fitted
A vector of fitted values returned by the model.
residuals
A vector of residuals between original and fitted values returned by the model.
x0
The value of the predictor variable indicating the start of kinetics.
coefs
A dataframe of model coefficients, including manually fixed parameters.
diagnostics
A dataframe of model goodness-of-fit metrics (
AIC
,BIC
,R2adj
RMSE
,MAE
,MAPE
).call
The model call.
Details
method %in% c("monoexponential", "sigmoidal")
use nls()
for nonlinear (weighted) least-squares estimates.
Model Parameterisation
"monoexponential": A four-parameter monoexponential association function
in the form: ifelse(x <= TD, A, A + (B - A) * (1 - exp((TD - x) / tau)))
.
A
: textB
: textTD
: texttau
: textMRT
: text
"sigmoidal": A four-parameter generalised logistic (sigmoidal) function
in the form: A + (B - A) / (1 + exp((xmid - x) / scal))
.
A
: textB
: textxmid
: textscal
: text
"half_time": A non-parametric estimate of the time to recover half of the total kinetics amplitude.
A
: textB
: texthalf_value
: text
"peak_slope": A non-parametric estimate of the peak dynamic rate of change in
the mNIRS signal, calculated from rolling local linear regression slope
within a window in units of x
defined by peak_slope_span
.
A
: textB
: textx_fitted
: textpeak_slope
: text
Examples
set.seed(13)
x1 <- seq(-10, 60, by = 2)
A <- 10; B <- 100; TD <- 5; tau <- 12
y1 <- monoexponential(x1, A, B, TD, tau) + rnorm(length(x1), 0, 3)
## monoexponential kinetics ===============================
model <- process_kinetics(y1, x1, method = "monoexponential")
model
#>
#> Monoexponential Association Nonlinear nls() Regression
#> model: y ~ SSmonoexp(x, A, B, TD, tau)
#> data: data.frame(x = x, y = y)
#> equation: y ~ A + (B - A) * (1 - exp((TD - x) / tau))
#>
#> Model Coefficients:
#> A B TD tau MRT MRT_fitted
#> 11.97 98.9 5.436 11.65 17.08 66.92
#>
#> Model Diagnostics:
#> AIC BIC RMSE MAE MAPE
#> 181 188.9 2.6 2.046 5.465
#>
if (FALSE) { # \dontrun{
## require(ggplot2)
plot(model, plot_coefs = TRUE, plot_diagnostics = TRUE, plot_residuals = TRUE)
} # }
## sigmoidal kinetics ===============================
model <- process_kinetics(y1, x1, method = "sigmoidal")
model
#>
#> Sigmoidal Generalised Logistic Nonlinear nls() Regression
#> model: y ~ SSfpl(x, A, B, xmid, scal)
#> data: data.frame(x = x, y = y)
#> equation: y ~ A + (B - A) / (1 + exp((xmid - x) / scal))
#>
#> Model Coefficients:
#> A B scal xmid xmid_fitted
#> 7.929 95.47 5.178 13.92 51.7
#>
#> Model Diagnostics:
#> AIC BIC RMSE MAE MAPE
#> 199.9 207.8 3.381 2.661 9.035
#>
if (FALSE) { # \dontrun{
## require(ggplot2)
plot(model, plot_coefs = TRUE, plot_diagnostics = TRUE, plot_residuals = TRUE)
} # }
## half recovery time ===============================
model <- process_kinetics(y1, x1, method = "half_time")
model
#>
#> Non-parametric Half Response Time
#> data: data.frame(x = x, y = y)
#> equation: half_x ~ interp_x(y = A + (B - A) / 2)
#>
#> Model Coefficients:
#> A B half_x half_value
#> 11.9 100.8 16 56.32
#>
if (FALSE) { # \dontrun{
## require(ggplot2)
plot(model, plot_coefs = TRUE, plot_diagnostics = TRUE, plot_residuals = TRUE)
} # }
## peak slope ===============================
model <- process_kinetics(y1, x1, method = "peak_slope", peak_slope_span = 10)
model
#>
#> Peak Linear Regression Rolling Slope
#> data: data.frame(x = x, y = y)
#> equation: peak_slope ~ max(rolling_slope(y ~ x, peak_slope_span = 10))
#>
#> Model Coefficients:
#> peak_slope_x peak_slope_fitted peak_slope
#> 8 33.22 5.11
#>
if (FALSE) { # \dontrun{
## require(ggplot2)
plot(model, plot_coefs = TRUE, plot_diagnostics = TRUE, plot_residuals = TRUE)
} # }