rlapplyancova

How to use stack() function correctly to extract marginal means of ANCOVA in R?

I'm a beginner of R, and I would like to extract marginal means from ANCOVA tests performed on over 200 outcome variables. It worked well when I only used stack() on only one outcome variable, but I got error when I use both stack() and lapply().

Here I use the built-in dataset "iris" to display the problem. The dataset "iris" has three levels at Species, and I use Petal.Width as covariate, Species as predictive variable, and the first three columns of variables as outcome variable.

My purpose is to extract multiple marginal means of the corresponding outcome variables at the same time rather than perform extraction one by one.

#load data and packages
data("iris")
library(car); library(compute.es); library(effects); library(ggplot2);
library(multcomp); library(pastecs); library(WRS)

#set contrasts for the following ANCOVA tests
contrasts(iris$Species) <- contr.poly(3)

#perform 
list2 <- lapply(colnames(iris)[1:3], function(x){
anova_fit = aov(reformulate(c("Petal.Width","Species"),x), data = iris)
summary(effect("Species",anova_fit, se=TRUE))
})

The above code worked well with the help of @StupidWolf after I raised the former question (How to extract marginal means of multiple variables with effect() function?). And then I got error when I perform the following code:

means.all <- stack(lapply(colnames(iris)[1:3], function(x){
anova_fit = aov(reformulate(c("Petal.Width","Species"),x), data = iris)
summary(effect("Species",anova_fit, se=TRUE))[[5]][1]
}))[2:1]

The error is Error in rep.int(factor(names(x), unique(names(x))), lengths(x)) : invalid 'times' value.

But when I extract marginal mean on only one outcome variable (take Sepal.Length as an example), I could extract the marginal mean with the code below:

anova_fit = aov(reformulate(c("Petal.Width","Species"),"Sepal.Length"), data = iris)
means1 <- summary(effect("Species",anova_fit, se=TRUE))[[5]][1]

I do not know how to correctly use both the stack() and lapply() to extract marginal means.

Many thanks!

Ella

Solution

  • I am not sure how you want your final expected output to look like.

    Probably, you can try this approach :

    do.call(rbind, lapply(list2, function(x) 
  data.frame(prop = c('effect', 'lower', 'upper'), 
            rbind(x$effect, x$lower, x$upper))))


#    prop setosa versicolor virginica
#1 effect   5.88       5.82      5.83
#2  lower   5.49       5.68      5.49
#3  upper   6.27       5.96      6.17
#4 effect   4.17       2.67      2.33
#5  lower   3.93       2.58      2.11
#6  upper   4.42       2.76      2.54
#7 effect   2.43       4.13      4.71
#8  lower   2.13       4.02      4.44
#9  upper   2.74       4.24      4.98

    You can also simplify this by replacing do.call + rbind with purrr's map_df :

    purrr::map_df(list2, function(x) data.frame(prop = c('effect', 'lower', 'upper'), 
                                 rbind(x$effect, x$lower, x$upper)))