1 |
data |
|
cls1 |
|
cls2 |
|
shrinkage |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 | ##---- Should be DIRECTLY executable !! ----
##-- ==> Define data, use random,
##-- or do help(data=index) for the standard data sets.
## The function is currently defined as
function (data, cls1, cls2, shrinkage = F)
{
meanS1 <- apply(data[, cls1], 1, mean)
meanS2 <- apply(data[, cls2], 1, mean)
covS1 <- covS2 <- NULL
if ((nrow(data) - 1) > (ncol(data))) {
shrinkage = T
}
if (shrinkage) {
covS1 <- cov.shrink(t(data[, cls1]), verbose = F)
covS2 <- cov.shrink(t(data[, cls2]), verbose = F)
}
else {
covS1 <- cov(t(data[, cls1]))
covS2 <- cov(t(data[, cls2]))
}
n1 = ncol(data[, cls1])
n2 = ncol(data[, cls2])
p <- nrow(data)
covS = ((n1 - 1) * covS1 + (n2 - 1) * covS2)/(n1 + n2 - 2)
covSinv <- solve(covS)
htsqx <- (t(meanS1 - meanS2) %*% covSinv %*% (meanS1 - meanS2)) *
(n1 * n2/(n1 + n2))
list(statistic = htsqx, df1 = p, df2 = n1 + n2 - p - 1, n1 = n1,
n2 = n2)
}
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.