knitr::opts_chunk$set(echo = TRUE)

Import of Data

We import the data and select the variables needed for the analysis.

load(paste0(PATH.OUT, "/CalibrationData.RData"))
tmp <- strsplit(FORMULA, "~")[[1]]
y.var <- names(unlist(sapply(colnames(CalibrationData), 
                             grep, x = tmp[1])))
x.vars <- names(unlist(sapply(colnames(CalibrationData), 
                              grep, x = tmp[2])))
if(SUBSET != ""){
  calData <- eval(call("subset", x =  CalibrationData,
                       subset = parse(text = SUBSET)))
}else{
  calData <- CalibrationData
}
calData <- calData[,c(y.var, x.vars)]
calData

Model

We will apply the following model.

as.formula(FORMULA)

Analysis

We now run a linear regression analysis including stepwise model selection by the Akaike information criterion (AIC).

Full Model

In the fist step, we compute the full model.

library(MASS)
full.fit <- lm(as.formula(FORMULA), data = calData)
summary(full.fit)

Selected Model

We perform a stepwise model selection starting with the full model.

AIC.fit <- stepAIC(full.fit, trace = 0)
AIC.fit
AIC.fit.sum <- summary(AIC.fit)
AIC.fit.sum

We plot the given concentrations against the fitted values.

library(ggplot2)
R2 <- round(AIC.fit.sum$r.squared, 3)
adj.R2 <- round(AIC.fit.sum$adj.r.squared, 3)
DF <- data.frame(Observed = AIC.fit$model[,1], 
                 Fitted = fitted(AIC.fit))
ggplot(DF, aes(x = Observed, y = Fitted)) +
  geom_point() + geom_abline(slope = 1, intercept = 0) +
  xlab("Observed values") + ylab("Fitted values") + 
  annotate("text",  x=min(DF$Observed), y = max(DF$Fitted), 
           label = substitute(paste(R^2, " = ", R2, ", adj. ", R^2, " = ", adj.R2), 
                              list(R2 = R2, adj.R2 = adj.R2)), 
           vjust=1, hjust=0, size = 5)

Computation of LOB, LOD and LOQ

We compute limit of blank (LOB), limit of detection (LOD) and limit of quantification (LOQ) using the linear regression fit. We get the LOB as the upper bound of the one-sided 95\% confidence interval of the intercept. Based on the LOB and the standard error of the residuals, we compute LOD and finally LOQ based on LOD.

LOB <- confint(AIC.fit, level = 0.90)[1,2]
names(LOB) <- "LOB"
LOB
LOD <- LOB + 1.645*summary(AIC.fit)$sigma
names(LOD) <- "LOD"
LOD
LOQ <- 3.3*LOD
names(LOQ) <- "LOQ"
LOQ

Save Results

We save the results.

save(AIC.fit, LOB, LOD, LOQ, file = paste0(PATH.OUT, "/CalibrationResults.RData"))

Software

sessionInfo()


stamats/MultiFlow documentation built on Sept. 7, 2020, 12:21 p.m.