CV from CI, pooling, sample size [🇷 for BE/BA]
Dear all,
following a question at LinkedIn’s non-public user-group “BA/BE professionals”, code which combines some of
Which gives:
Doing the math is just the first step. Before you pool CVs I would suggest to inspect whether the confidence intervals of the CVs overlap (columns 6–7). If not, try to find out why (different CROs, populations, bioanalytical methods, …). Use common sense to decide which CVs are reliable enough to pool.
Another question at LinkedIn was whether one should use the highest CV of several studies for conservatism. Maybe, maybe not. If the sample sizes were very close and the CVs not too different, why not. Otherwise (or if different designs were used) I prefer the upper CL of the pooled CV.
Note:
following a question at LinkedIn’s non-public user-group “BA/BE professionals”, code which combines some of
PowerTOST
’s functions to- calculate the CVs of several studies by the function
CVfromCI()
,
- calculate 95% CI of the CVs by the
function CVCL()
,
- pool the CVs and calculate the 80% (upper/lower) CL of the pooled CV by the function
CVpooled()
,
- check whether the studies were balanced, and
- estimate the sample size of a planned study (given its design, expected T/R-ratio, desired power, and regulatory body) by the functions
sampleN.TOST()
– crossover, paired, parallel;sampleN.scABEL()
– for the EMA or Health Canada;sampleN.RSABE()
– for the FDA. The FDA’s reference-scaling method for NTIDs bysampleN.NTIDFDA()
is not implemented yet.
###################
# your data below #
###################
design <- c("3x6x3", "2x2x2", "2x2x4")
N <- c(17, 34, 33) # (total) sample sizes
lower <- c(0.8666, 1.0211, 0.9750) # lower CL
upper <- c(1.2845, 1.2665, 1.1390) # upper CL
CI <- c(95, 90, 90) # level of the CI in %
theta0 <- 1.10 # expected T/R-ratio
target <- 0.80 # target power
des <- "2x2x4" # planned design
RSABE <- FALSE # reference-scaling desired?
reg <- "EMA" # any of "EMA", "FDA", "HC"
adjust <- FALSE # iteratively adjust alpha to
# prevent inflation of the TIE?
##########################################################
# don't change below unless you know what your are doing #
##########################################################
library(PowerTOST)
alpha <- (1-CI/100)/2
digits <- max(nchar(as.character(c(lower, upper))))-2
if (RSABE) RSABE.des <- "yes" else RSABE.des <- "no"
PE <- round(sqrt(lower*upper), digits)
CL <- df <- calc <- bal <- CVs <- CLlo <- CLhi <-vector()
# calculate the CVs and their confidence intervals
for (j in seq_along(N)) {
CL[j] <- paste0(CI[j], "%")
n <- N[j]
df[j] <- eval(parse(text=known.designs()[which(known.designs()[["design"]] == design[j]), "df"],
srcfile=NULL))
calc[j] <- known.designs()[which(known.designs()[["design"]] == design[j]), "df"]
steps <- known.designs()[which(known.designs()[["design"]] == design[j]), "steps"]
if (n %% steps == 0) bal[j] <- "yes" else bal[j] <- "no"
CVs[j] <- CVfromCI(upper=lower[j], lower=upper[j], n=N[j], design=design[j], alpha=alpha[j])
CLlo[j] <- CVCL(CV=CVs[j], df=df[j], side="2-sided", alpha=0.05)[[1]]
CLhi[j] <- CVCL(CV=CVs[j], df=df[j], side="2-sided", alpha=0.05)[[2]]
}
CVdata <- data.frame(PE, CL, lower, upper, CVs, CLlo, CLhi, N, design,
paste("study", 1:length(N)), bal, calc, paste("=", df))
names(CVdata) <- c("PE", "CI", "lo", "hi", "CV", "lower CL", "upper CL",
"n", "design", "source", "balanced", "df", "")
# calculate the pooled CV and the 80% upper one-sided CL
CVpooled <- CVpooled(CVdata[, c(5, 8:10)], alpha=0.2)
CVdata[, 5:7] <- round(CVdata[, 5:7], digits) # precision like input
# estimate sample size and expected power for the pooled CV
# and its upper CL
alpha.1 <- rep(0.05, 2)
if (RSABE) { # reference-scaling
CL.lo <- CVCL(CV=CVpooled[["CV"]], df=CVpooled[["df"]],
side="lower", alpha=0.2)[["lower CL"]]
CL.hi <- CVCL(CV=CVpooled[["CV"]], df=CVpooled[["df"]],
side="upper", alpha=0.2)[["upper CL"]]
txt <- paste("Pooled CV =", round(CVpooled[["CV"]], digits),
"with", CVpooled[["df"]], "degrees of freedom")
if (reg == "FDA") { # the FDA's RSABE
if (des == "2x3x3" | des == "2x2x4" | des == "2x2x3") {
tmp <- sampleN.RSABE(CV=CVpooled[["CV"]], theta0=theta0,
targetpower=target, design=des,
details=FALSE, print=FALSE)
n.est1 <- tmp[["Sample size"]]
pwr.est1 <- round(tmp[["Achieved power"]], digits)
L.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["lower"]]
U.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["upper"]]
tmp1 <- sampleN.RSABE(CV=CL.lo, theta0=theta0,
targetpower=target, design=des,
details=FALSE, print=FALSE)
tmp2 <- sampleN.RSABE(CV=CL.hi, theta0=theta0,
targetpower=target, design=des,
details=FALSE, print=FALSE)
if (tmp2[["Sample size"]] >= tmp1[["Sample size"]]) { # critical: upper CL
n.est2 <- tmp2[["Sample size"]]
pwr.est2 <- round(tmp2[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.hi, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.hi, regulator=reg)[["upper"]]
txt <- paste(txt, "\nUpper 80% confidence limit of the CV =",
round(CL.hi, digits), "\n")
} else { # critical: lower CL
n.est2 <- tmp1[["Sample size"]]
pwr.est2 <- round(tmp1[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.lo, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.lo, regulator=reg)[["upper"]]
txt <- paste(txt, "\nUpper 80% confidence limit of the CV =",
round(CL.lo, digits), "\n")
}
} else {
stop("Design not implemented.")
}
} else { # the EMA's and Health Canada's ABEL
if (des == "2x3x3" | des == "2x2x4" | des == "2x2x3") {
if (!adjust) { # according to the GL: alpha 0.05
tmp <- sampleN.scABEL(CV=CVpooled[["CV"]], theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
n.est1 <- tmp[["Sample size"]]
pwr.est1 <- round(tmp[["Achieved power"]], digits)
L.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["lower"]]
U.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["upper"]]
tmp1 <- sampleN.scABEL(CV=CL.lo, theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
tmp2 <- sampleN.scABEL(CV=CL.hi, theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
if (tmp2[["Sample size"]] >= tmp1[["Sample size"]]) { # critical: upper CL
n.est2 <- tmp2[["Sample size"]]
pwr.est2 <- round(tmp2[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.hi, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.hi, regulator=reg)[["upper"]]
txt <- paste(txt, "\nUpper 80% confidence limit of the CV =",
round(CL.hi, digits), "\n")
} else { # critical: lower CL
n.est2 <- tmp1[["Sample size"]]
pwr.est2 <- round(tmp1[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.lo, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.lo, regulator=reg)[["upper"]]
txt <- paste(txt, "\nUpper 80% confidence limit of the CV =",
round(CL.lo, digits), "\n")
}
} else { # iteratively adjusted alpha
tmp <- sampleN.scABEL.ad(CV=CVpooled[["CV"]], theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
n.est1 <- tmp[["Sample size"]]
pwr.est1 <- round(tmp[["Achieved power"]], digits)
L.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["lower"]]
U.1 <- scABEL(CV=CVpooled[["CV"]], regulator=reg)[["upper"]]
if (!is.na(tmp[["adj. alpha"]])) alpha.1[1] <- tmp[["adj. alpha"]]
tmp1 <- sampleN.scABEL.ad(CV=CL.lo, theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
tmp1 <- sampleN.scABEL.ad(CV=CL.hi, theta0=theta0,
targetpower=target, design=des,
regulator=reg, details=FALSE, print=FALSE)
if (tmp2[["Sample size"]] >= tmp1[["Sample size"]]) { # critical: upper CL
n.est2 <- tmp[["Sample size"]]
pwr.est2 <- round(tmp[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.hi, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.hi, regulator=reg)[["upper"]]
} else { # critical: lower CL
n.est2 <- tmp1[["Sample size"]]
pwr.est2 <- round(tmp1[["Achieved power"]], digits)
L.2 <- scABEL(CV=CL.lo, regulator=reg)[["lower"]]
U.2 <- scABEL(CV=CL.lo, regulator=reg)[["upper"]]
}
if (!is.na(tmp[["adj. alpha"]])) alpha.1[2] <- tmp[["adj. alpha"]]
}
} else {
stop("Design not implemented.")
}
}
} else { # conventional ABE
tmp <- sampleN.TOST(CV=CVpooled[["CV"]], theta0=theta0,
targetpower=target, design=des, print=FALSE)
n.est1 <- tmp[["Sample size"]]
pwr.est1 <- round(tmp[["Achieved power"]], digits)
tmp <- sampleN.TOST(CV=CVpooled[["CVupper"]], theta0=theta0,
targetpower=target, design=des, print=FALSE)
n.est2 <- tmp[["Sample size"]]
pwr.est2 <- round(tmp[["Achieved power"]], digits)
L.2 <- L.1 <- 0.8
U.2 <- U.1 <- 1.25
txt <- paste("Pooled CV =", round(CVpooled[["CV"]], digits), "with",
CVpooled[["df"]], "degrees of freedom",
"\nUpper 80% confidence limit of the CV =",
round(CVpooled[["CVupper"]], digits), "\n")
}
Nest <- data.frame(rep(des, 2), rep(RSABE.des, 2), rep(reg, 2),
rep(sprintf("%.4f", theta0), 2),
sprintf("%.4f", c(L.1, L.2)),
sprintf("%.4f", c(U.1, U.2)), alpha.1,
round(c(CVpooled[["CV"]], CVpooled[["CVupper"]]), digits),
c(n.est1, n.est2), c(pwr.est1, pwr.est2))
names(Nest) <- c("Design", "RSABE desired", "regulator", "theta0",
"theta1", "theta2", "alpha", "CV", "n", "Achieved power")
if (length(unique(design)) > 1 & "parallel" %in% design)
warning("Mixing crossover and parallel studies = nonsense!")
print(as.data.frame(CVdata), row.names=F);cat(txt);print(as.data.frame(Nest), row.names=FALSE)
Which gives:
PE CI lo hi CV lower CL upper CL n design source balanced df
1.0551 95% 0.8666 1.2845 0.2832 0.2247 0.3844 17 3x6x3 study 1 no 2*n-4 = 30
1.1372 90% 1.0211 1.2665 0.2667 0.2132 0.3574 34 2x2x2 study 2 yes n-2 = 32
1.0538 90% 0.9750 1.1390 0.2736 0.2386 0.3210 33 2x2x4 study 3 no 3*n-4 = 95
Pooled CV = 0.2741 with 157 degrees of freedom
Upper 80% confidence limit of CV = 0.2888
Design RSABE desired regulator theta0 theta1 theta2 alpha CV n Achieved power
2x2x4 no EMA 1.1000 0.8000 1.2500 0.05 0.2741 28 0.8015
2x2x4 no EMA 1.1000 0.8000 1.2500 0.05 0.2888 32 0.8135
Doing the math is just the first step. Before you pool CVs I would suggest to inspect whether the confidence intervals of the CVs overlap (columns 6–7). If not, try to find out why (different CROs, populations, bioanalytical methods, …). Use common sense to decide which CVs are reliable enough to pool.
Another question at LinkedIn was whether one should use the highest CV of several studies for conservatism. Maybe, maybe not. If the sample sizes were very close and the CVs not too different, why not. Otherwise (or if different designs were used) I prefer the upper CL of the pooled CV.
Note:
- Since in the example two studies are unbalanced, the function
CVfromCI()
throws these messages:
Unbalanced 3x6x3 design. n(i)= 3/3/3/3/3/2 assumed.
Unbalanced 2x2x4 design. n(i)= 17/16 assumed.
CVfromCI()
tries to keep sequences as balanced as possible – which gives a conservative estimate of the CV if the study was more unbalanced. If you know the subjects / sequences calculate the CV and useCVpooled()
“manually” (seehelp(CVfromCI)
andhelp(CVpooled)
how to do that).
- The CV calculated from the CI of a replicate design is always pooled from CVwR and CVwT. Since they are not necessarily equal, be aware.
- If you opt for reference-scaling (
RSABE <- TRUE
) the upper CL of the CV might not be the most conservative solution. It is possible that a lower CV results in less power (since less scaling is applicable) and therefore, requires a higher sample size. The script checks for both possibilities and shows the worst case (lower/upper CL). I suggest to explore how power behaves. After executing this script run:
plot(pa.scABE(CV=round(CVpooled[["CV"]], digits), theta0=theta0,
targetpower=target, design=des, regulator=reg),
pct=FALSE)
And:
if (tmp2[["Sample size"]] >= tmp1[["Sample size"]]) {
plot(pa.scABE(CV=round(CL.hi, digits), theta0=theta0,
targetpower=target, design=des, regulator=reg),
pct=FALSE)
} else {
plot(pa.scABE(CV=round(CL.lo, digits), theta0=theta0,
targetpower=target, design=des, regulator=reg),
pct=FALSE)
}
—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz
The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz
The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes
Complete thread:
- CV from CI, pooling, sample sizeHelmut 2016-07-27 14:09 [🇷 for BE/BA]