From 9e0dae397df8c18e7333d2604cae96b2a7927420 Mon Sep 17 00:00:00 2001 From: ranke Date: Fri, 23 Jun 2006 15:33:27 +0000 Subject: - inverse.predict now has a var.s argument instead of the never tested ss argument. This is documented in the updated vignette - loq() now has w.loq and var.loq arguments, and stops with a message if neither are specified and the model has weights. - calplot doesn't stop any more for weighted regression models, but only refrains from drawing prediction bands - Added method = "din" to lod(), now that I actually have it (DIN 32645) and was able to see which approximation is used therein. - removed the demos, as the examples and tests are already partially duplicated - The vignette is more of a collection of various notes, but should certainly be helpful for the user. - Version bump to 0.1-xxx git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/chemCal@16 5fad18fb-23f0-0310-ab10-e59a3bee62b4 --- man/calplot.lm.Rd | 12 +++++----- man/din32645.Rd | 45 +++++++++++++++++++++----------------- man/inverse.predict.Rd | 25 +++++++++++---------- man/lod.Rd | 20 ++++++++++++----- man/loq.Rd | 58 +++++++++++++++++++++++-------------------------- man/massart97ex3.Rd | 59 +++++++++++++++++++++++++++----------------------- 6 files changed, 118 insertions(+), 101 deletions(-) (limited to 'man') diff --git a/man/calplot.lm.Rd b/man/calplot.lm.Rd index de63022..bf3f616 100644 --- a/man/calplot.lm.Rd +++ b/man/calplot.lm.Rd @@ -5,12 +5,11 @@ \title{Plot calibration graphs from univariate linear models} \description{ Produce graphics of calibration data, the fitted model as well - as prediction and confidence bands. + as confidence, and, for unweighted regression, prediction bands. } \usage{ - calplot(object, xlim = c("auto","auto"), ylim = c("auto","auto"), - xlab = "Concentration", ylab = "Response", alpha=0.05, - varfunc = NULL) + calplot(object, xlim = c("auto", "auto"), ylim = c("auto", "auto"), + xlab = "Concentration", ylab = "Response", alpha=0.05, varfunc = NULL) } \arguments{ \item{object}{ @@ -40,7 +39,8 @@ } \value{ A plot of the calibration data, of your fitted model as well as lines showing - the confidence limits as well as the prediction limits. + the confidence limits. Prediction limits are only shown for models from + unweighted regression. } \note{ Prediction bands for models from weighted linear regression require weights @@ -51,7 +51,7 @@ } \examples{ data(massart97ex3) -m <- lm(y ~ x, data=massart97ex3) +m <- lm(y ~ x, data = massart97ex3) calplot(m) } \author{ diff --git a/man/din32645.Rd b/man/din32645.Rd index 94486c4..cacbf07 100644 --- a/man/din32645.Rd +++ b/man/din32645.Rd @@ -11,39 +11,44 @@ } \examples{ data(din32645) -m <- lm(y ~ x, data=din32645) +m <- lm(y ~ x, data = din32645) calplot(m) -(prediction <- inverse.predict(m,3500,alpha=0.01)) -# This should give 0.07434 according to Dintest test data, as -# collected from Procontrol 3.1 (isomehr GmbH) + +## Prediction of x with confidence interval +(prediction <- inverse.predict(m, 3500, alpha = 0.01)) + +# This should give 0.07434 according to test data from Dintest, which +# was collected from Procontrol 3.1 (isomehr GmbH) in this case round(prediction$Confidence,5) -# According to Dintest test data, we should get 0.0698 for the critical value +## Critical value: +(crit <- lod(m, alpha = 0.01, beta = 0.5)) + +# According to DIN 32645, we should get 0.07 for the critical value # (decision limit, "Nachweisgrenze") -(lod <- lod(m, alpha = 0.01, beta = 0.5)) -round(lod$x,4) +round(crit$x, 2) +# and according to Dintest test data, we should get 0.0698 from +round(crit$x, 4) +## Limit of detection (smallest detectable value given alpha and beta) # In German, the smallest detectable value is the "Erfassungsgrenze", and we -# should get 0.140 according to Dintest test data, but with chemCal, we can't -# reproduce this, -lod(m, alpha = 0.01, beta = 0.01) -# except by using an equivalent to the approximation -# xD = 2 * Sc / A (Currie 1999, p. 118, or Orange Book, Chapter 18.4.3.7) -lod.approx <- 2 * lod$x -round(lod.approx, digits=3) -# which seems to be the pragmatic definition in DIN 32645, as judging from -# the Dintest test data. - -# This accords to the test data from Dintest again, except for the last digit -# of the value cited for Procontrol 3.1 (0.2121) +# should get 0.14 according to DIN, which we achieve by using the method +# described in it: +lod.din <- lod(m, alpha = 0.01, beta = 0.01, method = "din") +round(lod.din$x, 2) + +## Limit of quantification +# This accords to the test data coming with the test data from Dintest again, +# except for the last digits of the value cited for Procontrol 3.1 (0.2121) (loq <- loq(m, alpha = 0.01)) round(loq$x,4) + # A similar value is obtained using the approximation # LQ = 3.04 * LC (Currie 1999, p. 120) 3.04 * lod(m,alpha = 0.01, beta = 0.5)$x } \references{ - DIN 32645 (equivalent to ISO 11843) + DIN 32645 (equivalent to ISO 11843), Beuth Verlag, Berlin, 1994 Dintest. Plugin for MS Excel for evaluations of calibration data. Written by Georg Schmitt, University of Heidelberg. diff --git a/man/inverse.predict.Rd b/man/inverse.predict.Rd index 5be0250..925f3e9 100644 --- a/man/inverse.predict.Rd +++ b/man/inverse.predict.Rd @@ -5,7 +5,7 @@ \alias{inverse.predict.default} \title{Predict x from y for a linear calibration} \usage{inverse.predict(object, newdata, \dots, - ws, alpha=0.05, ss = "auto") + ws, alpha=0.05, var.s = "auto") } \arguments{ \item{object}{ @@ -21,15 +21,17 @@ future implementations. } \item{ws}{ - The weight attributed to the sample. The default is to take the - mean of the weights in the model, if there are any. + The weight attributed to the sample. This argument is obligatory + if \code{object} has weights. } \item{alpha}{ The error tolerance level for the confidence interval to be reported. } - \item{ss}{ - The estimated standard error of the sample measurements. The - default is to take the residual standard error from the calibration. + \item{var.s}{ + The estimated variance of the sample measurements. The default is to take + the residual standard error from the calibration and to adjust it + using \code{ws}, if applicable. This means that \code{var.s} + overrides \code{ws}. } } \value{ @@ -59,12 +61,13 @@ \examples{ data(massart97ex3) attach(massart97ex3) - yx <- split(y,factor(x)) - s <- round(sapply(yx,sd),digits=2) - w <- round(1/(s^2),digits=3) + yx <- split(y, x) + ybar <- sapply(yx, mean) + s <- round(sapply(yx, sd), digits = 2) + w <- round(1 / (s^2), digits = 3) weights <- w[factor(x)] - m <- lm(y ~ x,w=weights) + m <- lm(y ~ x, w = weights) - inverse.predict(m,15,ws = 1.67) + inverse.predict(m, 15, ws = 1.67) # 5.9 +- 2.5 } \keyword{manip} diff --git a/man/lod.Rd b/man/lod.Rd index c2fe4e9..055074c 100644 --- a/man/lod.Rd +++ b/man/lod.Rd @@ -5,7 +5,7 @@ \alias{lod.default} \title{Estimate a limit of detection (LOD)} \usage{ - lod(object, \dots, alpha = 0.05, beta = 0.05) + lod(object, \dots, alpha = 0.05, beta = 0.05, method = "default") } \arguments{ \item{object}{ @@ -24,10 +24,18 @@ \item{beta}{ The error tolerance beta for the detection limit. } + \item{method}{ + The default method uses a prediction interval at the LOD + for the estimation of the LOD, which obviously requires + iteration. This is described for example in Massart, p. 432 ff. + The \dQuote{din} method uses the prediction interval at + x = 0 as an approximation. + } } \value{ A list containig the corresponding x and y values of the estimated limit of - detection of a model used for calibration. } + detection of a model used for calibration. +} \description{ The decision limit (German: Nachweisgrenze) is defined as the signal or analyte concentration that is significantly different from the blank signal @@ -39,7 +47,7 @@ one-sided significance test). } \note{ - - The default values for alpha and beta are recommended by IUPAC. + - The default values for alpha and beta are the ones recommended by IUPAC. - The estimation of the LOD in terms of the analyte amount/concentration xD from the LOD in the signal domain SD is done by simply inverting the calibration function (i.e. assuming a known calibration function). @@ -68,8 +76,8 @@ # The critical value (decision limit, German Nachweisgrenze) can be obtained # by using beta = 0.5: lod(m, alpha = 0.01, beta = 0.5) - # or approximated by - 2 * lod(m, alpha = 0.01, beta = 0.5)$x - # for the case of known, constant variance (homoscedastic data) +} +\seealso{ + Examples for \code{\link{din32645}} } \keyword{manip} diff --git a/man/loq.Rd b/man/loq.Rd index 4850487..0250098 100644 --- a/man/loq.Rd +++ b/man/loq.Rd @@ -5,14 +5,17 @@ \alias{loq.default} \title{Estimate a limit of quantification (LOQ)} \usage{ - loq(object, \dots, alpha = 0.05, k = 3, n = 1, w = "auto") + loq(object, \dots, alpha = 0.05, k = 3, n = 1, w.loq = "auto", + var.loq = "auto") } \arguments{ \item{object}{ A univariate model object of class \code{\link{lm}} or \code{\link[MASS:rlm]{rlm}} with model formula \code{y ~ x} or \code{y ~ x - 1}, - optionally from a weighted regression. + optionally from a weighted regression. If weights are specified + in the model, either \code{w.loq} or \code{var.loq} have to + be specified. } \item{alpha}{ The error tolerance for the prediction of x values in the calculation. @@ -29,53 +32,46 @@ The number of replicate measurements for which the LOQ should be specified. } - \item{w}{ + \item{w.loq}{ The weight that should be attributed to the LOQ. Defaults to one for unweighted regression, and to the mean of the weights for weighted regression. See \code{\link{massart97ex3}} for an example how to take advantage of knowledge about the variance function. } + \item{var.loq}{ + The approximate variance at the LOQ. The default value is + calculated from the model. + } } \value{ The estimated limit of quantification for a model used for calibration. } \description{ - A useful operationalisation of a limit of quantification is simply the - solution of the equation + The limit of quantification is the x value, where the relative error + of the quantification given the calibration model reaches a prespecified + value 1/k. Thus, it is the solution of the equation \deqn{L = k c(L)}{L = k * c(L)} - where c(L) is half of the length of the confidence interval at the limit L as - estimated by \code{\link{inverse.predict}}. By virtue of this formula, the - limit of detection is the x value, where the relative error - of the quantification with the given calibration model is 1/k. + where c(L) is half of the length of the confidence interval at the limit L + (DIN 32645, equivalent to ISO 11843). c(L) is internally estimated by + \code{\link{inverse.predict}}, and L is obtained by iteration. } \note{ - IUPAC recommends to base the LOQ on the standard deviation of the - signal where x = 0. The approach taken here is to my knowledge - original to the chemCal package. + - IUPAC recommends to base the LOQ on the standard deviation of the signal + where x = 0. + - The calculation of a LOQ based on weighted regression is non-standard + and therefore not tested. Feedback is welcome. } \examples{ data(massart97ex3) attach(massart97ex3) - m0 <- lm(y ~ x) - loq(m0) - - # Now we use a weighted regression - yx <- split(y,factor(x)) - s <- round(sapply(yx,sd),digits=2) - w <- round(1/(s^2),digits=3) - weights <- w[factor(x)] - mw <- lm(y ~ x,w=weights) - loq(mw) - - # In order to define the weight at the loq, we can use - # the variance function 1/y for the model - mwy <- lm(y ~ x, w = 1/y) + m <- lm(y ~ x) + loq(m) - # Let's do this with one iteration only - loq(mwy, w = 1 / predict(mwy,list(x = loq(mwy)$x))) - - # We can get better by doing replicate measurements - loq(mwy, n = 3, w = 1 / predict(mwy,list(x = loq(mwy)$x))) + # We can get better by using replicate measurements + loq(m, n = 3) +} +\seealso{ + Examples for \code{\link{din32645}} } \keyword{manip} diff --git a/man/massart97ex3.Rd b/man/massart97ex3.Rd index eb00e79..e7cd383 100644 --- a/man/massart97ex3.Rd +++ b/man/massart97ex3.Rd @@ -3,7 +3,7 @@ \alias{massart97ex3} \title{Calibration data from Massart et al. (1997), example 3} \description{ - Sample dataset to test the package. + Sample dataset from p. 188 to test the package. } \usage{data(massart97ex3)} \format{ @@ -11,37 +11,42 @@ observations of y for each level. } \examples{ -data(massart97ex3) -attach(massart97ex3) -yx <- split(y,x) -ybar <- sapply(yx,mean) -s <- round(sapply(yx,sd),digits=2) -w <- round(1/(s^2),digits=3) -weights <- w[factor(x)] -m <- lm(y ~ x,w=weights) -# The following concords with the book -inverse.predict(m, 15, ws = 1.67) -inverse.predict(m, 90, ws = 0.145) + data(massart97ex3) + attach(massart97ex3) + yx <- split(y, x) + ybar <- sapply(yx, mean) + s <- round(sapply(yx, sd), digits = 2) + w <- round(1 / (s^2), digits = 3) + weights <- w[factor(x)] + m <- lm(y ~ x, w = weights) + calplot(m) -# Some of the following examples are commented out, because the require -# prediction intervals from predict.lm for weighted models, which is not -# available in R at the moment. + # The following concords with the book p. 200 + inverse.predict(m, 15, ws = 1.67) # 5.9 +- 2.5 + inverse.predict(m, 90, ws = 0.145) # 44.1 +- 7.9 -#calplot(m) + # The LOD is only calculated for models from unweighted regression + # with this version of chemCal + m0 <- lm(y ~ x) + lod(m0) -m0 <- lm(y ~ x) -lod(m0) -#lod(m) + # Limit of quantification from unweighted regression + m0 <- lm(y ~ x) + loq(m0) -# Now we want to take advantage of the lower weights at lower y values -#m2 <- lm(y ~ x, w = 1/y) -# To get a reasonable weight for the lod, we need to estimate it and predict -# a y value for it -#yhat.lod <- predict(m,data.frame(x = lod(m2))) -#lod(m2,w=1/yhat.lod,k=3) + # For calculating the limit of quantification from a model from weighted + # regression, we need to supply weights, internally used for inverse.predict + # If we are not using a variance function, we can use the weight from + # the above example as a first approximation (x = 15 is close to our + # loq approx 14 from above). + loq(m, w.loq = 1.67) + # The weight for the loq should therefore be derived at x = 7.3 instead + # of 15, but the graphical procedure of Massart (p. 201) to derive the + # variances on which the weights are based is quite inaccurate anyway. } \source{ - Massart, L.M, Vandenginste, B.G.M., Buydens, L.M.C., De Jong, S., Lewi, P.J., - Smeyers-Verbeke, J. (1997) Handbook of Chemometrics and Qualimetrics: Part A, p. 188 + Massart, L.M, Vandenginste, B.G.M., Buydens, L.M.C., De Jong, S., Lewi, P.J., + Smeyers-Verbeke, J. (1997) Handbook of Chemometrics and Qualimetrics: Part A, + Chapter 8. } \keyword{datasets} -- cgit v1.2.1