Title: | Tools for Nonlinear Regression Analysis |
---|---|
Description: | Several tools for assessing the quality of fit of a gaussian nonlinear model are provided. |
Authors: | Florent Baty [aut], Marie-Laure Delignette-Muller [aut], Sandrine Charles [ctb], Jean-Pierre Flandrois [ctb], Christian Ritz [ctb], Aurelie Siberchicot [aut, cre] |
Maintainer: | Aurelie Siberchicot <[email protected]> |
License: | GPL-3 |
Version: | 2.1-0 |
Built: | 2024-11-16 05:20:17 UTC |
Source: | https://github.com/lbbe-software/nlstools |
Produces confidence intervals for the parameters in nonlinear regression model fit. The intervals can either be based large sample results or on profiling.
confint2(object, parm, level = 0.95, method = c("asymptotic", "profile"), ...)
confint2(object, parm, level = 0.95, method = c("asymptotic", "profile"), ...)
object |
object of class |
parm |
a vector character strings with names of the parameter for which to calculate confidence intervals (by default all parameters). |
level |
the confidence level required. |
method |
method to be used: "asympotic" for large sample and "profile" for profiling approach. |
... |
additional argument(s) to pass on the method doing the profiling. |
The profiling used is the method confint.nls.
A matrix with columns giving lower and upper confidence limits for each parameter.
Christian Ritz
L.minor.m1 <- nls(rate ~ Vm*conc/(K+conc), data = L.minor, start = list(K=20, Vm=120)) confint2(L.minor.m1) confint2(L.minor.m1, "K")
L.minor.m1 <- nls(rate ~ Vm*conc/(K+conc), data = L.minor, start = list(K=20, Vm=120)) confint2(L.minor.m1) confint2(L.minor.m1, "K")
Enzyme kinetics
data(L.minor)
data(L.minor)
A data frame with 8 observations on the following 2 variables.
conc
a numeric vector
rate
a numeric vector
Cedergreen, N. and Madsen, T. V. (2002) Nitrogen uptake by the floating macrophyte Lemna minor, New Phytologist, 155, 285–292.
Michaelis Menten data sets
data(vmkm) data(vmkmki)
data(vmkm) data(vmkmki)
vmkm
is a data frame with 2 columns (S: concentration of substrat, v: reaction rate)vmkmki
is a data frame with 3 columns (S: concentration of substrat, I: concentration of inhibitor, v: reaction rate)
These datasets were provided by the French research unit INRA UMR1233.
data(vmkm) data(vmkmki) plot(vmkm) plot(vmkmki)
data(vmkm) data(vmkmki) plot(vmkm) plot(vmkmki)
Formula of Michaelis-Menten model commonly used to describe enzyme kinetics, and derived formulas taking into account the effect of a competitive or a non-competitive inhibitor
michaelis compet_mich non_compet_mich
michaelis compet_mich non_compet_mich
These models describe the evolution of the reaction rate (v) as a function
of the concentration of substrate (S) and the concentration of inhibitor (I) for compet_mich
and non_compet_mich
.
michaelis
is the classical Michaelis-Menten model (Dixon, 1979) with two parameters (Km, Vmax) :
compet_mich
is the Michaelis-Menten derived model with three parameters (Km, Vmax, Ki), describing
a competitive inhibition :
non_compet_mich
is the Michaelis-Menten derived model with three parameters (Km, Vmax, Ki), describing
a non-competitive inhibition :
A formula
Florent Baty, Marie-Laure Delignette-Muller
Dixon M and Webb EC (1979) Enzymes, Academic Press, New York.
# Example 1 data(vmkm) nls1 <- nls(michaelis,vmkm,list(Km=1,Vmax=1)) plotfit(nls1, smooth = TRUE) # Example 2 data(vmkmki) def.par <- par(no.readonly = TRUE) par(mfrow = c(2,2)) nls2_c <- nls(compet_mich, vmkmki, list(Km=1,Vmax=20,Ki=0.5)) plotfit(nls2_c, variable=1) overview(nls2_c) res2_c <- nlsResiduals(nls2_c) plot(res2_c, which=1) nls2_nc <- nls(non_compet_mich, vmkmki, list(Km=1, Vmax=20, Ki=0.5)) plotfit(nls2_nc, variable=1) overview(nls2_nc) res2_nc <- nlsResiduals(nls2_nc) plot(res2_nc, which=1) par(def.par)
# Example 1 data(vmkm) nls1 <- nls(michaelis,vmkm,list(Km=1,Vmax=1)) plotfit(nls1, smooth = TRUE) # Example 2 data(vmkmki) def.par <- par(no.readonly = TRUE) par(mfrow = c(2,2)) nls2_c <- nls(compet_mich, vmkmki, list(Km=1,Vmax=20,Ki=0.5)) plotfit(nls2_c, variable=1) overview(nls2_c) res2_c <- nlsResiduals(nls2_c) plot(res2_c, which=1) nls2_nc <- nls(non_compet_mich, vmkmki, list(Km=1, Vmax=20, Ki=0.5)) plotfit(nls2_nc, variable=1) overview(nls2_nc) res2_nc <- nlsResiduals(nls2_nc) plot(res2_nc, which=1) par(def.par)
Bootstrap resampling
nlsBoot (nls, niter = 999) ## S3 method for class 'nlsBoot' plot(x, type = c("pairs", "boxplot"), mfr = c(ceiling(sqrt(ncol(x$coefboot))), ceiling(sqrt(ncol(x$coefboot)))), ask = FALSE, ...) ## S3 method for class 'nlsBoot' print(x, ...) ## S3 method for class 'nlsBoot' summary(object, ...)
nlsBoot (nls, niter = 999) ## S3 method for class 'nlsBoot' plot(x, type = c("pairs", "boxplot"), mfr = c(ceiling(sqrt(ncol(x$coefboot))), ceiling(sqrt(ncol(x$coefboot)))), ask = FALSE, ...) ## S3 method for class 'nlsBoot' print(x, ...) ## S3 method for class 'nlsBoot' summary(object, ...)
nls |
an object of class 'nls' |
niter |
number of iterations |
x , object
|
an object of class 'nlsBoot' |
type |
type of representation (options are "pairs" or "boxplot") |
mfr |
layout definition (number of rows and columns in the graphics device) |
ask |
if TRUE, draw plot interactively |
... |
further arguments passed to or from other methods |
Non-parametric bootstrapping is used. Mean centered residuals are bootstrapped. By default, 999 resampled data sets are created from which parameter estimates are obtained by fitting the model on each of these data sets. Whenever the fit fails to converge, a flag reports the number of non-convergences. If the fitting procedure fails to converge in more than 50% of the cases, the procedure is interrupted with a flag and no result is given. The function summary
returns the bootstrap estimates (mean and std. dev. of the bootstrapped estimates) and the median and 95 percent confidence intervals (50, 2.5, and 97.5 percentiles of the bootstrapped estimates). The bootstrapped estimate distributions can be visualized using the function plot.nlsBoot
either by plotting the bootstrapped sample for each pair of parameters or by displaying the boxplot representation of the bootstrapped sample for each parameter. Notice that nlsBoot
does not currently handle transformed dependent variables specified in the left side of the nls
formula.
nlsBoot
returns a list of 5 objects:
coefboot |
contains the bootstrap parameter estimates |
bootCI |
contains the bootstrap medians and the bootstrap 95% confidence intervals |
estiboot |
contains the means and std. errors of the bootstrap parameter estimates |
rse |
is the vector of bootstrap residual errors |
nls |
the object of class 'nls' given in input |
Florent Baty, Marie-Laure Delignette-Muller
Bates DM and Watts DG (1988) Nonlinear regression analysis and its applications. Wiley, Chichester, UK.
Huet S, Bouvier A, Poursat M-A, Jolivet E (2003) Statistical tools for nonlinear regression: a practical guide with S-PLUS and R examples. Springer, Berlin, Heidelberg, New York.
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2res + (t > 5.883) * (VO2res + (VO2peak - VO2res) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2res = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.boot1 <- nlsBoot(O2K.nls1, niter = 200) plot(O2K.boot1) plot(O2K.boot1, type = "boxplot", ask = FALSE) summary(O2K.boot1)
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2res + (t > 5.883) * (VO2res + (VO2peak - VO2res) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2res = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.boot1 <- nlsBoot(O2K.nls1, niter = 200) plot(O2K.boot1) plot(O2K.boot1, type = "boxplot", ask = FALSE) summary(O2K.boot1)
Computation of confidence intervals on predictions from Bootstrap resampling
nlsBootPredict(nlsBoot, newdata, interval = c("confidence", "prediction"))
nlsBootPredict(nlsBoot, newdata, interval = c("confidence", "prediction"))
nlsBoot |
An object of class 'nlsBoot'. |
newdata |
A data frame in which to look for values of independent variables for the predictions.If omitted, the data used for fitting are used. |
interval |
Type of interval to compute, "confidence", or "prediction". |
nlsBootPredict
produces confidence intervals on predicted values
that can be obtained using predict.nls
for values of the
independent variable(s) defined in the data frame newdata.
Non-parametric bootstrapping is used (results of nlsBoot
).
For confidence intervals the bootstrap sample of predictions is simply computed
from the bootstrap sample of estimations of the model parameters, by evaluating
the mean value of the model on each new data. For prediction intervals,
to take into account the residual errors, a residual error sampled
in the mean centered residuals is added to each mean predicted
value. In both cases, bootstrap predictions are summarized by the median and 95
percent confidence intervals (50, 2.5, and 97.5 percentiles of the bootstrapped values).
nlsBoot
returns a matrix of predictions with three columns
respectively corresponding to the 50, 2.5 and 97.5 percentiles of bootstrap predictions.
Florent Baty, Marie-Laure Delignette-Muller
Huet S, Bouvier A, Poursat M-A, Jolivet E (2003) Statistical tools for nonlinear regression: a practical guide with S-PLUS and R examples. Springer, Berlin, Heidelberg, New York.
See nlsBoot
and predict.nls
.
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2res + (t > 5.883) * (VO2res + (VO2peak - VO2res) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2res = 400, VO2peak = 1600, mu = 1), data = O2K) niter <- 200 ### To reach stable prediction intervals use far greater niter (>> 1000) O2K.boot1 <- nlsBoot(O2K.nls1, niter = niter) newdata <- data.frame(t = seq(0, 12, length.out = 50)) (pred.clim <- nlsBootPredict(O2K.boot1, newdata = newdata, interval = "confidence")) (pred.plim <- nlsBootPredict(O2K.boot1, newdata = newdata, interval = "prediction")) plotfit(O2K.nls1, smooth = TRUE, ylim = c(200, 1800)) lines(newdata$t, pred.clim[, 2], col = "red") lines(newdata$t, pred.clim[, 3], col = "red") lines(newdata$t, pred.plim[, 2], col = "blue") lines(newdata$t, pred.plim[, 3], col = "blue") ### An example without giving newdata # plot of data plot(O2K$t, O2K$VO2) # add of predictions computed using predict.nls() pred <- predict(O2K.nls1) points(O2K$t, pred, pch = 16) # add of prediction intervals using nlsBootPredict() (pred.plim <- nlsBootPredict(O2K.boot1, interval = "prediction")) segments(O2K$t, pred.plim[, 2], O2K$t, pred.plim[, 3], col = "blue")
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2res + (t > 5.883) * (VO2res + (VO2peak - VO2res) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2res = 400, VO2peak = 1600, mu = 1), data = O2K) niter <- 200 ### To reach stable prediction intervals use far greater niter (>> 1000) O2K.boot1 <- nlsBoot(O2K.nls1, niter = niter) newdata <- data.frame(t = seq(0, 12, length.out = 50)) (pred.clim <- nlsBootPredict(O2K.boot1, newdata = newdata, interval = "confidence")) (pred.plim <- nlsBootPredict(O2K.boot1, newdata = newdata, interval = "prediction")) plotfit(O2K.nls1, smooth = TRUE, ylim = c(200, 1800)) lines(newdata$t, pred.clim[, 2], col = "red") lines(newdata$t, pred.clim[, 3], col = "red") lines(newdata$t, pred.plim[, 2], col = "blue") lines(newdata$t, pred.plim[, 3], col = "blue") ### An example without giving newdata # plot of data plot(O2K$t, O2K$VO2) # add of predictions computed using predict.nls() pred <- predict(O2K.nls1) points(O2K$t, pred, pch = 16) # add of prediction intervals using nlsBootPredict() (pred.plim <- nlsBootPredict(O2K.boot1, interval = "prediction")) segments(O2K$t, pred.plim[, 2], O2K$t, pred.plim[, 3], col = "blue")
Draws parameter values in the Beale's 95 percent unlinearized confidence region
nlsConfRegions (nls, length = 1000, exp = 1.5) ## S3 method for class 'nlsConfRegions' plot(x, bounds = FALSE, ask = FALSE, ...) ## S3 method for class 'nlsConfRegions' print(x, ...)
nlsConfRegions (nls, length = 1000, exp = 1.5) ## S3 method for class 'nlsConfRegions' plot(x, bounds = FALSE, ask = FALSE, ...) ## S3 method for class 'nlsConfRegions' print(x, ...)
nls |
an object of class 'nls' |
length |
number of points to draw in the confidence region |
exp |
expansion factor of the hypercube in which random values of parameters are drawn |
x |
an object of class 'nlsConfRegions' |
bounds |
logical defining whether bounds of the drawing hypercube are plotted |
ask |
if TRUE, draw plot interactively |
... |
further arguments passed to or from other methods |
A sample of points in the 95 percent confidence region is computed according to Beale's criterion (Beale, 1960). This region is also named the joint parameter likelihood region (Bates and Watts, 1988). The method used consists in a random sampling of parameters values in a hypercube centered on the least squares estimate and rejecting the parameters values whose residual sum of squares do not verify the Beale criterion. The confidence region is plotted by projection of the sampled points in each plane defined by a couple of parameters. Bounds of the hypercube in which random values of parameters are drawn may be plotted in order to check if the confidence region was totally included in the hypercube defined by default. If not the hypercube should be expanded in order to obtain the full confidence region
nlsConfRegions
returns a list of four objects:
cr |
a data frame containing the sample drawn in the Beale's confidence region |
rss |
a vector containing the residual sums of squares corresponding to |
rss95 |
the 95 percent residual sum of squares threshold according to Beale (1960) |
bounds |
lower and upper bounds of the hypercube in which random values of parameters have been drawn |
Florent Baty, Marie-Laure Delignette-Muller
Beale EML (1960) Confidence regions in non-linear estimations. Journal of the Royal Statistical Society, 22B, 41-88.
Bates DM and Watts DG (1988) Nonlinear regression analysis and its applications. Wiley, Chichester, UK.
ellipse.nls
in the ellipse
library
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.conf1 <- nlsConfRegions(O2K.nls1, exp = 2, length = 200) plot(O2K.conf1, bounds = TRUE)
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.conf1 <- nlsConfRegions(O2K.nls1, exp = 2, length = 200) plot(O2K.conf1, bounds = TRUE)
Provides residual sum of squares (RSS) contours
nlsContourRSS (nls, lseq = 100, exp = 2) ## S3 method for class 'nlsContourRSS' plot(x, nlev = 0, col = TRUE, col.pal = terrain.colors(100), ask = FALSE, useRaster = TRUE, ...) ## S3 method for class 'nlsContourRSS' print(x, ...)
nlsContourRSS (nls, lseq = 100, exp = 2) ## S3 method for class 'nlsContourRSS' plot(x, nlev = 0, col = TRUE, col.pal = terrain.colors(100), ask = FALSE, useRaster = TRUE, ...) ## S3 method for class 'nlsContourRSS' print(x, ...)
nls |
an object of class 'nls' |
lseq |
length of the sequences of parameters |
exp |
expansion factor of the parameter intervals defining the grids |
nlev |
number of contour levels to add to the likelihood contour at level 95 percent |
col |
logical. Contours are plotted with colors if |
col.pal |
Palette of colors. Colors to be used as background (default is terrain.colors(100); unused if col is FALSE) |
x |
an object of class 'nlsContourRSS' |
ask |
if TRUE, draw plot interactively (default is FALSE) |
useRaster |
a bitmap raster is used to plot the image instead of polygons (default is TRUE) |
... |
further arguments passed to or from other methods |
The aim of these functions is to plot the residual sum of squares (RSS) contours which correspond to likelihood contours for a Gaussian model. For each pair of parameters the RSS is calculated on a grid centered on the least squares estimates of both parameters, the other parameters being fixed to their least square estimates. The contours of RSS values are then plotted for each pair of parameters. For each pair of parameters, one of this contour corresponds to a section of the 95 percent Beale's confidence region in the plane of these parameters. This contour is plotted in a different color.
nlsContourRSS
returns a list of three objects:
seqPara |
a matrix with the sequence of grid values for each parameter |
lrss |
a list of matrices with logarithm values of RSS in the grid for each pair of parameters |
lrss95 |
the logarithm of the 95 percent residual sum of squares threshold according to Beale (1960) |
Florent Baty, Marie-Laure Delignette-Muller
Beale EML (1960) Confidence regions in non-linear estimations. Journal of the Royal Statistical Society, 22B, 41-88.
Bates DM and Watts DG (1988) Nonlinear regression analysis and its applications. Wiley, Chichester, UK.
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.cont1 <- nlsContourRSS(O2K.nls1) plot(O2K.cont1)
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.cont1 <- nlsContourRSS(O2K.nls1) plot(O2K.cont1)
Jackknife resampling
nlsJack (nls) ## S3 method for class 'nlsJack' plot(x, mfr = c(nrow(x$reldif),1), ask = FALSE, ...) ## S3 method for class 'nlsJack' print(x, ...) ## S3 method for class 'nlsJack' summary(object, ...)
nlsJack (nls) ## S3 method for class 'nlsJack' plot(x, mfr = c(nrow(x$reldif),1), ask = FALSE, ...) ## S3 method for class 'nlsJack' print(x, ...) ## S3 method for class 'nlsJack' summary(object, ...)
nls |
an object of class 'nls' |
x , object
|
an object of class 'nlsJack' |
mfr |
layout definition, default is k rows (k: number of parameters) and 1 column |
ask |
if TRUE, draw plot interactively |
... |
further arguments passed to or from other methods |
A jackknife resampling procedure is performed. Each observation is sequentially removed from the initial data set using a leave-one-out strategy. A data set with n observations provides thus n resampled data sets of n-1 observations. The jackknife estimates with confidence intervals are calculated as described by Seber and Wild (1989) from the results of n new fits of the model on the n jackknife resampled data sets. The leave-one-out procedure is also employed to assess the influence of each observation on each parameter estimate. An observation is empirically defined as influential for one parameter if the difference between the estimate of this parameter with and without the observation exceeds twice the standard error of the estimate divided by sqrt(n). This empirical method assumes a small curvature of the nonlinear model. For each parameter, the absolute relative difference (in percent of the estimate) of the estimates with and without each observation is plotted. An asterisk is plotted for each influential observation.
nlsJack
returns a list with 7 objects:
estijack |
a data frame with jackknife estimates and bias |
coefjack |
a data frame with the parameter estimates for each jackknife sample |
reldif |
a data frame with the absolute relative difference (in percent of the estimate) of the estimates with and without each observation |
dfb |
a data frame with dfbetas for each parameter and each observation |
jackCI |
a data frame with jackknife confidence intervals |
rse |
a vector with residual standard error for each jackknife sample |
rss |
residual a vector with residual sum of squares for each jackknife sample |
Florent Baty, Marie-Laure Delignette-Muller
Seber GAF, Wild CJ (1989) Nonlinear regression. Wiley, New York.
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.jack1 <- nlsJack(O2K.nls1) plot(O2K.jack1) summary(O2K.jack1)
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.jack1 <- nlsJack(O2K.nls1) plot(O2K.jack1) summary(O2K.jack1)
Provides several plots and tests for the analysis of residuals
nlsResiduals (nls) ## S3 method for class 'nlsResiduals' plot(x, which = 0, ...) test.nlsResiduals (x) ## S3 method for class 'nlsResiduals' print(x, ...)
nlsResiduals (nls) ## S3 method for class 'nlsResiduals' plot(x, which = 0, ...) test.nlsResiduals (x) ## S3 method for class 'nlsResiduals' print(x, ...)
nls |
an object of class 'nls' |
x |
an object of class 'nlsResiduals' |
which |
an integer: |
... |
further arguments passed to or from other methods |
Several plots and tests are proposed to check the validity of the assumptions of the error model based on the analysis of residuals.
The function plot.nlsResiduals
proposes several plots of residuals from the nonlinear fit: plot of non-transformed residuals against fitted values, plot of standardized residuals against fitted values, plot of square root of absolute value of standardized residuals against fitted values, auto-correlation plot of residuals (i+1th residual against ith residual), histogram of the non-transformed residuals and normal Q-Q plot of standardized residuals.test.nlsResiduals
tests the normality of the residuals with the Shapiro-Wilk test (shapiro.test in package stats) and the randomness of residuals with the runs test (Siegel and Castellan, 1988). The runs.test function used in nlstools
is the one implemented in the package tseries
.
nlsResiduals
returns a list of five objects:
std95 |
the Student value for alpha=0.05 (bilateral) and the degree of freedom of the model |
resi1 |
a matrix with fitted values vs. non-transformed residuals |
resi2 |
a matrix with fitted values vs. standardized residuals |
resi3 |
a matrix with fitted values vs. sqrt(abs(standardized residuals)) |
resi4 |
a matrix with ith residuals vs. i+1th residuals |
Florent Baty, Marie-Laure Delignette-Muller
Bates DM and Watts DG (1988) Nonlinear regression analysis and its applications. Wiley, Chichester, UK.
Siegel S and Castellan NJ (1988) Non parametric statistics for behavioral sciences. McGraw-Hill international, New York.
# Plots of residuals formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.res1 <- nlsResiduals(O2K.nls1) plot(O2K.res1, which = 0) # Histogram and qq-plot plot(O2K.res1, which = 5) plot(O2K.res1, which = 6) # Tests test.nlsResiduals(O2K.res1)
# Plots of residuals formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) O2K.res1 <- nlsResiduals(O2K.nls1) plot(O2K.res1, which = 0) # Histogram and qq-plot plot(O2K.res1, which = 5) plot(O2K.res1, which = 6) # Tests test.nlsResiduals(O2K.res1)
Tools to help the fit of nonlinear models with nls
preview (formula, data, start, variable = 1) plotfit (x, smooth = FALSE, variable = 1, xlab = NULL, ylab = NULL, pch.obs = 1, pch.fit = "+", lty = 1, lwd = 1, col.obs = "black", col.fit = "red", ...) overview (x)
preview (formula, data, start, variable = 1) plotfit (x, smooth = FALSE, variable = 1, xlab = NULL, ylab = NULL, pch.obs = 1, pch.fit = "+", lty = 1, lwd = 1, col.obs = "black", col.fit = "red", ...) overview (x)
formula |
formula of a non-linear model |
data |
a data frame with header matching the variables given in the formula |
start |
a list of parameter starting values which names match the parameters given in the formula |
variable |
index of the variable to be plotted against the predicted values; default is the first independent variable as it appears in the orginal dataset |
x |
an object of class 'nls' |
smooth |
a logical value, default is FALSE. If smooth is TRUE, a plot of observed values is plotted as a function of 1000 values continuously taken in the range interval [min(variable),max(variable)]. This option can only be used if the number of controlled variables is 1. |
xlab |
X-label |
ylab |
Y-label |
pch.obs |
type of point of the observed values |
pch.fit |
type of point of the fitted values (not applicable if smooth=TRUE) |
lty |
type of line of the smoothed fitted values (if smooth=TRUE) |
lwd |
thickness of line of the smoothed fitted values (if smooth=TRUE) |
col.obs |
color of the observed points |
col.fit |
color of the fitted values |
... |
further arguments passed to or from other methods |
The function preview
helps defining the parameter starting values prior fitting the model. It provides a superimposed plot of observed (circles) and predicted (crosses) values of the dependent variable versus one of the independent variables with the model evaluated at the starting values of the parameters. The function overview
returns the parameters estimates, their standard errors as well as their asymptotic confidence intervals and the correlation matrix (alternately, the function confint
provides better confidence interval estimates whenever it converges). plotfit
displays a superimposed plot of the dependent variable versus one the independent variables together with the fitted model.
Florent Baty, Marie-Laure Delignette-Muller
Baty F, Ritz C, Charles S, Brutsche M, Flandrois J-P, Delignette-Muller M-L (2015). A Toolbox for Nonlinear Regression in R: The Package nlstools. Journal of Statistical Software, 66(5), 1-21.
Bates DM and Watts DG (1988) Nonlinear regression analysis and its applications. Wiley, Chichester, UK.
nls
in the stats
library and confint.nls
in the package MASS
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) preview(formulaExp, O2K, list(VO2rest = 400, VO2peak = 1600, mu = 1)) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) overview(O2K.nls1) plotfit(O2K.nls1, smooth = TRUE)
formulaExp <- as.formula(VO2 ~ (t <= 5.883) * VO2rest + (t > 5.883) * (VO2rest + (VO2peak - VO2rest) * (1 - exp(-(t - 5.883) / mu)))) preview(formulaExp, O2K, list(VO2rest = 400, VO2peak = 1600, mu = 1)) O2K.nls1 <- nls(formulaExp, start = list(VO2rest = 400, VO2peak = 1600, mu = 1), data = O2K) overview(O2K.nls1) plotfit(O2K.nls1, smooth = TRUE)
The models or data sets listed here are no longer part of package nlstools. In order to access these models and data set in the future, please load the additional package nlsMicrobio.
Defunct functions are:geeraerd
geeraerd_without_Nres
geeraerd_without_Sl
mafart
albert
trilinear
bilinear_without_Nres
bilinear_without_Sl
baranyi
baranyi_without_Nmax
baranyi_without_lag
buchanan
buchanan_without_Nmax
buchanan_without_lag
gompertzm
jameson_buchanan
jameson_baranyi
jameson_without_lag
cpm_T
cpm_pH_4p
cpm_pH_3p
cpm_aw_3p
cpm_aw_2p
cpm_T_pH_aw
competition1
competition2
growthcurve1
growthcurve2
growthcurve3
growthcurve4
ross
survivalcurve1
survivalcurve2
survivalcurve3
Oxygen uptake kinetics during a 6-minute walking test in a patient with pulmonary disease. The first 5.83 minutes correspond to the resting phase prior to exercise.
data(O2K)
data(O2K)
O2K
is a data frame with 2 columns (t: time, VO2: oxygen uptake)
This data set was provided by the Cantonal Hospital St. Gallen, Switzerland.
data(O2K) plot(O2K)
data(O2K) plot(O2K)