R code for the p-value for td obtained using Dunnett's test
The below gives R code for computing the p-value for Dunnett's test given the number of treatments, number of subjects per group and the error degrees of freedom for the group main effect.
You may need to install the mvtnorm package using
install.packages("mvtnorm")The function for computing the p-value for Dunnett's test is then read in as below
library('mvtnorm')
dunnett_p <- function(td, k, n, df)
{
# td: Dunnett's t-statistic
# k: number of treatments, including the reference
# n: number of subject (assuming equal in all treatments)
# df: degrees of freedom
cMat <- rbind(-1,diag(k-1))
corMat <- t(cMat)%*%(cMat/n)
den <- sqrt(crossprod(t(colSums(cMat^2/n))))
corMat <- corMat / den
p <- pmvt(lower=rep(-td, k-1), upper=rep(td, k-1), delta=rep(0, k-1), df=df, corr=corMat)
1-p[[1]]
}The above function may then be run, for example, as below to give the p-value
dunnett_p(2.3,3,5,23)
[1] 0.05614953