Hipótesis: Están las frecuencuencias Dominante de F12 Khz explicadas por la estructura vertical, horizontal y área basal.
Cargamos paquetes:
library(lme4)
library(nlme)
library(MCMCglmm)
library(AICcmodavg)
library(mgcv)
library(MuMIn)
library(MASS)
library(lattice)
library(lavaan)
library(piecewiseSEM)
library(car)
Cargamos la base de datos:
load("F12.Rda")
Estructura de variables:
str(F12)
## 'data.frame': 34 obs. of 7 variables:
## $ Site : Factor w/ 36 levels "CE01","CE02",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ mean : num 1310 189 476 363 683 ...
## $ sd : num 1102 436 797 636 731 ...
## $ Frec_Dom : chr "F12" "F12" "F12" "F12" ...
## $ prom.vert: num 0.1421 0.091 0.1065 0.0828 0.0906 ...
## $ prom.hor : num 0.0542 0.0601 0.0596 0.0586 0.0583 ...
## $ AB : num 20296 10729 11029 7827 15654 ...
## - attr(*, ".internal.selfref")=<externalptr>
## - attr(*, "sorted")= chr "Site"
Estandarización de variables:
Site <- as.factor(F12$Site)
FDm <- as.numeric(F12$mean)
FD <- as.factor(F12$Frec_Dom)
Fv <- as.numeric(((F12$prom.vert)-mean(F12$prom.vert))/sd(F12$prom.vert))
Fh <- as.numeric(((F12$prom.hor)-mean(F12$prom.hor))/sd(F12$prom.hor))
AB <- as.numeric(((F12$AB)-mean(F12$AB))/sd(F12$AB))
Agrupación de variables transformadas:
data<- data.frame(Site,FDm, FD, Fv, Fh,AB )
attach(data)
## The following objects are masked _by_ .GlobalEnv:
##
## AB, FD, FDm, Fh, Fv, Site
mod1 <- lm(FDm~Fv+Fh+AB,data = F12)
mod2 <- lm(FDm~Fh,data = F12)
mod3 <- lm(FDm~Fv,data = F12)
mod4 <- lm(FDm~Fv+Fh+Fv:Fh ,data = F12)
mod5 <- lm(FDm~AB)
mod0 <- lm(FDm~1,data = F12)
AIC(mod1, mod2, mod3, mod4, mod5)
## df AIC
## mod1 5 467.8711
## mod2 3 475.6497
## mod3 3 476.1118
## mod4 5 479.0901
## mod5 3 468.0821
AIC(mod1, mod0)
## df AIC
## mod1 5 467.8711
## mod0 2 474.2632
No hay fuertes evidencias de que existe un mejor modelo.
#Comprobación de Supuestos
par(mfrow=c(2,2))
plot(mod2)
#Revizamos la ortogonalidad de las variables independientes
cor_PCA <- (princomp(F12[,5:6]))
biplot(cor_PCA)
#Seleccion de variables a traves del método “backward”
step <- stepAIC(mod1, direction="backward")
## Start: AIC=369.38
## FDm ~ Fv + Fh + AB
##
## Df Sum of Sq RSS AIC
## - Fh 1 74058 1478685 369.13
## <none> 1404626 369.38
## - Fv 1 96015 1500641 369.63
## - AB 1 576525 1981151 379.08
##
## Step: AIC=369.13
## FDm ~ Fv + AB
##
## Df Sum of Sq RSS AIC
## <none> 1478685 369.13
## - Fv 1 111140 1589825 369.59
## - AB 1 534648 2013332 377.62
No se muetsra una efecto importante para eliminar alguna variable
#Para la importancia relativa de cada variable del modelo
library(leaps)
attach(data)
## The following objects are masked _by_ .GlobalEnv:
##
## AB, FD, FDm, Fh, Fv, Site
## The following objects are masked from data (pos = 4):
##
## AB, FD, FDm, Fh, Fv, Site
leaps <- regsubsets(FDm~Fv+Fh,data = data, nbest = 10)
summary(leaps)
## Subset selection object
## Call: regsubsets.formula(FDm ~ Fv + Fh, data = data, nbest = 10)
## 2 Variables (and intercept)
## Forced in Forced out
## Fv FALSE FALSE
## Fh FALSE FALSE
## 10 subsets of each size up to 2
## Selection Algorithm: exhaustive
## Fv Fh
## 1 ( 1 ) " " "*"
## 1 ( 2 ) "*" " "
## 2 ( 1 ) "*" "*"
plot(leaps, scale = "r2")
#intercepto y Fh explica 0.0056 de la varianza
#Autocorrelación
vif(mod1) #Si es menor a 2-4 esta bien
## Fv Fh AB
## 1.133011 1.032976 1.142911
sqrt(vif(mod1)) #Muetras la inflacion de de cada variable (entre 2 a 4)
## Fv Fh AB
## 1.064430 1.016354 1.069070
max(vif(mod1)) #Si es menor a 2-4 esta bien
## [1] 1.142911
#Relación entre variables
par(mfrow=c(1,3))
library(visreg)
visreg(mod1)
Modelo lineal generalizado:
glm1 <- glm((FDm+1)~Fv+Fh+AB, family=poisson, data=data)
summary(glm1)
##
## Call:
## glm(formula = (FDm + 1) ~ Fv + Fh + AB, family = poisson, data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -22.8911 -9.0653 -0.8553 5.2278 22.4665
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.845081 0.009441 619.12 <2e-16 ***
## Fv 0.171149 0.010058 17.02 <2e-16 ***
## Fh 0.136155 0.008813 15.45 <2e-16 ***
## AB 0.326198 0.008185 39.85 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 4999.7 on 33 degrees of freedom
## Residual deviance: 3471.3 on 30 degrees of freedom
## AIC: Inf
##
## Number of Fisher Scoring iterations: 4
glm1quasi<- glm((FDm+1)~Fv+Fh+AB, family=quasipoisson, data) #Modelo de quasi, por la sobredispersión
summary(glm1quasi)
##
## Call:
## glm(formula = (FDm + 1) ~ Fv + Fh + AB, family = quasipoisson,
## data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -22.8911 -9.0653 -0.8553 5.2278 22.4665
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.84508 0.09859 59.285 < 2e-16 ***
## Fv 0.17115 0.10504 1.629 0.113692
## Fh 0.13616 0.09204 1.479 0.149481
## AB 0.32620 0.08548 3.816 0.000631 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for quasipoisson family taken to be 109.0583)
##
## Null deviance: 4999.7 on 33 degrees of freedom
## Residual deviance: 3471.3 on 30 degrees of freedom
## AIC: NA
##
## Number of Fisher Scoring iterations: 4