© Ricardo Solar/UFMG - compartilhamento e utilização não-comercial livres. Não reproduzir sem autorização
rm(list=ls())
dados<-read.table("Formigas.txt",header=T)
summary(dados)
## Canopy_cover SolInv CreTen
## Min. :0.0000 Min. :0.00 Min. :0.0000
## 1st Qu.:0.7999 1st Qu.:0.00 1st Qu.:0.0000
## Median :0.8468 Median :0.00 Median :0.0000
## Mean :0.8054 Mean :0.25 Mean :0.3793
## 3rd Qu.:0.8878 3rd Qu.:0.25 3rd Qu.:1.0000
## Max. :0.9520 Max. :1.00 Max. :1.0000
plot(dados)
names(dados)
## [1] "Canopy_cover" "SolInv" "CreTen"
Solenopsis invicta (SolInv) Crematogaster tenuicola (CreTen)
Dados de presenca/ausencia (Bernoulli, binarios) : funcao de ligacao ainda eh logit = log (p/1-p)
m.completoSI = glm(SolInv~Canopy_cover, family=binomial, data=dados)
anova(m.completoSI,test="Chisq")
summary(m.completoSI)
##
## Call:
## glm(formula = SolInv ~ Canopy_cover, family = binomial, data = dados)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.8229 -0.6465 -0.5603 -0.2981 2.1202
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.636 1.520 3.051 0.002281 **
## Canopy_cover -7.253 1.892 -3.834 0.000126 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 130.46 on 115 degrees of freedom
## Residual deviance: 109.30 on 114 degrees of freedom
## AIC: 113.3
##
## Number of Fisher Scoring iterations: 4
with(m.completoSI, deviance/df.residual)
## [1] 0.9587829
##Dispersao OK
library(RT4Bio)
## Loading required package: MASS
## Loading required package: survival
rdiagnostic(m.completoSI)
DHARMa::simulateResiduals(m.completoSI, plot=T)
## Object of Class DHARMa with simulated residuals based on 250 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help.
##
## Scaled residual values: 0.1474335 0.8900209 0.9375736 0.4623027 0.4761113 0.4115664 0.9587229 0.08897467 0.4367838 0.5749464 0.7218649 0.352572 0.3414428 0.1908792 0.7419757 0.9394594 0.5293834 0.07718971 0.3055862 0.1929606 ...
m.completoCT = glm(CreTen~Canopy_cover, family=binomial, data=dados)
anova(m.completoCT,test="Chisq")
with(m.completoCT, deviance/df.residual)
## [1] 1.194863
##Dispersao OK
rdiagnostic(m.completoCT)
DHARMa::simulateResiduals(m.completoCT, plot=T)
## Object of Class DHARMa with simulated residuals based on 250 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help.
##
## Scaled residual values: 0.3997087 0.8388238 0.4180393 0.4868932 0.7773231 0.7790711 0.3279211 0.06762075 0.7504804 0.3531814 0.488426 0.7947877 0.1064237 0.05854699 0.9313812 0.743161 0.8135357 0.6945322 0.2143449 0.1376952 ...
Uma forma um pouco diferente de representar os resultados, vamos plotar ambas as especies no mesmo grafico
SI <- coef(m.completoSI)
CT <- coef(m.completoCT)
plot(SolInv~Canopy_cover, data=dados, axes=F, type="n", xlab="", ylab="Ocorrencia da especie")
axis(1, lwd=2)
axis(2,lwd=2, las=1, at=c(0,0.5,1))
box(bty="l", lwd=3)
##Pontos e curve para Solenopsis invicta
points(SolInv~Canopy_cover, data=dados, pch=21, bg=2)
points(SolInv~Canopy_cover, data=dados, pch=21, bg="#CC000050")
curve(plogis(SI[1]+SI[2]*x), add=T, lwd=2, col=2)
##Pontos e curve para Crematogaster tenuicola
##Vou dar uma melhorada aqui, e deslocar de leve os pontos de CreTen pra poder visualizar melhor
CreTen2 <- ifelse(dados$CreTen==1, dados$CreTen-0.03, dados$CreTen+0.03)
points(CreTen2~Canopy_cover, data=dados, pch=21, bg=1)
curve(plogis(CT[1]+CT[2]*x), add=T, lwd=2, col=1)
legend(0,.3, c("Cre ten", "Sol inv"), col=1:2, lty=1, lwd=2, pch=21, pt.bg = 1:2, bty="n")