Michael Friendly
2007-Aug-28 20:04 UTC
[R] help with aggregate(): tables of means for terms in an mlm
I'm trying to extend some work in the car and heplots packages that requires getting a table of multivariate means for one (or later, more) terms in an mlm object. I can do this for concrete examples, using aggregate(), but can't figure out how to generalize it. I want to return a result that has the factor-level combinations as rownames, and the means as the body of the table (aggregate returns the factors as initial columns). # Examples: m1 & m2 are desired results > library(car) > soils.mod <- lm(cbind(pH,N,Dens,P,Ca,Mg,K,Na,Conduc) ~ Block + Contour*Depth, data=Soils) > term.names(soils.mod) [1] "(Intercept)" "Block" "Contour" "Depth" [5] "Contour:Depth" > > # response variables > resp<- model.response(model.frame(soils.mod)) > # 1-factor means: term="Contour" > m1<-aggregate(resp, list(Soils$Contour), mean) > rownames(m1) <- m1[,1] > ( m1 <- m1[,-1] ) pH N Dens P Ca Mg K Na Conduc Depression 4.692 0.08731 1.343 188.2 7.101 8.986 0.3781 5.823 6.946 Slope 4.746 0.10594 1.333 159.4 8.109 8.320 0.4156 6.103 6.964 Top 4.570 0.11256 1.272 150.9 8.877 8.088 0.6050 4.872 5.856 > > # 2-factor means: term="Contour:Depth" > m2<-aggregate(resp, list(Soils$Contour, Soils$Depth), mean) > rownames(m2) <- paste(m2[,1], m2[,2],sep=":") > ( m2 <- m2[,-(1:2)] ) pH N Dens P Ca Mg K Na Conduc Depression:0-10 5.353 0.17825 0.9775 333.0 10.685 7.235 0.6250 1.5125 1.473 Slope:0-10 5.508 0.21900 1.0500 258.0 12.248 7.232 0.6350 1.9900 2.050 Top:0-10 5.332 0.19550 1.0025 242.8 13.385 6.590 0.8000 0.9225 1.373 Depression:10-30 4.880 0.08025 1.3575 187.5 7.548 9.635 0.4500 4.6400 5.480 Slope:10-30 5.283 0.10100 1.3475 160.2 9.515 8.980 0.4800 4.9350 4.910 Top:10-30 4.850 0.11750 1.3325 147.5 10.238 8.090 0.6500 2.9800 3.583 Depression:30-60 4.362 0.05050 1.5350 124.2 5.402 9.918 0.2400 7.5875 9.393 Slope:30-60 4.268 0.06075 1.5100 114.5 5.877 8.968 0.3000 7.6300 8.925 Top:30-60 4.205 0.07950 1.3225 116.2 6.620 8.742 0.5450 6.2975 7.440 Depression:60-90 4.173 0.04025 1.5025 108.0 4.770 9.157 0.1975 9.5525 11.438 Slope:60-90 3.927 0.04300 1.4225 105.0 4.798 8.100 0.2475 9.8575 11.970 Top:60-90 3.893 0.05775 1.4300 97.0 5.268 8.928 0.4250 9.2900 11.030 > Here is the current version of a function that doesn't work, because I can't supply the factor names to aggregate in the proper way. Can someone help me make it work? termMeans.mlm <- function( object, term ) { resp<- model.response(model.frame(object)) terms <- term.names(soils.mod) terms <- terms[terms != "(Intercept)"] factors <- strsplit(term, ":") # browser() means <- aggregate(resp, factors, mean) # rownames(means) <- ... # means <- means[, -(1:length(factors)] } > termMeans.mlm(soils.mod, "Contour") Error in FUN(X[[1L]], ...) : arguments must have same length thanks, -Michael -- Michael Friendly Email: friendly AT yorku DOT ca Professor, Psychology Dept. York University Voice: 416 736-5115 x66249 Fax: 416 736-5814 4700 Keele Street http://www.math.yorku.ca/SCS/friendly.html Toronto, ONT M3J 1P3 CANADA