[R] tukey.1

[Prof Brian Ripley]

> From: Kjetil Halvorsen <kjetilh at umsanet.edu.bo>
> 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?

Yes, the terms structure is different in R. Try, for example,

> tt <- attr(npk.aov$terms, "variables")
list(yield, block, N, P, K)
> as.character(tt)[-1]
[1] "yield" "block" "N"     "P"     "K"    

so you need (as.character(tt)[-1])[attr(aov.obj$terms, "response")]
(but I think the response is always first in R).

>     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]]

The assign structure is even more different, but maybe you have already
altered this.

A good place to look is aov.R or glm.R which have

        mf <- eval(lmcall, sys.frame(sys.parent()))
        xvars <- as.character(attr(Terms, "variables"))[-1]
        if ((yvar <- attr(Terms, "response")) > 0)
            xvars <- xvars[-yvar]
....



-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272860 (secr)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._