R function parallelismTest

The function parallelismTest tests the parallelism of a set of standard curves in an immunoassay. It takes two arguments: x, which is a vector of dilution factors, and y, which is a matrix of observed values of the assay. The function first calculates the slope and standard error of each standard sample using linear regression. Then it computes the difference and standard error of the difference between every pair of slopes. It also calculates the upper and lower bounds of the confidence interval for the slope differences using the t-distribution. 

The function then plots the reference curves and the confidence interval for the slope differences using barplots. The function also finds the maximum absolute value of the confidence interval bounds, which is used as a limit for parallelism.

The function is useful for testing whether the standard samples have parallel slopes, which is an assumption for immunoassay validation. If the confidence interval for any pair of slopes includes zero, then the slopes are not significantly different and parallelism is not violated.

#Parallelism Test for immunoassay

parallelismTest <- function(x,y,alpha=0.05) {

#x:dilution factor

#y: matrix of the observed value of the assay 

nD <- length(x) #number of dilutions

nS <- dim(y)[2] #number of Standard samples 

Slope <- sapply(1:nS, function(i) summary(mod <- lm(y[,i] ~ x))$coef[2])

SE <- sapply(1:nS, function(i) summary(mod <- lm(y[,i] ~ x))$coef[[4]])#St.Err of slope

dS <- matrix(0,ncol=nS,nrow=nS,byrow=TRUE)

SEdS <- matrix(0,ncol=nS,nrow=nS,byrow=TRUE)

text1 <- matrix("",ncol=nS,nrow=nS,byrow=TRUE)

ym <- rowMeans(y)

mod0 <- lm(ym ~ x)

for (i in 1:(nS-1)) {

for (j in 2:nS) {

dS[i,j] <- Slope[i]-Slope[j]

SEdS[i,j] <- sqrt(SE[i]^2 + SE[j]^2)

text1[i,j] <- paste0(i,"-",j,sep="")

}

}

ut <- upper.tri(dS,diag=FALSE)

dS <- dS[ut] 

SEdS <- SEdS[ut]

text1 <- text1[ut]

Low <- dS-qt(1-alpha/2,2*nD -4)*SEdS

Up <- dS+qt(1-alpha/2,2*nD -4)*SEdS

D <- max(abs(c(Low,Up)))

par(las=0)

leg.text2 <- sapply(1:nS,function(i) paste0("Standard sample ",i,sep=""))

plot(x,ym,xaxt="n",type="o",pch=19,lwd=2,lty=1,main="reference curves",ylim=range(y),xlab=expression(bold("fold dilution")),ylab=expression(bold("Standard Sample readout")))

axis(side=1,at=x)

abline(reg=mod0,lwd=2,col="red",lty=2)

cf <- round(coef(mod0), 2) 

eq <- paste0("meanStandardValue = ", cf[1],ifelse(sign(cf[2])==1, " + ", " - "), abs(cf[2]), " FD ")

mtext(eq, 3, line=-2)

sapply(1:nS,function(i) lines(x,y[,i],type="o",pch=1,col=i,lty=3))

legend("bottomleft",legend=c("mean Standard readout","mean Standard regr.line",leg.text2),bty="n",cex=0.8,

col=c("black","red",1:nS),lty=c(1,2,rep(3,times=nS)),lwd=c(2,2,rep(1,times=nS)),pch=c(19,rep(1,times=nS)))

par(las=2)

barplot(dS,horiz = TRUE,xaxt="n",main="Standard samples:\nconfidence interval",sub=expression(bold("slope differences")),

xlim=c(-1.2*abs(min(c(range(dS,-D)))),1.2*abs(max(c(range(dS,D))))),names.arg=text1)

axis(side=1,at=sort(c(signif(dS,3),signif(-D,3),signif(D,3),0)))

abline(v=c(-D,D),lty=3,col="red")

return(list(Slope=Slope,MeanStandard=ym,diffSlope=dS,limit=D))

}

# EXAMPLE 

dilutions <- c(64,32,16,8,4,2,1)

x <- log(1/dilutions,base=2)

Input <- c(

6.913, 8.360,  7.175, 9.256,

6.018, 7.087, 5.973, 7.813,

5.505, 6.524, 5.262, 7.000,

4.966, 5.200,     5.456, 5.698,

3.440,     4.429, 4.952, 4.616,

3.154, 3.643, 3.657, 3.021,

1.802, 2.256, 3.268, 1.879) 

y<- matrix(Input,ncol=4,nrow=7,byrow=TRUE) 

# Apply the function to the data

result <- parallelismTest(x, y)

# The slope vector shows the slope of each standard sample regression line

print(result$Slope)

[1] -0.8259286 -0.9748214 -0.5951071 -1.2178214

# The mean standard value vector shows the mean of each standard sample readout

print(result$mStd)

[1] 7.92600 6.72275 6.07275 5.33000 4.35925 3.36875 2.30125

# The slope difference matrix shows the difference between each pair of slopes

print(result$dSlope)

[1]  0.1488929 -0.2308214 -0.3797143  0.3918929  0.2430000  0.6227143

# The limit value shows the maximum absolute value of the confidence interval for the slope difference

print(result$limit)

[1] 0.8048761

We get these plots

In this example, all the pairwise contrasts between the slopes are within the confidence limits. This finding indicates that the hypothesis of parallelism is not rejected.