[R-sig-eco] mvabund-interprating the uni.test

Thu Oct 23 00:07:38 CEST 2014

I'm using mvabund uni.test to select species that are changing according to a treatment gradient. I've pulled out species that have a significantly different intercept from the null but I'd also like to indicate if they are positive or negative with treatment and extract the "sum-of-LR" as Warton does in the paper showing the usefulness of this approach (http://onlinelibrary.wiley.com/doi/10.1111/j.2041-210X.2011.00127.x/abstract).  I've run manyglm followed by anova on the spider data to try to match the mvabund documentation but I can't figure out where the results beyond the p-values are stored, here are the str for the manyglm and anova output.  (the dput are huge)

thanks

Kendra

> str(glm.spid)
List of 41
 $ coefficients       : num [1:3, 1:12] 2.2353 0.0582 -1.3425 2.1552 -0.4229 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ var.coefficients   : num [1:3, 1:12] 0.1068 0.0185 0.1351 0.252 0.0466 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ fitted.values      : num [1:28, 1:12] 9.349 0.843 9.349 9.349 9.349 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ linear.predictor   : num [1:28, 1:12] 2.24 -0.17 2.24 2.24 2.24 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ residuals          : num [1:28, 1:12] 1.767 -0.778 0.638 -0.83 -0.943 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ PIT.residuals      : num [1:28, 1:12] 0.934 0.534 0.79 0.238 0.101 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ sqrt.1_Hii         : num [1:28, 1:12] 0.945 0.89 0.945 0.945 0.945 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ var.estimator      : num [1:28, 1:12] 87.79 1.48 87.79 87.79 87.79 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ sqrt.weight        : num [1:28, 1:12] 0.998 0.693 0.998 0.998 0.998 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ theta              : Named num [1:12] 1.114 0.428 1.301 0.156 0.715 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ two.loglike        : Named num [1:12] -121.5 -142.1 -78.9 -60.6 -44.5 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ deviance           : Named num [1:12] 22.2 29.4 18.7 16 10.5 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ aic                : Named num [1:12] 129.5 150.1 86.9 68.6 52.5 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ iter               : Named num [1:12] 3 3 6 4 5 4 4 5 4 4 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ data               :'data.frame':	28 obs. of  3 variables:
  ..$ spiddat      : int [1:28, 1:12] 25 0 15 2 1 0 2 0 1 3 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
  .. ..- attr(*, "class")= chr [1:2] "mvabund" "matrix"
  ..$ bare.sand    : num [1:28] 0 0 0 0 0 ...
  ..$ fallen.leaves: num [1:28] 0 1.79 0 0 0 ...
  ..- attr(*, "terms")=Classes 'terms', 'formula' length 3 spiddat ~ bare.sand + fallen.leaves
  .. .. ..- attr(*, "variables")= language list(spiddat, bare.sand, fallen.leaves)
  .. .. ..- attr(*, "factors")= int [1:3, 1:2] 0 1 0 0 0 1
  .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. ..$ : chr [1:3] "spiddat" "bare.sand" "fallen.leaves"
  .. .. .. .. ..$ : chr [1:2] "bare.sand" "fallen.leaves"
  .. .. ..- attr(*, "term.labels")= chr [1:2] "bare.sand" "fallen.leaves"
  .. .. ..- attr(*, "order")= int [1:2] 1 1
  .. .. ..- attr(*, "intercept")= int 1
  .. .. ..- attr(*, "response")= int 1
  .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. .. ..- attr(*, "predvars")= language list(spiddat, bare.sand, fallen.leaves)
  .. .. ..- attr(*, "dataClasses")= Named chr [1:3] "nmatrix.12" "numeric" "numeric"
  .. .. .. ..- attr(*, "names")= chr [1:3] "spiddat" "bare.sand" "fallen.leaves"
 $ stderr.coefficients: num [1:3, 1:12] 0.327 0.136 0.368 0.502 0.216 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ phi                : Named num [1:12] 0.897 2.338 0.768 6.43 1.399 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ tol                : num 1.49e-08
 $ maxiter            : num 25
 $ maxiter2           : num 10
 $ AIC                : Named num [1:12] 129.5 150.1 86.9 68.6 52.5 ...
  ..- attr(*, "names")= chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ AICsum             : num 1658
 $ family             : chr "negative.binomial"
 $ K                  : num 1
 $ theta.method       : chr "PHI"
 $ cor.type           : chr "I"
 $ shrink.param       : num 0
 $ call               : language manyglm(formula = spiddat ~ bare.sand + fallen.leaves, family = "negative.binomial", data = X)
 $ terms              :Classes 'terms', 'formula' length 3 spiddat ~ bare.sand + fallen.leaves
  .. ..- attr(*, "variables")= language list(spiddat, bare.sand, fallen.leaves)
  .. ..- attr(*, "factors")= int [1:3, 1:2] 0 1 0 0 0 1
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:3] "spiddat" "bare.sand" "fallen.leaves"
  .. .. .. ..$ : chr [1:2] "bare.sand" "fallen.leaves"
  .. ..- attr(*, "term.labels")= chr [1:2] "bare.sand" "fallen.leaves"
  .. ..- attr(*, "order")= int [1:2] 1 1
  .. ..- attr(*, "intercept")= int 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. ..- attr(*, "predvars")= language list(spiddat, bare.sand, fallen.leaves)
  .. ..- attr(*, "dataClasses")= Named chr [1:3] "nmatrix.12" "numeric" "numeric"
  .. .. ..- attr(*, "names")= chr [1:3] "spiddat" "bare.sand" "fallen.leaves"
 $ assign             : int [1:3] 0 1 2
 $ formula            :Class 'formula' length 3 spiddat ~ bare.sand + fallen.leaves
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
 $ rank               : int 3
 $ xlevels            : Named list()
 $ df.residual        : int 25
 $ x                  : num [1:28, 1:3] 1 1 1 1 1 1 1 1 1 1 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
 $ y                  : num [1:28, 1:12] 25 0 15 2 1 0 2 0 1 3 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ qr                 :List of 4
  ..$ qr   : num [1:28, 1:3] -5.292 0.189 0.189 0.189 0.189 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:28] "1" "2" "3" "4" ...
  .. .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  ..$ rank : int 3
  ..$ qraux: num [1:3] 1.19 1.11 1.18
  ..$ pivot: int [1:3] 1 2 3
  ..- attr(*, "class")= chr "qr"
 $ show.coef          : logi FALSE
 $ show.fitted        : logi FALSE
 $ show.residuals     : logi FALSE
 $ offset             : num [1:28, 1:12] 0 0 0 0 0 0 0 0 0 0 ...
 - attr(*, "class")= chr [1:2] "manyglm" "mglm"
> str(an.spid)
List of 11
 $ family      : chr "negative.binomial"
 $ p.uni       : chr "adjusted"
 $ test        : chr "Dev"
 $ cor.type    : chr "I"
 $ resamp      : chr "pit.trap"
 $ nBoot       : num 1000
 $ shrink.param: num [1:3] 0 0 0
 $ n.bootsdone : num 999
 $ table       :'data.frame':	3 obs. of  4 variables:
  ..$ Res.Df  : int [1:3] 27 26 25
  ..$ Df.diff : num [1:3] NA 1 1
  ..$ Dev     : num [1:3] NA 70.4 70.3
  ..$ Pr(>Dev): num [1:3] NA 0.002 0.007
  ..- attr(*, "heading")= chr [1:3] "Analysis of Deviance Table\n" "Model: manyglm(formula = spiddat ~ bare.sand + fallen.leaves, family = \"negative.binomial\", " "Model:     data = X)"
  ..- attr(*, "title")= chr "\nMultivariate test:\n"
 $ uni.p       : num [1:3, 1:12] NA 0.568 0.001 NA 0.485 0.988 NA 0.001 0.074 NA ...
  ..- attr(*, "title")= chr "\nUnivariate Tests:\n"
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 $ uni.test    : num [1:3, 1:12] NA 1.22 27.64 NA 2.42 ...
  ..- attr(*, "title")= chr "\nUnivariate Tests:\n"
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "(Intercept)" "bare.sand" "fallen.leaves"
  .. ..$ : chr [1:12] "Alopacce" "Alopcune" "Alopfabr" "Arctlute" ...
 - attr(*, "class")= chr "anova.manyglm"

--
Kendra Maas, Ph.D.
Post Doctoral Research Fellow
University of British Columbia
604-822-5646