Hola!
I have been trying to get thetukey.1 function from one of the s books
to function
under R, but have still problems. Below is my slightly changed version
given.
> tukey.1
function (aov.obj, data)
{
vnames <- names(aov.obj$contrasts)
if (length(vnames) != 2)
stop("The model must be two-way.")
vara <- data[, vnames[1]]
varb <- data[, vnames[2]]
na <- length(levels(vara))
nb <- length(levels(varb))
where.resp <- as.character(attr(aov.obj$terms,
"variables")[attr(aov.obj$terms,
"response")])
# This gets the value "list"!!, so resp below becomes null I cannot
find
a workaround
3 for this. Any suggestions?
resp <- data[, where.resp]
cfs <- coef(aov.obj)
alpha.A <- aov.obj$contrasts[[vnames[1]]] %*%
cfs[aov.obj$assign[aov.obj$assign ==
1]]
alpha.B <- aov.obj$contrasts[[vnames[2]]] %*%
cfs[aov.obj$assign[aov.obj$assign ==
2]]
r.mat <- matrix(0, nb, na)
r.mat[cbind(as.vector(unclass(varb)), as.vector(unclass(vara)))] <-
resp
SS.theta.num <- sum((alpha.B %*% t(alpha.A)) * r.mat)^2
SS.theta.den <- sum(alpha.A^2) * sum(alpha.B^2)
SS.theta <- SS.theta.num/SS.theta.den
SS.res <- sum(resid(aov.obj)^2)
SS.res.1 <- SS.res - SS.theta
T.1df <- ((na * nb - na - nb) * SS.theta)/SS.res.1
p.value <- 1 - pf(T.1df, 1, na * nb - na - nb)
list(T.1df = T.1df, p.value = p.value)
}
Kjetil Halvorsen
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._