Hi,
my name is Marco Bressan i'm working to improve ADati package. I study
psicology ad Padua University (Italy). I have this problem: why bartlett.test
function running good and my anova.welch function no?
Ciao,
il mio nome è Marco Bressan e sto lavorando per migliorare il pacchetto ADati.
Studio psicologia all'università di Padova. Non capisco come mai la funzione
che ho fatto mi dia quell'errore
this is anova.welch:
anova.welch <- function(x, ...) UseMethod("anova.welch")
anova.welch.default <-
## this is the algoritm, I think it's ok (I copy this from Welch() inside ADati)
function (x, y = NULL, nu = c(0,0) ,...)
{
mx <- tapply(x,y,mean)
s2x <- tapply(x,y,var)
k <- length(nu)
w <- nu/s2x
Xp <- sum(w * mx)/sum(w)
Fnum <- sum(w * (mx - Xp)^2)/(k - 1)
Fden <- 0
for (h in 1:k) {
a <- 1/(nu[h] - 1)
b <- (1 - (w[h]/sum(w)))^2
Fden <- Fden + (a * b)
}
gl2.den <- Fden * 3
Fden <- Fden * ((2 * (k - 2))/(k^2 - 1)) + 1
Fw <- Fnum/Fden
STATISTIC <- Fw
gl1 <- k - 1
gl2.num <- k^2 - 1
gl2 <- gl2.num/gl2.den
PARAMETER <- c(gl1, gl2)
PVAL <- pf(Fw, gl1, gl2, lower.tail = FALSE)
METHOD <- "Welch ANOVA"
DNAME <- NA
names(STATISTIC) <- "F"
names(PARAMETER) <- c("num df", "denom df")
RVAL <- list(statistic = STATISTIC, parameter = PARAMETER,
p.value = PVAL, method = METHOD)
class(RVAL) <- "htest"
return(RVAL)
}
anova.welch.formula <-
## I copy this from bartlett.test
function(formula, data, subset, na.action, ...)
{
if(missing(formula) || (length(formula) != 3L))
stop("'formula' mancante o incorretta")
m <- match.call(expand.dots = FALSE)
if(is.matrix(eval(m$data, parent.frame())))
m$data <- as.data.frame(data)
m[[1L]] <- as.name("model.frame")
mf <- eval(m, parent.frame())
DNAME <- paste(names(mf), collapse = " by ")
names(mf) <- NULL
y <- do.call("anova.welch", as.list(mf))
y$data.name <- DNAME
y
}
> n1=10
> ## this is the test
> n2=15
> n3=20
> y=c(rnorm((n1+n2),5,2),rnorm(n3,7,8))
> A=factor(c(rep(1,n1),rep(2,n2),rep(3,n3)))
> anova.welch(y,A,c(n1,n2,n3))
Welch ANOVA
data:
F = 2.3025, num df = 2.000, denom df = 27.384, p-value = 0.1191
> anova.welch(y~A,nu=c(n1,n2,n3))
Errore in model.frame.default(formula = y ~ A, ... = list(nu = c(n1, n2, :
invalid type (pairlist) for variable '(...)'
Sorry for my english, tanks you if you can help me :)
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
