## Generate some p values and compare the three possibilities of
## adjusting for multiple testing
n <- 10000
x <- pmin(rexp(n,rate =1/0.01), 1) # Generate some p values
x[sample(c(F,T), n/10, TRUE)] <- NA # delete some observations
# Three different methods of p value calculation
x1 <- p.adjust(x, method = 'holm')
x2 <- p.adjust(x, method = 'holm', n = length(x))
x3 <- multtest::mt.rawp2adjp(x, proc = 'Holm')
x3 <- x3$adjp[order(x3$index),"Holm"]
# Compare p.adjust and mt.rawp2adjp
par(mfrow=c(1,2))
plot(x1, x3); title('p.adj(x) n vs. mt.rawp2adjp')
plot(x2, x3); title('p.adj(x, len = length(x)) n vs. mt.rawp2adjp')

function (p, method = p.adjust.methods, n = length(p))
{
    method <- match.arg(method)
    if (method == "fdr")
        method <- "BH"
    nm <- names(p)
    p <- as.numeric(p)
    p0 <- setNames(p, nm)
    if (all(nna <- !is.na(p)))
        nna <- TRUE

p <- p[nna]
    lp <- length(p)

stopifnot(n >= lp)
[ remaining source code omitted ]
}