fit.vfp {VFP}

Estimate (Im)Precision-Profiles Modeling the Relationship 'Var ~ Mean'.

Description

This function fits one out of ten or any subset of ten non-linear functions at once. This is done on precision data consisting of information about the variance, concentration at which this variance was observed and the respective degrees of freedom. Provided data must contain at least three columns with this information. There are following variance-functions covered:

• constant variance sigma^2

• constant CV sigma^2=beta_1*u^2

• mixed constant, proportional variance sigma^2=beta_1+beta_2*u^2

• constrained power model, constant exponent sigma^2=(beta_1+beta_2*u)^K

• alternative constrained power model sigma^2=beta_1+beta_2*u^K

• alternative unconstrained power model for VF's with a minimum sigma^2=beta_1+beta_2*u+beta_3*u^J

• alternative unconstrained power model sigma^2=beta_1+beta_2*u^J

• unconstrained power model (default model of Sadler) sigma^2=(beta_1 + beta_2 * u)^J

• CLSI EP17 similar model sigma^2=beta_1 * u^J

• Exact CLSI EP17 model (fit by linear regression on logarithmic scale) cv=beta_1 * u^J

Fitting all ten models is somehow redundant if constant 'C' is chosen to be equal to 2, since models 3 and 5 are equivalent and these are constrained versions of model 7 where the exponent is also estimated. The latter also applies to model 4 which is a constrained version of model 8. Note that model 10 fits the same predictor function as model 9 using simple linear regression on a logarithmic scale. Therefore, regression statistics like AIC, deviance etc. is not directly comparable to that of models 1 to 9. Despite these possible redundancies, as computation time is not critical here for typical precision-profiles (of in-vitro diagnostics precision experiments) we chose to offer batch-processing as well. During computation, all models are internally reparameterized so as to guarantee that the variance function is positive in the range 'u' from 0 to 'u_max'. In models 7 and 8, 'J' is restricted to 0.1<J<10 to avoid the appearance of sharp hooks. Occasionally, these restrictions may lead to a failure of convergence. This is then a sign that the model parameters are on the boundary and that the model fits the data very badly. This should not be taken as reason for concern. It occurs frequently for model 6 when the variance function has no minimum, which is normally the case.

Usage

fit.vfp(
  Data,
  model.no = 7,
  K = 2,
  startvals = NULL,
  quiet = T,
  col.mean = "Mean",
  col.var = "VC",
  col.df = "DF",
  col.sd = NULL,
  col.cv = NULL,
  minVC = NA,
  ...
)

Arguments

Data

(data.frame, matrix) containing mean-values, estimated variances and degrees of freedom for each sample

model.no

(integer) in 1:10, may be any combination of these integer-values specifying models 1 to 10 which are consequently fitted to the data; defaults to 7

K

(numeric) value defining the constant used in models 4 and 5; defaults to 2

startvals

(numeric) vector of start values for the optimization algorithm, only used if a single model is specified by the user, otherwise, start values will be sought automatically

quiet

(logical) TRUE = all output messages (warnings, errors) and regular console output is suppressed and can be found in elements "stderr" and "stdout" of the returned object

col.mean

(character) string specifying a column/variable in 'Data' containing mean-values; defaults to "Mean"

col.var

(character) string specifying a column/variable in 'Data' containing the variances; Defaults to "VC"

col.df

(character) string specifying a column/variable in 'Data' containing the degrees of freedom; defaults to "DF"

col.sd

(character) string (optional) specifying a column/variable in 'Data' containing the SD-values, used for model 10 to derive more precise CV-values compared to derivation from 'col.var' in case 'col.cv' is not specified directly. Note that column "col.var" must nevertheless be set and existing

col.cv

(character) string (optional) specifying a column/variable in 'Data' containing the CV-values, if missing, first 'col.sd' is checked, if missing 'col.var' used to derive per-sample CV-values. Note that column "col.var" must nevertheless be set and existing

minVC

(numeric) value assigned to variances being equal to zero, defaults to NA, which results in removing these observations; could also be set to the smallest possible positive double (.Machine$double.eps)

...

additional parameters passed forward, e.g. 'vc' of function getMat.VCA for selecting a specific variance component in case of 'Data' being a list of 'VCA'-objects (see examples)

Details

Functions predict.VFP and predictMean are useful to make inference based on a fitted model. It is possible to derive concentrations at which a predefined variability is reached, which is sometimes referred to as "functional sensitivity" and/or "limit of quantitation" (LoQ). Funtion predictMean returns the fitted value at which a user-defined variance ("vc"), SD or CV is reached with its corresponding 100(1-alpha)% CI derived from the CI of the fitted model. The plotting method for objects of class 'VFP' can automatically add this information to a plot using arguments 'Prediction' and 'Pred.CI' (see plot.VFP for details. Function predict.VFP makes predictions for specified mean-values based on fitted models.

Note, that in cases where a valid solution was found in the re-parameterized space but the final fit with 'gnm' in the original parameter-space did not converge no variance-covariance matrix can be estimated. Therefore, no confidence-intervals will be available downstream.

Value

(object) of class 'VFP' with elements:

Models

objects of class 'gnm' representing the fitted model

RSS

residual sum of squares for each fitted model

AIC

the Akaike information criterion for each fitted model

Deviance

the deviance for each fitted model

Formulas

formula(s) as expression for each fitted model

Mu.Func

function(s) to transform mean value to re-parameterized scale

Mu

list of transformed mean values

Data

the original data set provided to fit model(s)

K

constant as specified by the user

startvals

start values as specified by the user

errors

messages of all errors caught

output

if 'quiet=TRUE' all output that was redirected to a file

warning

messages of all warnings caught

notes

all other notes caught during execution

Author(s)

Florian Dufey [e-mail], Andre Schuetzenmeister [e-mail]

See Also

plot.VFP, predict.VFP, predictMean

Examples

# load VCA-package and data
library(VCA)
data(VCAdata1)
# perform VCA-anaylsis
lst <- anovaVCA(y~(device+lot)/day/run, VCAdata1, by="sample")
# transform list of VCA-objects into required matrix
mat <- getMat.VCA(lst)		# automatically selects "total"
mat

# fit all 9 models batch-wise
res <- fit.vfp(model.no=1:10, Data=mat)

# if 'mat' is not required for later usage, following works
# equally well
res2 <- fit.vfp(lst, 1:10)

# plot best-fitting model
plot(res)
plot(res, type="cv")
plot(res, type="cv", ci.type="lines", ci.col="red",
		Grid=list(col="wheat"), Points=list(pch=2, lwd=2, col="black"))

# now derive concentation at which a specific reproducibility-
# imprecision of 10\% is reached and add this to the plot
pred <- plot(res, type="cv", ci.type="band", 
				ci.col=as.rgb("red", .25), Grid=list(col="orange"), 
				Points=list(pch=2, lwd=2, col="black"),
				Prediction=list(y=10, col="red"), Pred.CI=TRUE)

# (invisibly) returned object contains all relevant information
pred

# same for repeatability
mat.err <- getMat.VCA(lst, "error")
res.err <- fit.vfp(1:10, Data=mat.err)

# without extracting 'mat.err'
res.err2 <- fit.vfp(lst, 1:10, vc="error")

plot(res.err)

#######################################################################
# another example using CA19_9 data from CLSI EP05-A3

data(CA19_9)

# fit reproducibility model to data
fits.CA19_9 <- anovaVCA(result~site/day, CA19_9, by="sample")

# fit within-laboratory-model treating site as fixed effect
fits.ip.CA19_9 <- anovaMM(result~site/(day), CA19_9, by="sample")

# the variability "between-site" is not part of "total"
fits.ip.CA19_9[[1]]
fits.CA19_9[[1]]

# extract repeatability
rep.CA19_9 <- getMat.VCA(fits.CA19_9, "error")

# extract reproducibility
repro.CA19_9 <- getMat.VCA(fits.CA19_9, "total")

# extract intermediate-precision (within-lab)
ip.CA19_9 <- getMat.VCA(fits.ip.CA19_9, "total")

# fit model (a+bX)^C (model 8) to all three matrices

mod8.repro 	<- fit.vfp(repro.CA19_9, 8)
mod8.ip		<- fit.vfp(ip.CA19_9, 8)
mod8.rep		<- fit.vfp(rep.CA19_9, 8)

# plot reproducibility precision profile first
# leave enough space in right margin for a legend
plot(mod8.repro, mar=c(5.1, 7, 4.1, 15), 
		type="cv", ci.type="none", Model=FALSE,
		Line=list(col="blue", lwd=3), 
		Points=list(pch=15, col="blue", cex=1.5),  
		xlim=c(10, 450), ylim=c(0,10),
		Xlabel=list(text="CA19-9, kU/L (LogScale) - 3 Patient Pools, 3 QC Materials",
				cex=1.5), Title=NULL,
		Ylabel=list(text="% CV", cex=1.5),
		Grid=NULL, Crit=NULL, log="x")

# add intermediate precision profile
plot	(mod8.ip, type="cv", add=TRUE, ci.type="none",
		Points=list(pch=16, col="deepskyblue", cex=1.5),
		Line=list(col="deepskyblue", lwd=3), log="x")

# add repeatability precision profile
plot(	mod8.rep, type="cv", add=TRUE, ci.type="none",
		Points=list(pch=17, col="darkorchid3", cex=1.5),
		Line=list(col="darkorchid3", lwd=3), log="x")

# add legend to right margin
legend.rm( x="center", pch=15:17, col=c("blue", "deepskyblue", "darkorchid3"),
		cex=1.5, legend=c("Reproducibility", "Within-Lab Precision", "Repeatability"),
		box.lty=0)
#'