[R-sig-eco] plyr and mvabund, conceptual issue

[Eduard Szöcs]

Hai Kendra,

i've used a simple for-loop to do this in the past.

Something along these lines:


###-----------------------------------------------------------------
mymv <- function(response, env, zone) {
  df <- data.frame(env = env, zone = zone)
  out <- NULL
  for (i in levels(zone)) {
    rsp <- mvabund(response[zone == i, ])
    out[[i]]$mod <- manyglm(rsp ~ env, data = df[zone == i, ])
    out[[i]]$anova <- anova(out[[i]]$mod, p.uni = "adjusted", resamp =
"perm.resid", nBoot = 10)
  }
  return(out)
}

# try it
require(mvabund)
# some data
env <- runif(100)
zone <- gl(2, 50)
response <- data.frame(A = rpois(100, 5), B = rpois(100, 1))
# run it
mymv(response, env, zone)
###-----------------------------------------------------------------

Hope this helps,

Eduard


-- 
Eduard Szöcs
Quantitative Landscape Ecology
Institute for Environmental Sciences
University Koblenz-Landau
Fortstrasse 7
76829 Landau
Germany
http://www.uni-koblenz-landau.de/campus-landau/faculty7/environmental-sciences/landscape-ecology/Staff/eduardszoecs


On 29/10/14 00:31, Maas, Kendra wrote:
> I'm trying to run mvabund (generate glms for each species and do univariate anova to determine "indicator species" that respond to my treatments) on a lot of subsets of my data.  I'm having theoretical difficulty with how to use plyr on multiple dataframes or lists and outputting lists.  Previously I've run this series of commands using text editor to change the selected zone, I know that this is what plyr is designed for but I'm getting stuck
> 
> "b.red.otu" is a sample by species dataframe, "b.env" is a sample by factor dataframe containing the variables "zone" and "om"
> 
>> b.oJP.mva<-mvabund(subset(b.red.otu,zone=="oJP"))
>> b.oJP.nb<-manyglm(b.oJP.mva~subset(b.env, b.env$zone=="oJP")$om, family="negative.binomial")
>> b.oJP.nb.anova<-anova(b.oJP.nb,p.uni="adjusted", resamp="perm.resid")
> 
> 
> This code works, it's just really ugly and requires a lot of copy and paste/find and replace for every possible zone (I have tens of subsets that I want to look at)
> 
> 
> Here is how I'm attempting to work out my code with the Tasmania data packaged with mvabund.  I convert it to dataframes because I much more comfortable with them.
> 
> 
>> tas.env <- data.frame(Tasmania$treatment, Tasmania$block)
> 
>> tas.abund <- data.frame(Tasmania$abund)
> 
>> mva.out <- dlply(tas.abund, ~tas.env$block, function(x) {
>       mvabund(x~tas.env$treatment) 
>        })
> 
> Which returns an empty list[0]
> 
> I tried creating vectors for treatment and block
> 
>> block <- Tasmania$block
>> treatment <- Tasmania$treatment
> 
>> mva.out <- dlply(tas.abund, ~block, function(x) {
>      mvabund(x~treatment)
>      })
> 
> Error in as.data.frame.default(x[[i]], optional = TRUE) : 
>   cannot coerce class ""formula"" to a data.frame
> 
> I tried llply since mvabund puts out a list and Tasmania is already a list
> 
>> mva.out <- llply(Tasmania$abund, ~Tasmania$block, function(x) {
>       mvabund(x~Tasmania$treatment)
>       })
> 
> Error in llply(Tasmania$abund, ~Tasmania$block, function(x) { : 
>   .fun is not a function.
> 
> 
> I'm sure this is possible to do with plyr, I just can't figure out how.  Suggestions please.
> 
> thanks
> 
> Kendra
> _______________________________________________
> R-sig-ecology mailing list
> R-sig-ecology at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
> .
>