On Aug 26, 2010, at 8:47 AM, Bryan Hanson wrote:
Hello Again Gurus and Lurkers:
I’m trying to build a very user-friendly function which does aov
without
having the user type in a formula (which would be tedious in this
case).
The idea is to take the response from a PCA score matrix, and the
factors
from a list. A simple example is the function given below, along
with test
data and a sample call to the function.
I'm certainly having trouble understanding the proper ways to work
with
formulas and related items, but I think what I do in the function
should
work (it's built on ideas dug out of the archives). However, when
the data
is passed to aov (directly or via manova), something in the bowels
of aov
complains with the following error:
Error in model.frame.default(formula = form, drop.unused.levels =
TRUE) :
object is not a matrix
Actually I got a different error (also on a Mac. but with a version
that is two month later than yours and with quite a smaller number of
packages loaded):
> hypTestScores(mylist = td2, score.matrix = td1,
+ fac = c("f1", "f2"))
Error in eval(expr, envir, enclos) : object 'scores' not found
Which went away if I removed the "env=" argument. My guess is that you
were telling manova to look further up the lexical tree than it should
have (or that attach() operation messed up the environment somehow):
> hypTestScores(mylist = td2, score.matrix = td1,
+ fac = c("f1", "f2"))
Df Pillai approx F num Df den Df Pr(>F)
f1 1 0.31663 0.61777 3 4 0.6392
f2 1 0.36652 0.77145 3 4 0.5673
f1:f2 1 0.15687 0.24808 3 4 0.8593
Residuals 6
> detach("mylist") # needed if there is an error
Error in detach("mylist") : invalid 'name' argument
You might want to do inside the function something like
res =try( your function )
if (res="try-error") {stop() } else{
<process>}
To me, the formula looks legitimate, and the variables in the
formula are
all in the environment (I think: The way I am doing this is
basically that
described in ?summary.manova where only a formula is passed, no data
argument). Based upon reading the archives, the problem might arise
with
one of the deparse statements in aov, but I can't resolve it. It
might also
be one of scoping/environment, but again, this is only an idea.
TIA for any assistance. Bryan
*************
Bryan Hanson
Professor of Chemistry & Biochemistry
DePauw University, Greencastle IN USA
hypTestScores <-
function(mylist, score.matrix, pcs = 1:3, fac = NULL, ...) {
scores <- score.matrix[,pcs]
# str(scores) # looks correct to me
form <- as.formula(paste("scores", "~", paste(fac, collapse =
"*")),
env = parent.frame())
# str(form) # looks correct to me
attach(mylist)
if (length(pcs) > 1) out <- manova(formula = form, ...)
if (length(pcs) == 1) out <- aov(formula = form, ...)
print(summary(out))
detach(mylist)
invisible(out)
}
# test data
td1 <- matrix(rnorm(50), ncol = 5) # like PCA scores
td2 <- list(
f1 = as.factor(sample(c("A", "B"), 10, replace = TRUE)),
f2 = as.factor(sample(c("C", "D"), 10, replace = TRUE)))
# test call
hypTestScores(mylist = td2, score.matrix = td1,
fac = c("f1", "f2"))
detach("mylist") # needed if there is an error
sessionInfo()
R version 2.11.0 (2010-04-22)
x86_64-apple-darwin9.8.0
locale:
[1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] datasets tools grid graphics grDevices utils stats
[8] methods base
other attached packages:
[1] faraway_1.0.4 GGally_0.2 xtable_1.5-6
[4] mvbutils_2.5.1 ggplot2_0.8.8 proto_0.3-8
[7] reshape_0.8.3 ChemoSpec_1.45 R.utils_1.4.0
[10] R.oo_1.7.2 R.methodsS3_1.2.0 rgl_0.91
[13] lattice_0.18-5 mvoutlier_1.4 plyr_1.0.3
[16] RColorBrewer_1.0-2 chemometrics_0.8 som_0.3-5
[19] robustbase_0.5-0-1 rpart_3.1-46 pls_2.1-0
[22] pcaPP_1.8-1 mvtnorm_0.9-9 nnet_7.3-1
[25] mclust_3.4.4 MASS_7.3-5 lars_0.9-7
[28] e1071_1.5-23 class_7.3-2
______________________________________________
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.
David Winsemius, MD
West Hartford, CT
______________________________________________
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.