Práctica

library(readr)
## Warning: package 'readr' was built under R version 4.5.2
library(vegan)
## Warning: package 'vegan' was built under R version 4.5.2
## Warning: package 'permute' was built under R version 4.5.2
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.5.2
# Matriz biológica: abundancia porcentual de líquenes
Matriz_liquenes <- read_csv("~/Ecología del Paisaje/Semana 2/Matriz_liquenes.csv")

# Matriz ambiental
Matriz_amb_liq <- read_csv("~/Ecología del Paisaje/Semana 2/Matriz_amb_liq.csv")
head(Matriz_liquenes)
## # A tibble: 6 × 21
##   Sitio   Sp1   Sp2   Sp3   Sp4   Sp5   Sp6   Sp7   Sp8   Sp9  Sp10  Sp11  Sp12
##   <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 S1        5     0     8     2     6     1     9     4     7     3    10     2
## 2 S2        3     1     7     3     5     0     8     3     6     2     9     3
## 3 S3        6     2     9     1     7     2     7     5     8     4    11     1
## 4 S4        2     1     6     2     4     1     8     4     7     3     8     2
## 5 S5        7     3     8     3     8     3     9     6     9     5    10     3
## 6 S6        8     2     9     4     7     2     7     5     8     4    12     2
## # ℹ 8 more variables: Sp13 <dbl>, Sp14 <dbl>, Sp15 <dbl>, Sp16 <dbl>,
## #   Sp17 <dbl>, Sp18 <dbl>, Sp19 <dbl>, Sp20 <dbl>
head(Matriz_amb_liq)
## # A tibble: 6 × 16
##   Sitio  Temp HumRel Radiacion Elevacion CobDosel Abertura  NDVI BiomHerb
##   <chr> <dbl>  <dbl>     <dbl>     <dbl>    <dbl>    <dbl> <dbl>    <dbl>
## 1 S1     18.3     72       105      1870       85       12  0.71      230
## 2 S2     17.8     74       110      1890       80       15  0.74      220
## 3 S3     18.5     70       120      1905       78       18  0.78      240
## 4 S4     19.1     68       125      1920       82       10  0.7       210
## 5 S5     17.6     75       130      1950       76       20  0.82      260
## 6 S6     18       71       118      1885       79       14  0.72      240
## # ℹ 7 more variables: DensArboles <dbl>, DensArbustos <dbl>, AreaBasal <dbl>,
## #   Rugosidad <dbl>, DistAgua <dbl>, Ruido <dbl>, TempCorteza <dbl>

Transformación de los datos

library(vegan)
library(dplyr)
Matriz_liquenes <- Matriz_liquenes[,-1] #Elimino la primera columna
Matriz_amb_liq <- Matriz_amb_liq[,-1] #Elimino la primera columna

Matriz_liquenes_hell <- decostand(Matriz_liquenes, "total") #Hellinger
Matriz_amb_liq <- decostand(Matriz_amb_liq, "standardize")

Modelo

modelo_rda<- rda(Matriz_liquenes_hell ~., data=Matriz_amb_liq) # El punto luego de la virgulilla significa que está usando todas las columnas de la Matriz ambiental  
r2adj<-RsquareAdj(modelo_rda)$adj.r.squared 
anova.cca(modelo_rda, step= 1000)
## Permutation test for rda under reduced model
## Permutation: free
## Number of permutations: 999
## 
## Model: rda(formula = Matriz_liquenes_hell ~ Temp + HumRel + Radiacion + Elevacion + CobDosel + Abertura + NDVI + BiomHerb + DensArboles + DensArbustos + AreaBasal + Rugosidad + DistAgua + Ruido + TempCorteza, data = Matriz_amb_liq)
##          Df   Variance      F Pr(>F)
## Model    15 0.00067109 1.0255  0.478
## Residual 14 0.00061077

Afinación del modelo

mod.sel <- ordiR2step (rda(Matriz_liquenes_hell ~ 1, data=Matriz_amb_liq), # Corre el modelo sin utilizar columnas de la ambiental                          
           rda(Matriz_liquenes_hell ~ ., data=Matriz_amb_liq), # Corre el modelo utilizando todas las columnas                            
           direction = "both")  # Va añadiendo parámetros uno por uno y mide la significancia del modelo en cada paso (forward regression). Luego, agarra todos los parámetros y los va quitando uno por uno (backward regression).
## Step: R2.adj= 0 
## Call: Matriz_liquenes_hell ~ 1 
##  
##                   R2.adjusted
## + HumRel         0.0372131522
## <All variables>  0.0130231810
## + AreaBasal      0.0051961406
## + TempCorteza    0.0016879330
## + Temp           0.0010094345
## <model>          0.0000000000
## + Rugosidad     -0.0003531162
## + DensArbustos  -0.0056694274
## + BiomHerb      -0.0076732947
## + DensArboles   -0.0085270793
## + Abertura      -0.0092063597
## + CobDosel      -0.0113869066
## + Elevacion     -0.0113945506
## + NDVI          -0.0140102544
## + Radiacion     -0.0179338310
## + Ruido         -0.0182783116
## + DistAgua      -0.0183333021

Ver la significancia de cada variable (p-value) para el modelo.

  • Valores positivos significa que aporta al modelo

  • Valores negativos significan que deben desecharse porque no empeoran el modelo.

set.seed(1234)
rda1<- rda(Matriz_liquenes_hell ~1, data=Matriz_amb_liq) 
rda2<- rda(Matriz_liquenes_hell ~., data=Matriz_amb_liq) 
ordistep(rda1, scop= formula(rda2), direction= "forward", pstep=1000)
## 
## Start: Matriz_liquenes_hell ~ 1 
## 
##                Df     AIC      F Pr(>F)  
## + HumRel        1 -198.99 2.1209  0.020 *
## + AreaBasal     1 -198.01 1.1515  0.245  
## + TempCorteza   1 -197.90 1.0490  0.355  
## + Temp          1 -197.88 1.0293  0.385  
## + Rugosidad     1 -197.84 0.9898  0.420  
## + DensArbustos  1 -197.68 0.8365  0.560  
## + BiomHerb      1 -197.62 0.7792  0.640  
## + DensArboles   1 -197.60 0.7548  0.660  
## + Abertura      1 -197.58 0.7355  0.695  
## + CobDosel      1 -197.51 0.6735  0.700  
## + Elevacion     1 -197.51 0.6733  0.750  
## + NDVI          1 -197.44 0.5993  0.810  
## + Radiacion     1 -197.32 0.4891  0.895  
## + Ruido         1 -197.31 0.4794  0.940  
## + DistAgua      1 -197.31 0.4779  0.950  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Step: Matriz_liquenes_hell ~ HumRel 
## 
##                Df     AIC      F Pr(>F)
## + TempCorteza   1 -198.63 1.5130  0.120
## + Temp          1 -198.65 1.5329  0.155
## + AreaBasal     1 -198.37 1.2736  0.230
## + Rugosidad     1 -198.11 1.0232  0.395
## + DensArboles   1 -197.91 0.8446  0.560
## + DensArbustos  1 -197.89 0.8241  0.580
## + BiomHerb      1 -197.91 0.8373  0.585
## + Abertura      1 -197.89 0.8223  0.625
## + NDVI          1 -197.75 0.6957  0.715
## + CobDosel      1 -197.74 0.6822  0.755
## + Elevacion     1 -197.61 0.5669  0.840
## + Ruido         1 -197.60 0.5533  0.880
## + Radiacion     1 -197.55 0.5079  0.890
## + DistAgua      1 -197.51 0.4706  0.895
## 
## Call: rda(formula = Matriz_liquenes_hell ~ HumRel, data = Matriz_amb_liq)
## 
##                 Inertia Proportion Rank
## Total         1.282e-03  1.000e+00     
## Constrained   9.026e-05  7.041e-02    1
## Unconstrained 1.192e-03  9.296e-01   19
## 
## Inertia is variance
## 
## Eigenvalues for constrained axes:
##      RDA1 
## 9.026e-05 
## 
## Eigenvalues for unconstrained axes:
##        PC1        PC2        PC3        PC4        PC5        PC6        PC7 
## 0.00027355 0.00025831 0.00015498 0.00013120 0.00007835 0.00007145 0.00005446 
##        PC8 
## 0.00003879 
## (Showing 8 of 19 unconstrained eigenvalues)

Visualización

plot(modelo_rda, type = "n", scaling = 2)

# Sitios
points(modelo_rda, display = "sites", pch = 21, bg = "gray90", col = "black")

# Especies
text(modelo_rda, display = "species", col = "blue", cex = 0.8)

# Variables ambientales
bp <- scores(modelo_rda, display = "bp", scaling = 2)

arrows(0, 0, bp[,1], bp[,2], col = "red", length = 0.1)
text(bp[,1], bp[,2], labels = rownames(bp), col = "red", pos = 4)