Description Usage Arguments Author(s) Examples
This function performs Tukey's single degree of freedom test for interaction in an unreplicated two-factor design
1 |
data |
input - this is a data frame with three variables, the first variable is a numeric response and next two variables are factors. There should be ab lines in the data frame where a is the number of levels of the first factor, and b is the number of levels of the second factor. |
John Lawson
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 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 | library(daewr)
data(virus)
Tukey1df(virus)
## The function is currently defined as
function (data)
{
y <- data[, 1]
Afactor <- data[, 2]
Bfactor <- data[, 3]
tst1 <- is.factor(Afactor)
tst2 <- is.factor(Bfactor)
tst3 <- is.numeric(y)
if (tst1 & tst2 & tst3) {
a <- nlevels(Afactor)
b <- nlevels(Bfactor)
}
else {
stop("The first column of the data frame is the numeric response,
the 2nd and 3rd columns should be coded as factors")
}
tst4 <- max(a, b) > 2
tst5 <- length(y) == a * b
if (tst4 & tst5) {
ybb <- with(data, tapply(y, Bfactor, mean))
yba <- with(data, tapply(y, Afactor, mean))
sbb <- with(data, tapply(y, Bfactor, sum))
sba <- with(data, tapply(y, Afactor, sum))
ybardd <- mean(y)
CT <- (sum(y)^2)/(a * b)
ssA <- sum(sba^2/b) - CT
ssB <- sum(sbb^2/a) - CT
ssE <- sum(y^2) - CT - ssA - ssB
ybdj <- rep(ybb, 6)
prody <- y * ybdj
sumprod <- tapply(prody, Afactor, sum)
leftsum <- sum(sumprod * yba)
ssAB <- (a * b * (leftsum - (ssA + ssB + a * b * ybardd^2) *
ybardd)^2/(ssA * ssB))
ssR <- ssE - ssAB
F <- ssAB/(ssR/((a - 1) * (b - 1) - 1))
Pval <- 1 - pf(F, 1, ((a - 1) * (b - 1) - 1))
cat("Source df SS MS F Pr>F",
"\n")
cat("A ", paste(format(a - 1, width = 6),
" ", format(round(ssA, 4), justify = "right"), " ",
format(round(ssA/(a - 1), 4), justify = "right"),
"\n"), sep = "")
cat("B ", paste(format(b - 1, width = 6),
" ", format(round(ssB, 4), justify = "right"), " ",
format(round(ssB/(b - 1), 4), justify = "right"),
"\n"), sep = "")
cat("Error ", paste(format((b - 1) * (a - 1),
width = 6), " ", format(round(ssE, 4), justify = "right"),
" ", format(round(ssE/(a - 1) * (b - 1), 4), justify = "right"),
"\n"), sep = "")
cat("NonAdditivity", paste(format(1, width = 6), " ",
format(round(ssAB, 4), justify = "right"), " ",
format(round(ssAB, 4), justify = "right"), " ",
format(round(F, 2), justify = "right"), " ", format(round(Pval,
4), justify = "right"), "\n"), sep = "")
cat("Residual ", paste(format((b - 1) * (a - 1) -
1, width = 6), " ", format(round(ssR, 4), justify = "right"),
" ", format(round(ssR/((a - 1) * (b - 1) - 1), 4),
justify = "right"), "\n"), sep = "")
}
else {
stop("This function only works for unreplicated 2-factor
factorials with >2 levels for one of the factors")
}
}
|
Source df SS MS F Pr>F
A 5 0.1948 0.039
B 2 3.1664 1.5832
Error 10 0.1283 0.0513
NonAdditivity 1 0.0069 0.0069 0.51 0.4932
Residual 9 0.1214 0.0135
function (data)
{
y <- data[, 1]
Afactor <- data[, 2]
Bfactor <- data[, 3]
tst1 <- is.factor(Afactor)
tst2 <- is.factor(Bfactor)
tst3 <- is.numeric(y)
if (tst1 & tst2 & tst3) {
a <- nlevels(Afactor)
b <- nlevels(Bfactor)
}
else {
stop("The first column of the data frame is the numeric response, \n\t\tthe 2nd and 3rd columns should be coded as factors")
}
tst4 <- max(a, b) > 2
tst5 <- length(y) == a * b
if (tst4 & tst5) {
ybb <- with(data, tapply(y, Bfactor, mean))
yba <- with(data, tapply(y, Afactor, mean))
sbb <- with(data, tapply(y, Bfactor, sum))
sba <- with(data, tapply(y, Afactor, sum))
ybardd <- mean(y)
CT <- (sum(y)^2)/(a * b)
ssA <- sum(sba^2/b) - CT
ssB <- sum(sbb^2/a) - CT
ssE <- sum(y^2) - CT - ssA - ssB
ybdj <- rep(ybb, 6)
prody <- y * ybdj
sumprod <- tapply(prody, Afactor, sum)
leftsum <- sum(sumprod * yba)
ssAB <- (a * b * (leftsum - (ssA + ssB + a * b * ybardd^2) *
ybardd)^2/(ssA * ssB))
ssR <- ssE - ssAB
F <- ssAB/(ssR/((a - 1) * (b - 1) - 1))
Pval <- 1 - pf(F, 1, ((a - 1) * (b - 1) - 1))
cat("Source df SS MS F Pr>F",
"\n")
cat("A ", paste(format(a - 1, width = 6),
" ", format(round(ssA, 4), justify = "right"), " ",
format(round(ssA/(a - 1), 4), justify = "right"),
"\n"), sep = "")
cat("B ", paste(format(b - 1, width = 6),
" ", format(round(ssB, 4), justify = "right"), " ",
format(round(ssB/(b - 1), 4), justify = "right"),
"\n"), sep = "")
cat("Error ", paste(format((b - 1) * (a - 1),
width = 6), " ", format(round(ssE, 4), justify = "right"),
" ", format(round(ssE/(a - 1) * (b - 1), 4), justify = "right"),
"\n"), sep = "")
cat("NonAdditivity", paste(format(1, width = 6), " ",
format(round(ssAB, 4), justify = "right"), " ",
format(round(ssAB, 4), justify = "right"), " ",
format(round(F, 2), justify = "right"), " ", format(round(Pval,
4), justify = "right"), "\n"), sep = "")
cat("Residual ", paste(format((b - 1) * (a - 1) -
1, width = 6), " ", format(round(ssR, 4), justify = "right"),
" ", format(round(ssR/((a - 1) * (b - 1) - 1), 4),
justify = "right"), "\n"), sep = "")
}
else {
stop("This function only works for unreplicated 2-factor \n\t\tfactorials with >2 levels for one of the factors")
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.