29Abril
Sindy Sandoval
30 4 2021
\(Clase\)
preparando librerias
library(spdep)## Loading required package: sp## Loading required package: spData## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`## Loading required package: sf## Linking to GEOS 3.8.1, GDAL 3.1.4, PROJ 6.3.1library(ape)## Registered S3 method overwritten by 'ape':
## method from
## plot.mst spdeplibrary(sp)
library(MVA)## Loading required package: HSAUR2## Loading required package: toolslibrary(Hmisc)## Loading required package: lattice## Loading required package: survival## Loading required package: Formula## Loading required package: ggplot2##
## Attaching package: 'Hmisc'## The following object is masked from 'package:ape':
##
## zoom## The following objects are masked from 'package:base':
##
## format.pval, unitslibrary(normtest)
library(nortest)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:Hmisc':
##
## src, summarize## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(spatialreg)## Loading required package: Matrix## Registered S3 methods overwritten by 'spatialreg':
## method from
## residuals.stsls spdep
## deviance.stsls spdep
## coef.stsls spdep
## print.stsls spdep
## summary.stsls spdep
## print.summary.stsls spdep
## residuals.gmsar spdep
## deviance.gmsar spdep
## coef.gmsar spdep
## fitted.gmsar spdep
## print.gmsar spdep
## summary.gmsar spdep
## print.summary.gmsar spdep
## print.lagmess spdep
## summary.lagmess spdep
## print.summary.lagmess spdep
## residuals.lagmess spdep
## deviance.lagmess spdep
## coef.lagmess spdep
## fitted.lagmess spdep
## logLik.lagmess spdep
## fitted.SFResult spdep
## print.SFResult spdep
## fitted.ME_res spdep
## print.ME_res spdep
## print.lagImpact spdep
## plot.lagImpact spdep
## summary.lagImpact spdep
## HPDinterval.lagImpact spdep
## print.summary.lagImpact spdep
## print.sarlm spdep
## summary.sarlm spdep
## residuals.sarlm spdep
## deviance.sarlm spdep
## coef.sarlm spdep
## vcov.sarlm spdep
## fitted.sarlm spdep
## logLik.sarlm spdep
## anova.sarlm spdep
## predict.sarlm spdep
## print.summary.sarlm spdep
## print.sarlm.pred spdep
## as.data.frame.sarlm.pred spdep
## residuals.spautolm spdep
## deviance.spautolm spdep
## coef.spautolm spdep
## fitted.spautolm spdep
## print.spautolm spdep
## summary.spautolm spdep
## logLik.spautolm spdep
## print.summary.spautolm spdep
## print.WXImpact spdep
## summary.WXImpact spdep
## print.summary.WXImpact spdep
## predict.SLX spdep##
## Attaching package: 'spatialreg'## The following objects are masked from 'package:spdep':
##
## anova.sarlm, as_dgRMatrix_listw, as_dsCMatrix_I, as_dsCMatrix_IrW,
## as_dsTMatrix_listw, as.spam.listw, bptest.sarlm, can.be.simmed,
## cheb_setup, coef.gmsar, coef.sarlm, coef.spautolm, coef.stsls,
## create_WX, deviance.gmsar, deviance.sarlm, deviance.spautolm,
## deviance.stsls, do_ldet, eigen_pre_setup, eigen_setup, eigenw,
## errorsarlm, fitted.gmsar, fitted.ME_res, fitted.sarlm,
## fitted.SFResult, fitted.spautolm, get.ClusterOption,
## get.coresOption, get.mcOption, get.VerboseOption,
## get.ZeroPolicyOption, GMargminImage, GMerrorsar, griffith_sone,
## gstsls, Hausman.test, HPDinterval.lagImpact, impacts, intImpacts,
## Jacobian_W, jacobianSetup, l_max, lagmess, lagsarlm, lextrB,
## lextrS, lextrW, lmSLX, logLik.sarlm, logLik.spautolm, LR.sarlm,
## LR1.sarlm, LR1.spautolm, LU_prepermutate_setup, LU_setup,
## Matrix_J_setup, Matrix_setup, mcdet_setup, MCMCsamp, ME, mom_calc,
## mom_calc_int2, moments_setup, powerWeights, predict.sarlm,
## predict.SLX, print.gmsar, print.ME_res, print.sarlm,
## print.sarlm.pred, print.SFResult, print.spautolm, print.stsls,
## print.summary.gmsar, print.summary.sarlm, print.summary.spautolm,
## print.summary.stsls, residuals.gmsar, residuals.sarlm,
## residuals.spautolm, residuals.stsls, sacsarlm, SE_classic_setup,
## SE_interp_setup, SE_whichMin_setup, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOption, similar.listw, spam_setup, spam_update_setup,
## SpatialFiltering, spautolm, spBreg_err, spBreg_lag, spBreg_sac,
## stsls, subgraph_eigenw, summary.gmsar, summary.sarlm,
## summary.spautolm, summary.stsls, trW, vcov.sarlm, Wald1.sarlmModelo SACSAR–>Estima Rho y Lambda
Llamando datos
library(readxl)
df= read_excel('/Users/sindyluh/Downloads/Computación estadística/BD_MODELADO.xlsx')Separando información georeferenciada de la otra
df_xy=df[,c(1,2)] # Coords-->cogo de df los datos de columna 1 y 2 (las cuales son coordenadas)
X= df[,-c(1,2)] # Explicativas, resto de df la columna 1 y 2Creando matriz de pesos con una función específica que requiere la libreria para poder usar modelos espaciales
contnb=dnearneigh(coordinates(df_xy),0,380000,longlat = F)#380000 abarca todas las mediciones
contnb## Neighbour list object:
## Number of regions: 313
## Number of nonzero links: 97656
## Percentage nonzero weights: 99.68051
## Average number of links: 312class(contnb)## [1] "nb"df_xy=as.matrix(df_xy)
dlist <- nbdists(contnb, df_xy)
dlist <- lapply(dlist, function(x) 1/x)
Wve=nb2listw(contnb,glist=dlist,style = "W")Corriendo modelo SARAR
library(spdep)
map2= sacsarlm(Avg_CEa_07~1, data= X, listw= Wve)## Warning in sacsarlm(Avg_CEa_07 ~ 1, data = X, listw = Wve): inversion of asymptotic covariance matrix failed for tol.solve = 2.22044604925031e-16
## reciprocal condition number = 5.77785e-20 - using numerical Hessian.summary(map2)##
## Call:sacsarlm(formula = Avg_CEa_07 ~ 1, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.405867 -0.620449 -0.028054 0.640915 2.891325
##
## Type: sac
## Coefficients: (numerical Hessian approximate standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.7267 3.2577 -0.53 0.5961
##
## Rho: 0.97863
## Approximate (numerical Hessian) standard error: 0.021305
## z-value: 45.935, p-value: < 2.22e-16
## Lambda: 0.97863
## Approximate (numerical Hessian) standard error: 0.021299
## z-value: 45.946, p-value: < 2.22e-16
##
## LR test value: 254.56, p-value: < 2.22e-16
##
## Log likelihood: -448.7933 for sac model
## ML residual variance (sigma squared): 0.98538, (sigma: 0.99267)
## Number of observations: 313
## Number of parameters estimated: 4
## AIC: 905.59, (AIC for lm: 1156.2)Comprobando si el modelo map2 es adecuado–>Residuales:normalidad e independencia
#####Normalidad
residuales_map2 =map2$residuals
shapiro.test(residuales_map2)##
## Shapiro-Wilk normality test
##
## data: residuales_map2
## W = 0.99417, p-value = 0.2746#####Independecia
library(ape)
# Matriz de distancias
df.dists <- as.matrix(dist(cbind(df$Avg_X_MCB, df$Avg_Y_MCE)))
# Inversa de las matriz
df.dists.inv <- 1/df.dists
# Asignar ceros a la diagonal
diag(df.dists.inv) <- 0
# Redondear
df.dists.inv <- round(df.dists.inv,3)
# Matriz estandarizada
We<-df.dists.inv/rowSums(df.dists.inv)
# Moran
Moran.I(residuales_map2,We)## $observed
## [1] 0.1041319
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004650578
##
## $p.value
## [1] 0Modelo SLM Spatial Lag Model
library(spdep)
colnames(df)# Mirar el nombre de las demás variables## [1] "Avg_X_MCB" "Avg_Y_MCE" "Avg_CEa_07" "Avg_CEa_15" "NDVI"
## [6] "DEM" "SLOPE" "Avg_z"map3=errorsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map3)##
## Call:errorsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.019160 -0.540466 -0.045367 0.513314 2.592838
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -64.737579 5.752902 -11.2530 < 2.2e-16
## Avg_CEa_15 0.859898 0.083054 10.3535 < 2.2e-16
## NDVI -2.395368 1.907913 -1.2555 0.209301
## DEM 0.036792 0.020974 1.7542 0.079402
## SLOPE -0.073067 0.024760 -2.9510 0.003168
## Avg_z 0.257034 0.028465 9.0299 < 2.2e-16
##
## Lambda: 0.9825, LR test value: 99.359, p-value: < 2.22e-16
## Asymptotic standard error: 0.012342
## z-value: 79.604, p-value: < 2.22e-16
## Wald statistic: 6336.8, p-value: < 2.22e-16
##
## Log likelihood: -406.1005 for error model
## ML residual variance (sigma squared): 0.76603, (sigma: 0.87523)
## Number of observations: 313
## Number of parameters estimated: 8
## AIC: 828.2, (AIC for lm: 925.56)Comprobando si el modelo map3 es adecuado–>Residuales:normalidad e independencia
#####Normalidad
residuales_map3 =map3$residuals
shapiro.test(residuales_map3)##
## Shapiro-Wilk normality test
##
## data: residuales_map3
## W = 0.99166, p-value = 0.07491#####Independencia
Moran.I(residuales_map3,We)## $observed
## [1] 0.1298137
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004653655
##
## $p.value
## [1] 0###Modelo map3 pero descartando NDVI(p value>0,05)
map3b=errorsarlm(formula=Avg_CEa_07~Avg_CEa_15+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map3b)##
## Call:errorsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + DEM + SLOPE +
## Avg_z, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.068942 -0.573110 -0.041672 0.535538 2.620533
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -66.334323 5.621358 -11.8004 < 2.2e-16
## Avg_CEa_15 0.871288 0.082765 10.5273 < 2.2e-16
## DEM 0.039380 0.020925 1.8819 0.059845
## SLOPE -0.074849 0.024782 -3.0203 0.002525
## Avg_z 0.251732 0.028220 8.9203 < 2.2e-16
##
## Lambda: 0.98246, LR test value: 98.998, p-value: < 2.22e-16
## Asymptotic standard error: 0.012369
## z-value: 79.427, p-value: < 2.22e-16
## Wald statistic: 6308.6, p-value: < 2.22e-16
##
## Log likelihood: -406.8867 for error model
## ML residual variance (sigma squared): 0.76989, (sigma: 0.87744)
## Number of observations: 313
## Number of parameters estimated: 7
## AIC: 827.77, (AIC for lm: 924.77)Comprobando si el modelo map3b es adecuado–>Residuales:normalidad e independencia
#####Normalidad
residuales_map3b =map3b$residuals
shapiro.test(residuales_map3b)##
## Shapiro-Wilk normality test
##
## data: residuales_map3b
## W = 0.99235, p-value = 0.1078#####Independencia
Moran.I(residuales_map3b,We)## $observed
## [1] 0.1295797
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004653871
##
## $p.value
## [1] 0Modelo map3 pero descartando NDVI y DEM
map3C=errorsarlm(formula=Avg_CEa_07~Avg_CEa_15+SLOPE+Avg_z,data= X, listw= Wve)
summary(map3C)##
## Call:errorsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + SLOPE + Avg_z,
## data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.150527 -0.558459 -0.045187 0.540349 2.578564
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -65.325177 5.620712 -11.622 < 2.2e-16
## Avg_CEa_15 0.874324 0.083217 10.507 < 2.2e-16
## SLOPE -0.079881 0.024777 -3.224 0.001264
## Avg_z 0.286926 0.021256 13.498 < 2.2e-16
##
## Lambda: 0.98237, LR test value: 97.514, p-value: < 2.22e-16
## Asymptotic standard error: 0.012433
## z-value: 79.011, p-value: < 2.22e-16
## Wald statistic: 6242.7, p-value: < 2.22e-16
##
## Log likelihood: -408.6476 for error model
## ML residual variance (sigma squared): 0.77863, (sigma: 0.8824)
## Number of observations: 313
## Number of parameters estimated: 6
## AIC: 829.3, (AIC for lm: 924.81)Comprobando si el modelo map3C es adecuado–>Residuales:normalidad e independencia
#####Normalidad
residuales_map3C =map3C$residuals
shapiro.test(residuales_map3C)##
## Shapiro-Wilk normality test
##
## data: residuales_map3C
## W = 0.99348, p-value = 0.1948#####Independencia
Moran.I(residuales_map3C,We)## $observed
## [1] 0.1282906
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.00465405
##
## $p.value
## [1] 0Modelo SEM
map4=lagsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map4)##
## Call:lagsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.144807 -0.520913 -0.027436 0.549523 2.438704
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -60.552051 4.021615 -15.0567 < 2.2e-16
## Avg_CEa_15 0.831698 0.075608 11.0001 < 2.2e-16
## NDVI -2.108385 1.721660 -1.2246 0.2207171
## DEM 0.021807 0.017524 1.2445 0.2133291
## SLOPE -0.085893 0.023033 -3.7291 0.0001921
## Avg_z 0.212051 0.022919 9.2520 < 2.2e-16
##
## Rho: 0.98171, LR test value: 117.07, p-value: < 2.22e-16
## Asymptotic standard error: 0.012881
## z-value: 76.216, p-value: < 2.22e-16
## Wald statistic: 5808.9, p-value: < 2.22e-16
##
## Log likelihood: -397.2464 for lag model
## ML residual variance (sigma squared): 0.7241, (sigma: 0.85094)
## Number of observations: 313
## Number of parameters estimated: 8
## AIC: 810.49, (AIC for lm: 925.56)
## LM test for residual autocorrelation
## test value: 594.27, p-value: < 2.22e-16#Comprobando si el modelo map4 es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map4 =map4$residuals
shapiro.test(residuales_map4)##
## Shapiro-Wilk normality test
##
## data: residuales_map4
## W = 0.99368, p-value = 0.2154#####Independencia
Moran.I(residuales_map4,We)## $observed
## [1] 0.1120975
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004653229
##
## $p.value
## [1] 0###Modelo map4 sin NDVI
map4B=lagsarlm(formula=Avg_CEa_07~Avg_CEa_15+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map4B)##
## Call:lagsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + DEM + SLOPE + Avg_z,
## data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.138949 -0.529797 -0.020674 0.558870 2.458700
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -61.963159 3.862441 -16.0425 < 2.2e-16
## Avg_CEa_15 0.842203 0.075301 11.1845 < 2.2e-16
## DEM 0.023958 0.017477 1.3708 0.170425
## SLOPE -0.089126 0.022936 -3.8859 0.000102
## Avg_z 0.207258 0.022636 9.1560 < 2.2e-16
##
## Rho: 0.98167, LR test value: 116.78, p-value: < 2.22e-16
## Asymptotic standard error: 0.012907
## z-value: 76.055, p-value: < 2.22e-16
## Wald statistic: 5784.4, p-value: < 2.22e-16
##
## Log likelihood: -397.9945 for lag model
## ML residual variance (sigma squared): 0.72757, (sigma: 0.85298)
## Number of observations: 313
## Number of parameters estimated: 7
## AIC: 809.99, (AIC for lm: 924.77)
## LM test for residual autocorrelation
## test value: 589.69, p-value: < 2.22e-16#Comprobando si el modelo map4B es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map4B =map4B$residuals
shapiro.test(residuales_map4B)##
## Shapiro-Wilk normality test
##
## data: residuales_map4B
## W = 0.99448, p-value = 0.3199#####Independencia
Moran.I(residuales_map4B,We)## $observed
## [1] 0.1115907
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004653691
##
## $p.value
## [1] 0###Modelo map4 sin NDVI ni DEM
map4C=lagsarlm(formula=Avg_CEa_07~Avg_CEa_15+SLOPE+Avg_z,data= X, listw= Wve)
summary(map4C)##
## Call:lagsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + SLOPE + Avg_z, data = X,
## listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.201902 -0.551581 -0.027206 0.555763 2.423502
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -61.740674 3.870404 -15.9520 < 2.2e-16
## Avg_CEa_15 0.841509 0.075524 11.1423 < 2.2e-16
## SLOPE -0.093277 0.022803 -4.0905 4.304e-05
## Avg_z 0.230575 0.014985 15.3874 < 2.2e-16
##
## Rho: 0.98174, LR test value: 116.95, p-value: < 2.22e-16
## Asymptotic standard error: 0.012863
## z-value: 76.324, p-value: < 2.22e-16
## Wald statistic: 5825.3, p-value: < 2.22e-16
##
## Log likelihood: -398.9313 for lag model
## ML residual variance (sigma squared): 0.73193, (sigma: 0.85553)
## Number of observations: 313
## Number of parameters estimated: 6
## AIC: 809.86, (AIC for lm: 924.81)
## LM test for residual autocorrelation
## test value: 576.52, p-value: < 2.22e-16#Comprobando si el modelo map4C es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map4C =map4C$residuals
shapiro.test(residuales_map4C)##
## Shapiro-Wilk normality test
##
## data: residuales_map4C
## W = 0.9954, p-value = 0.4834#####Independencia
Moran.I(residuales_map4C,We)## $observed
## [1] 0.1103959
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.00465375
##
## $p.value
## [1] 0#Modelo SDE (variables explicativas con dependencia espacial-->matriz)
map5=lagsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve, type="mixed")
summary(map5)##
## Call:lagsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.602120 -0.556661 0.050368 0.549685 2.255787
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -100.419904 32.629529 -3.0776 0.002087
## Avg_CEa_15 0.931130 0.089892 10.3583 < 2.2e-16
## NDVI -2.604775 2.023958 -1.2870 0.198104
## DEM 0.018410 0.027985 0.6579 0.510626
## SLOPE -0.023314 0.024255 -0.9612 0.336464
## Avg_z 0.221732 0.036028 6.1545 7.531e-10
## lag.Avg_CEa_15 0.941605 0.826775 1.1389 0.254749
## lag.NDVI 55.886306 17.759652 3.1468 0.001651
## lag.DEM -0.061709 0.140317 -0.4398 0.660095
## lag.SLOPE -1.660298 0.269954 -6.1503 7.734e-10
## lag.Avg_z -0.037005 0.170947 -0.2165 0.828619
##
## Rho: 0.96661, LR test value: 54.061, p-value: 1.944e-13
## Asymptotic standard error: 0.023511
## z-value: 41.113, p-value: < 2.22e-16
## Wald statistic: 1690.3, p-value: < 2.22e-16
##
## Log likelihood: -375.6066 for mixed model
## ML residual variance (sigma squared): 0.63304, (sigma: 0.79564)
## Number of observations: 313
## Number of parameters estimated: 13
## AIC: 777.21, (AIC for lm: 829.27)
## LM test for residual autocorrelation
## test value: 336.46, p-value: < 2.22e-16#Modelo map5 sin DEM
map5B=lagsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+SLOPE+Avg_z,data= X, listw= Wve, type="mixed")
summary(map5B)##
## Call:lagsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + SLOPE + Avg_z,
## data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.621263 -0.558569 0.043779 0.550417 2.303183
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -100.733916 32.614023 -3.0887 0.002011
## Avg_CEa_15 0.940314 0.088664 10.6053 < 2.2e-16
## NDVI -2.704756 2.020080 -1.3389 0.180592
## SLOPE -0.023843 0.024251 -0.9832 0.325518
## Avg_z 0.231045 0.033233 6.9524 3.592e-12
## lag.Avg_CEa_15 0.864600 0.818549 1.0563 0.290850
## lag.NDVI 55.739319 17.768007 3.1371 0.001707
## lag.SLOPE -1.646244 0.267055 -6.1644 7.073e-10
## lag.Avg_z -0.081667 0.128877 -0.6337 0.526287
##
## Rho: 0.96609, LR test value: 53.69, p-value: 2.347e-13
## Asymptotic standard error: 0.023872
## z-value: 40.47, p-value: < 2.22e-16
## Wald statistic: 1637.8, p-value: < 2.22e-16
##
## Log likelihood: -375.833 for mixed model
## ML residual variance (sigma squared): 0.63402, (sigma: 0.79626)
## Number of observations: 313
## Number of parameters estimated: 11
## AIC: 773.67, (AIC for lm: 825.36)
## LM test for residual autocorrelation
## test value: 335.86, p-value: < 2.22e-16#Comprobando si el modelo map5B es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map5B =map5B$residuals
shapiro.test(residuales_map5B)##
## Shapiro-Wilk normality test
##
## data: residuales_map5B
## W = 0.99678, p-value = 0.7883#####Independencia
Moran.I(residuales_map5B,We)## $observed
## [1] 0.08279494
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004652996
##
## $p.value
## [1] 0##Modelo GNS= Modelo SARAR con type=mixed y variables epxlicativas
map6= sacsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+DEM+SLOPE+Avg_z, data= X, listw= Wve, type="mixed")
summary(map6)##
## Call:sacsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.380845 -0.510556 0.011811 0.459913 2.070208
##
## Type: sacmixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -95.095972 90.029432 -1.0563 0.290842
## Avg_CEa_15 0.943705 0.083668 11.2792 < 2.2e-16
## NDVI -2.153544 1.866171 -1.1540 0.248504
## DEM 0.020193 0.027389 0.7373 0.460957
## SLOPE -0.025529 0.022301 -1.1447 0.252323
## Avg_z 0.182539 0.038180 4.7810 1.744e-06
## lag.Avg_CEa_15 0.427790 1.224360 0.3494 0.726790
## lag.NDVI 39.077754 19.721890 1.9814 0.047542
## lag.DEM -0.110948 0.154355 -0.7188 0.472272
## lag.SLOPE -1.355530 0.435598 -3.1119 0.001859
## lag.Avg_z 0.127292 0.362028 0.3516 0.725131
##
## Rho: 0.9625
## Asymptotic standard error: 0.57593
## z-value: 1.6712, p-value: 0.09468
## Lambda: 0.96455
## Asymptotic standard error: 0.54509
## z-value: 1.7695, p-value: 0.076807
##
## LR test value: 206.43, p-value: < 2.22e-16
##
## Log likelihood: -352.5646 for sacmixed model
## ML residual variance (sigma squared): 0.5365, (sigma: 0.73246)
## Number of observations: 313
## Number of parameters estimated: 14
## AIC: 733.13, (AIC for lm: 925.56)#Comprobando si el modelo map6 es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map6 =map6$residuals
shapiro.test(residuales_map6)##
## Shapiro-Wilk normality test
##
## data: residuales_map6
## W = 0.99564, p-value = 0.5343#####Independencia
Moran.I(residuales_map6,We)## $observed
## [1] 0.07793079
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004652667
##
## $p.value
## [1] 0#Modelo 6 sin DEM
map6B= sacsarlm(formula=Avg_CEa_07~Avg_CEa_15+NDVI+SLOPE+Avg_z, data= X, listw= Wve, type="mixed")
summary(map6B)##
## Call:sacsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + SLOPE + Avg_z,
## data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.398592 -0.490370 0.012676 0.501434 2.113961
##
## Type: sacmixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -98.280498 85.714775 -1.1466 0.251547
## Avg_CEa_15 0.954131 0.082172 11.6114 < 2.2e-16
## NDVI -2.236429 1.863287 -1.2003 0.230039
## SLOPE -0.025516 0.022313 -1.1435 0.252813
## Avg_z 0.189962 0.037625 5.0488 4.446e-07
## lag.Avg_CEa_15 0.388549 1.213126 0.3203 0.748751
## lag.NDVI 39.332985 19.688965 1.9977 0.045747
## lag.SLOPE -1.337654 0.415882 -3.2164 0.001298
## lag.Avg_z 0.045157 0.308110 0.1466 0.883479
##
## Rho: 0.96203
## Asymptotic standard error: 0.53778
## z-value: 1.7889, p-value: 0.073632
## Lambda: 0.96424
## Asymptotic standard error: 0.50719
## z-value: 1.9011, p-value: 0.057283
##
## LR test value: 207.52, p-value: < 2.22e-16
##
## Log likelihood: -352.8831 for sacmixed model
## ML residual variance (sigma squared): 0.53766, (sigma: 0.73325)
## Number of observations: 313
## Number of parameters estimated: 12
## AIC: 729.77, (AIC for lm: 925.28)#Comprobando si el modelo map6B es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map6B =map6B$residuals
shapiro.test(residuales_map6B)##
## Shapiro-Wilk normality test
##
## data: residuales_map6B
## W = 0.99569, p-value = 0.5456#####Independencia
Moran.I(residuales_map6B,We)## $observed
## [1] 0.07741073
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004652274
##
## $p.value
## [1] 0#Modelo 6 sin DEM ni NDVI
map6C= sacsarlm(formula=Avg_CEa_07~Avg_CEa_15+SLOPE+Avg_z, data= X, listw= Wve, type="mixed")
summary(map6C)##
## Call:sacsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + SLOPE + Avg_z, data = X,
## listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.345410 -0.459164 0.023555 0.494094 2.080937
##
## Type: sacmixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -77.287189 86.625783 -0.8922 0.37229
## Avg_CEa_15 0.949660 0.082066 11.5719 < 2.2e-16
## SLOPE -0.024447 0.022450 -1.0890 0.27617
## Avg_z 0.177499 0.037337 4.7540 1.994e-06
## lag.Avg_CEa_15 0.107889 1.195788 0.0902 0.92811
## lag.SLOPE -1.045625 0.361782 -2.8902 0.00385
## lag.Avg_z 0.126374 0.320796 0.3939 0.69363
##
## Rho: 0.96407
## Asymptotic standard error: 0.57562
## z-value: 1.6748, p-value: 0.093966
## Lambda: 0.96668
## Asymptotic standard error: 0.53452
## z-value: 1.8085, p-value: 0.070529
##
## LR test value: 204.71, p-value: < 2.22e-16
##
## Log likelihood: -355.0471 for sacmixed model
## ML residual variance (sigma squared): 0.5447, (sigma: 0.73804)
## Number of observations: 313
## Number of parameters estimated: 10
## AIC: 730.09, (AIC for lm: 924.81)#Comprobando si el modelo map6B es adecuado-->Residuales:normalidad e independencia
#####Normalidad
residuales_map6C =map6C$residuals
shapiro.test(residuales_map6C)##
## Shapiro-Wilk normality test
##
## data: residuales_map6C
## W = 0.99456, p-value = 0.3319#####Independencia
Moran.I(residuales_map6C,We)## $observed
## [1] 0.08250214
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004651312
##
## $p.value
## [1] 0\(Tarea\)
#Regresion multiple
plot(df[, 3:8], main = "Matriz de correlación")
colnames(df)## [1] "Avg_X_MCB" "Avg_Y_MCE" "Avg_CEa_07" "Avg_CEa_15" "NDVI"
## [6] "DEM" "SLOPE" "Avg_z"mode=lm(formula=df$Avg_CEa_15~df$Avg_X_MCB+df$Avg_Y_MCE+df$Avg_CEa_07+df$NDVI+df$DEM+df$SLOPE+df$Avg_z)
summary(mode)##
## Call:
## lm(formula = df$Avg_CEa_15 ~ df$Avg_X_MCB + df$Avg_Y_MCE + df$Avg_CEa_07 +
## df$NDVI + df$DEM + df$SLOPE + df$Avg_z)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.5009 -0.3661 -0.0082 0.3716 2.2698
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.018e+02 3.806e+02 1.581 0.1149
## df$Avg_X_MCB -5.540e-04 2.646e-04 -2.094 0.0371 *
## df$Avg_Y_MCE -9.298e-05 4.870e-04 -0.191 0.8487
## df$Avg_CEa_07 2.713e-01 2.832e-02 9.579 < 2e-16 ***
## df$NDVI -1.458e+00 1.140e+00 -1.280 0.2016
## df$DEM -1.998e-02 1.624e-02 -1.231 0.2194
## df$SLOPE 7.070e-02 1.522e-02 4.646 5.04e-06 ***
## df$Avg_z -1.219e-01 1.802e-02 -6.763 6.89e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5636 on 305 degrees of freedom
## Multiple R-squared: 0.4285, Adjusted R-squared: 0.4154
## F-statistic: 32.67 on 7 and 305 DF, p-value: < 2.2e-16# Al analizar el p-value se ve que nuestras variables AVg_X_MCB,Avg_CEa_07,SLOPE y Avg_z si aportan a nuestro modelo
#Modelo SAR
mapp= spautolm(Avg_CEa_15~1, data= X, listw= Wve, family="SAR")
summary(mapp)#AIC=615.58##
## Call: spautolm(formula = Avg_CEa_15 ~ 1, data = X, listw = Wve, family = "SAR")
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.453255 -0.397645 -0.042934 0.322283 2.953512
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 19.3151 1.5515 12.45 < 2.2e-16
##
## Lambda: 0.97691 LR test value: 86.774 p-value: < 2.22e-16
## Numerical Hessian standard error of lambda: 0.023
##
## Log likelihood: -304.7918
## ML residual variance (sigma squared): 0.40168, (sigma: 0.63378)
## Number of observations: 313
## Number of parameters estimated: 3
## AIC: 615.58#Modelo SARAR
colnames(df)## [1] "Avg_X_MCB" "Avg_Y_MCE" "Avg_CEa_07" "Avg_CEa_15" "NDVI"
## [6] "DEM" "SLOPE" "Avg_z"library(spdep)
map22= sacsarlm(Avg_CEa_15~1, data= X, listw= Wve)## Warning in sacsarlm(Avg_CEa_15 ~ 1, data = X, listw = Wve): inversion of asymptotic covariance matrix failed for tol.solve = 2.22044604925031e-16
## reciprocal condition number = 1.77227e-18 - using numerical Hessian.summary(map22)#AIC=570.68#pvalueRho=2.22 e^-16##
## Call:sacsarlm(formula = Avg_CEa_15 ~ 1, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.364908 -0.357504 -0.063033 0.289858 2.878760
##
## Type: sac
## Coefficients: (numerical Hessian approximate standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.1253 1.1713 0.9607 0.3367
##
## Rho: 0.95895
## Approximate (numerical Hessian) standard error: 0.040792
## z-value: 23.508, p-value: < 2.22e-16
## Lambda: 0.95895
## Approximate (numerical Hessian) standard error: 0.040755
## z-value: 23.53, p-value: < 2.22e-16
##
## LR test value: 133.68, p-value: < 2.22e-16
##
## Log likelihood: -281.3411 for sac model
## ML residual variance (sigma squared): 0.34087, (sigma: 0.58384)
## Number of observations: 313
## Number of parameters estimated: 4
## AIC: 570.68, (AIC for lm: 700.36)residuales_map22 =map22$residuals
shapiro.test(residuales_map22)##
## Shapiro-Wilk normality test
##
## data: residuales_map22
## W = 0.94498, p-value = 2.076e-09Moran.I(residuales_map22,We)## $observed
## [1] 0.0563745
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004629181
##
## $p.value
## [1] 0#Modelo GNS
map66= sacsarlm(formula=Avg_CEa_15~Avg_CEa_07+NDVI+DEM+SLOPE+Avg_z, data= X, listw= Wve, type="mixed")
summary(map66)#AIC=407,pvalueRho=0.12487##
## Call:sacsarlm(formula = Avg_CEa_15 ~ Avg_CEa_07 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.259920 -0.257064 -0.024628 0.256195 1.934415
##
## Type: sacmixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 26.084548 36.251116 0.7196 0.47180
## Avg_CEa_07 0.362195 0.034530 10.4893 < 2.2e-16
## NDVI -0.433013 1.115198 -0.3883 0.69781
## DEM 0.025591 0.016294 1.5706 0.11628
## SLOPE 0.010971 0.013341 0.8224 0.41087
## Avg_z -0.112125 0.022705 -4.9384 7.875e-07
## lag.Avg_CEa_07 -0.486187 0.219517 -2.2148 0.02677
## lag.NDVI -29.360420 11.602286 -2.5306 0.01139
## lag.DEM -0.125137 0.091592 -1.3663 0.17186
## lag.SLOPE 0.883938 0.212106 4.1674 3.080e-05
## lag.Avg_z 0.205119 0.155215 1.3215 0.18633
##
## Rho: 0.9037
## Asymptotic standard error: 0.58886
## z-value: 1.5347, p-value: 0.12487
## Lambda: 0.94686
## Asymptotic standard error: 0.32846
## z-value: 2.8827, p-value: 0.0039429
##
## LR test value: 149.65, p-value: < 2.22e-16
##
## Log likelihood: -189.5009 for sacmixed model
## ML residual variance (sigma squared): 0.19093, (sigma: 0.43695)
## Number of observations: 313
## Number of parameters estimated: 14
## AIC: 407, (AIC for lm: 542.65)residuales_map66 =map66$residuals
shapiro.test(residuales_map66)##
## Shapiro-Wilk normality test
##
## data: residuales_map66
## W = 0.97364, p-value = 1.68e-05Moran.I(residuales_map66,We)## $observed
## [1] 0.04788866
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004639083
##
## $p.value
## [1] 0##MODELO GNS sin DEM
map66b= sacsarlm(formula=Avg_CEa_15~Avg_CEa_07+NDVI+SLOPE+Avg_z, data= X, listw= Wve, type="mixed")
summary(map66b)#AIC=405.6,pvalueRho=0.10807##
## Call:sacsarlm(formula = Avg_CEa_15 ~ Avg_CEa_07 + NDVI + SLOPE + Avg_z,
## data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.242977 -0.273753 -0.027046 0.262976 1.942185
##
## Type: sacmixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 17.956728 33.982429 0.5284 0.597213
## Avg_CEa_07 0.373281 0.034028 10.9697 < 2.2e-16
## NDVI -0.543373 1.116766 -0.4866 0.626571
## SLOPE 0.011352 0.013389 0.8479 0.396517
## Avg_z -0.105276 0.022422 -4.6952 2.664e-06
## lag.Avg_CEa_07 -0.559944 0.212632 -2.6334 0.008454
## lag.NDVI -28.742233 11.599813 -2.4778 0.013219
## lag.SLOPE 0.884282 0.203207 4.3516 1.351e-05
## lag.Avg_z 0.138273 0.128293 1.0778 0.281128
##
## Rho: 0.90567
## Asymptotic standard error: 0.5636
## z-value: 1.6069, p-value: 0.10807
## Lambda: 0.94735
## Asymptotic standard error: 0.31799
## z-value: 2.9792, p-value: 0.0028904
##
## LR test value: 147.85, p-value: < 2.22e-16
##
## Log likelihood: -190.7997 for sacmixed model
## ML residual variance (sigma squared): 0.19248, (sigma: 0.43873)
## Number of observations: 313
## Number of parameters estimated: 12
## AIC: 405.6, (AIC for lm: 541.44)residuales_map66b =map66b$residuals
shapiro.test(residuales_map66b)##
## Shapiro-Wilk normality test
##
## data: residuales_map66b
## W = 0.97534, p-value = 3.293e-05Moran.I(residuales_map66b,We)## $observed
## [1] 0.04841028
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004640051
##
## $p.value
## [1] 0#Modelo SEM
map44=lagsarlm(formula=Avg_CEa_15~Avg_CEa_07+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map44)#AIC=493.35, pvalueRho=7.9403 e^-13##
## Call:lagsarlm(formula = Avg_CEa_15 ~ Avg_CEa_07 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.4682926 -0.3245699 0.0049751 0.3215294 1.9926400
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 22.7597381 2.0758704 10.9639 < 2.2e-16
## Avg_CEa_07 0.2479677 0.0247839 10.0052 < 2.2e-16
## NDVI -1.2527443 1.0391615 -1.2055 0.2280
## DEM -0.0035254 0.0106002 -0.3326 0.7394
## SLOPE 0.0603502 0.0137383 4.3928 1.119e-05
## Avg_z -0.1133123 0.0145926 -7.7650 8.216e-15
##
## Rho: 0.96209, LR test value: 51.297, p-value: 7.9403e-13
## Asymptotic standard error: 0.02667
## z-value: 36.074, p-value: < 2.22e-16
## Wald statistic: 1301.3, p-value: < 2.22e-16
##
## Log likelihood: -238.6759 for lag model
## ML residual variance (sigma squared): 0.26412, (sigma: 0.51393)
## Number of observations: 313
## Number of parameters estimated: 8
## AIC: 493.35, (AIC for lm: 542.65)
## LM test for residual autocorrelation
## test value: 197.83, p-value: < 2.22e-16residuales_map44 =map44$residuals
shapiro.test(residuales_map44)##
## Shapiro-Wilk normality test
##
## data: residuales_map44
## W = 0.98401, p-value = 0.001482Moran.I(residuales_map44,We)## $observed
## [1] 0.06272672
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004645596
##
## $p.value
## [1] 0#Modelo SLM
map33=errorsarlm(formula=Avg_CEa_15~Avg_CEa_07+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve)
summary(map3)#AIC=828.2##
## Call:errorsarlm(formula = Avg_CEa_07 ~ Avg_CEa_15 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.019160 -0.540466 -0.045367 0.513314 2.592838
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -64.737579 5.752902 -11.2530 < 2.2e-16
## Avg_CEa_15 0.859898 0.083054 10.3535 < 2.2e-16
## NDVI -2.395368 1.907913 -1.2555 0.209301
## DEM 0.036792 0.020974 1.7542 0.079402
## SLOPE -0.073067 0.024760 -2.9510 0.003168
## Avg_z 0.257034 0.028465 9.0299 < 2.2e-16
##
## Lambda: 0.9825, LR test value: 99.359, p-value: < 2.22e-16
## Asymptotic standard error: 0.012342
## z-value: 79.604, p-value: < 2.22e-16
## Wald statistic: 6336.8, p-value: < 2.22e-16
##
## Log likelihood: -406.1005 for error model
## ML residual variance (sigma squared): 0.76603, (sigma: 0.87523)
## Number of observations: 313
## Number of parameters estimated: 8
## AIC: 828.2, (AIC for lm: 925.56)residuales_map33 =map33$residuals
shapiro.test(residuales_map33)##
## Shapiro-Wilk normality test
##
## data: residuales_map33
## W = 0.98846, p-value = 0.01377Moran.I(residuales_map33,We)## $observed
## [1] 0.07550844
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004647281
##
## $p.value
## [1] 0#Modelo SDE
map55=lagsarlm(formula=Avg_CEa_15~Avg_CEa_07+NDVI+DEM+SLOPE+Avg_z,data= X, listw= Wve, type="mixed")
summary(map55)#AIC=431.12, pvalueRho=5.8806e^06##
## Call:lagsarlm(formula = Avg_CEa_15 ~ Avg_CEa_07 + NDVI + DEM + SLOPE +
## Avg_z, data = X, listw = Wve, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.356577 -0.269279 -0.016311 0.269347 1.954889
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 28.265444 13.304190 2.1246 0.033624
## Avg_CEa_07 0.318816 0.033502 9.5163 < 2.2e-16
## NDVI -0.463244 1.170570 -0.3957 0.692295
## DEM 0.026816 0.016575 1.6179 0.105686
## SLOPE 0.010237 0.014082 0.7269 0.467260
## Avg_z -0.127758 0.021328 -5.9901 2.098e-09
## lag.Avg_CEa_07 -0.148588 0.197032 -0.7541 0.450771
## lag.NDVI -33.518708 10.304158 -3.2529 0.001142
## lag.DEM -0.135244 0.080945 -1.6708 0.094759
## lag.SLOPE 1.004251 0.157352 6.3822 1.746e-10
## lag.Avg_z 0.216661 0.099651 2.1742 0.029690
##
## Rho: 0.91953, LR test value: 20.527, p-value: 5.8806e-06
## Asymptotic standard error: 0.056295
## z-value: 16.334, p-value: < 2.22e-16
## Wald statistic: 266.81, p-value: < 2.22e-16
##
## Log likelihood: -202.5621 for mixed model
## ML residual variance (sigma squared): 0.21072, (sigma: 0.45905)
## Number of observations: 313
## Number of parameters estimated: 13
## AIC: 431.12, (AIC for lm: 449.65)
## LM test for residual autocorrelation
## test value: 128.64, p-value: < 2.22e-16residuales_map55 =map55$residuals
shapiro.test(residuales_map55)##
## Shapiro-Wilk normality test
##
## data: residuales_map55
## W = 0.97531, p-value = 3.25e-05Moran.I(residuales_map55,We)## $observed
## [1] 0.04750286
##
## $expected
## [1] -0.003205128
##
## $sd
## [1] 0.004640154
##
## $p.value
## [1] 0