setwd("~/Dropbox/Proyecto_posdoc_PUJ/Proyecto_Amazonia/data")
abu=read.table("means.txt", header=T)
attach(abu)
library(mvabund)
gfit1 = manyglm (abun.abs ~ log(hmax)+log(ldmc)+log(lth)+log(wd)+log(lma)+paisaje+paisaje:log(hmax)+paisaje:log(ldmc)+paisaje:log(lth)+paisaje:log(wd)+paisaje:log(lma), family="negative.binomial")
anova(gfit1)
## Time elapsed: 0 hr 0 min 22 sec
## Analysis of Deviance Table
##
## Model: manyglm(formula = abun.abs ~ log(hmax) + log(ldmc) + log(lth) +
## Model: log(wd) + log(lma) + paisaje + paisaje:log(hmax) + paisaje:log(ldmc) +
## Model: paisaje:log(lth) + paisaje:log(wd) + paisaje:log(lma), family = "negative.binomial")
##
## Multivariate test:
## Res.Df Df.diff Dev Pr(>Dev)
## (Intercept) 100
## log(hmax) 99 1 18.303 0.001 ***
## log(ldmc) 98 1 7.538 0.013 *
## log(lth) 97 1 0.487 0.516
## log(wd) 96 1 7.659 0.008 **
## log(lma) 95 1 0.876 0.395
## paisaje 94 1 5.533 0.037 *
## log(hmax):paisaje 93 1 5.715 0.037 *
## log(ldmc):paisaje 92 1 0.068 0.807
## log(lth):paisaje 91 1 1.165 0.361
## log(wd):paisaje 90 1 3.140 0.139
## log(lma):paisaje 89 1 3.491 0.094 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Arguments: P-value calculated using 999 resampling iterations via PIT-trap resampling (to account for correlation in testing).
coef(gfit1)
## abun.abs
## (Intercept) 0.410249525
## log(hmax) 0.751075492
## log(ldmc) -0.307379076
## log(lth) -0.553138665
## log(wd) 0.206927227
## log(lma) -0.739496190
## paisajeMontana 16.290685556
## log(hmax):paisajeMontana -1.807780041
## log(ldmc):paisajeMontana -1.463188133
## log(lth):paisajeMontana 0.001696549
## log(wd):paisajeMontana -0.885885961
## log(lma):paisajeMontana 1.668397160
par(mfrow=c(2,1),mgp=c(1.75,0.75,0),mar=c(3,3,2,1))
res=residuals(gfit1)
plot(res~fitted(gfit1), log="x", xlab="Fitted values [log scale]", ylab = "Dunn-Smyth residuals")
abline(h=0, col="red")
qqnorm(res,main = "")
abline(c(0,1), col="red")
library("dplyr", lib.loc="/Library/Frameworks/R.framework/Versions/3.5/Resources/library")
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
abu.lom <- filter(abu, paisaje == "Lomerio")
detach(abu)
attach(abu.lom)
gfit2 = manyglm (abun.abs ~ log(hmax)+log(ldmc)+log(lth)+log(wd)+log(lma), family="negative.binomial")
anova(gfit2)
## Time elapsed: 0 hr 0 min 3 sec
## Analysis of Deviance Table
##
## Model: manyglm(formula = abun.abs ~ log(hmax) + log(ldmc) + log(lth) +
## Model: log(wd) + log(lma), family = "negative.binomial")
##
## Multivariate test:
## Res.Df Df.diff Dev Pr(>Dev)
## (Intercept) 34
## log(hmax) 33 1 0.072 0.793
## log(ldmc) 32 1 4.885 0.050 *
## log(lth) 31 1 2.690 0.172
## log(wd) 30 1 0.714 0.478
## log(lma) 29 1 1.957 0.235
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Arguments: P-value calculated using 999 resampling iterations via PIT-trap resampling (to account for correlation in testing).
coef(gfit2)
## abun.abs
## (Intercept) 0.3779867
## log(hmax) 0.7330626
## log(ldmc) -0.2971001
## log(lth) -0.5517289
## log(wd) 0.1952171
## log(lma) -0.7443197
par(mfrow=c(2,1),mgp=c(1.75,0.75,0),mar=c(3,3,2,1))
res=residuals(gfit2)
plot(res~fitted(gfit2), log="x", xlab="Fitted values [log scale]", ylab = "Dunn-Smyth residuals")
abline(h=0, col="red")
qqnorm(res,main = "")
abline(c(0,1), col="red")
abu.mont <- filter(abu, paisaje == "Montana")
detach(abu.lom)
attach(abu.mont)
gfit3 = manyglm (abun.abs ~ log(hmax)+log(ldmc)+log(lth)+log(wd)+log(lma), family="negative.binomial")
anova(gfit3)
## Time elapsed: 0 hr 0 min 5 sec
## Analysis of Deviance Table
##
## Model: manyglm(formula = abun.abs ~ log(hmax) + log(ldmc) + log(lth) +
## Model: log(wd) + log(lma), family = "negative.binomial")
##
## Multivariate test:
## Res.Df Df.diff Dev Pr(>Dev)
## (Intercept) 65
## log(hmax) 64 1 19.430 0.001 ***
## log(ldmc) 63 1 8.247 0.010 **
## log(lth) 62 1 0.103 0.769
## log(wd) 61 1 6.521 0.018 *
## log(lma) 60 1 1.737 0.220
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Arguments: P-value calculated using 999 resampling iterations via PIT-trap resampling (to account for correlation in testing).
coef(gfit3)
## abun.abs
## (Intercept) 16.5757835
## log(hmax) -1.0493544
## log(ldmc) -1.7572646
## log(lth) -0.5470706
## log(wd) -0.6830834
## log(lma) 0.9162533
par(mfrow=c(2,1),mgp=c(1.75,0.75,0),mar=c(3,3,2,1))
res=residuals(gfit3)
plot(res~fitted(gfit3), log="x", xlab="Fitted values [log scale]", ylab = "Dunn-Smyth residuals")
abline(h=0, col="red")
qqnorm(res,main = "")
abline(c(0,1), col="red")
