Spatial Crossregressive
Hani
2022-10-23
LIBRARY
library(tmap)
library(sp)
library(raster)
library(spdep)## 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.9.1, GDAL 3.4.3, PROJ 7.2.1; sf_use_s2() is TRUElibrary(spData)
library(spatial)
library(spatialreg)## Loading required package: Matrix##
## Attaching package: 'spatialreg'## The following objects are masked from 'package:spdep':
##
## get.ClusterOption, get.coresOption, get.mcOption,
## get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOptionlibrary(rgdal)## Please note that rgdal will be retired by the end of 2023,
## plan transition to sf/stars/terra functions using GDAL and PROJ
## at your earliest convenience.
##
## rgdal: version: 1.5-32, (SVN revision 1176)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.4.3, released 2022/04/22
## Path to GDAL shared files: C:/Users/hanik/AppData/Local/R/win-library/4.2/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 7.2.1, January 1st, 2021, [PJ_VERSION: 721]
## Path to PROJ shared files: C:/Users/hanik/AppData/Local/R/win-library/4.2/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.5-0
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal.PETA JATIM
#Input Shp File
spJatim = shapefile("jatim baru.shp")## Warning in .local(x, ...): .prj file is missingnames(spJatim)## [1] "POLY_ID" "KODE_KAB" "NAMA_KAB" "KODE_PROP" "NAMA_PROP"
## [6] "KODE" "Y" "X1" "X2" "X3"
## [11] "X4" "X5" "X6" "X7" "X8"
## [16] "X9" "ZX1" "ZX2" "ZX3" "ZX4"
## [21] "ZX5" "ZX6" "ZX7" "ZX8" "ZX9"
## [26] "ZY" "COVID" "PNEUMONIA" "KEMISKINAN" "IPM"
## [31] "IG" "NOMOR"#Menampilkan Peta Jawa Timur
tm_shape(spJatim) + tm_polygons()## Warning: Currect projection of shape spJatim unknown. Long-lat (WGS84) is
## assumed.
Membuat Pemmbobot Queen
queen.nb=poly2nb(spJatim) #Pembobot queen
queen.listw=nb2listw(queen.nb) #convert nb to listw type
queen.jatim= queen.listw#PERSAMAAN REGRESI
reg.eq=Y~X1+X2+X3+X4+X5+X6+X7+X8#OLS
reg.OLS=lm(reg.eq,data=spJatim)
summary(reg.OLS)##
## Call:
## lm(formula = reg.eq, data = spJatim)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.056820 -0.013321 -0.002441 0.020950 0.058723
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.196e+00 7.576e-01 2.899 0.00707 **
## X1 -1.335e-06 4.420e-06 -0.302 0.76485
## X2 2.684e-01 1.359e-01 1.975 0.05780 .
## X3 -4.031e-02 1.525e-01 -0.264 0.79345
## X4 -2.639e-01 3.166e-01 -0.834 0.41128
## X5 1.871e-01 3.979e-01 0.470 0.64168
## X6 1.017e+01 4.677e+00 2.175 0.03791 *
## X7 2.583e+01 8.658e+00 2.983 0.00574 **
## X8 -2.314e+00 7.685e-01 -3.011 0.00535 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02866 on 29 degrees of freedom
## Multiple R-squared: 0.8558, Adjusted R-squared: 0.816
## F-statistic: 21.51 on 8 and 29 DF, p-value: 3.093e-10#SLX (SCR)
reg.SLX<- lmSLX(reg.eq, data=spJatim, listw=queen.jatim)
summary(reg.SLX)##
## Call:
## lm(formula = formula(paste("y ~ ", paste(colnames(x)[-1], collapse = "+"))),
## data = as.data.frame(x), weights = weights)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.032528 -0.016371 0.001995 0.012736 0.046155
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.508e+00 1.691e+00 2.666 0.01445 *
## X1 1.459e-06 5.540e-06 0.263 0.79483
## X2 1.149e-01 2.010e-01 0.572 0.57354
## X3 9.172e-02 1.577e-01 0.582 0.56692
## X4 -6.493e-01 3.463e-01 -1.875 0.07476 .
## X5 -3.516e-01 5.928e-01 -0.593 0.55938
## X6 1.559e+01 6.142e+00 2.539 0.01911 *
## X7 3.818e+01 1.116e+01 3.421 0.00257 **
## X8 -2.718e+00 1.097e+00 -2.477 0.02182 *
## lag.X1 -7.246e-06 1.483e-05 -0.488 0.63028
## lag.X2 -3.082e-02 2.932e-01 -0.105 0.91728
## lag.X3 2.523e-01 3.114e-01 0.810 0.42684
## lag.X4 -1.211e+00 9.401e-01 -1.288 0.21171
## lag.X5 7.302e-01 1.189e+00 0.614 0.54573
## lag.X6 1.953e+01 1.199e+01 1.630 0.11807
## lag.X7 -2.457e+01 2.071e+01 -1.186 0.24875
## lag.X8 -1.984e+00 1.445e+00 -1.373 0.18409
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02698 on 21 degrees of freedom
## Multiple R-squared: 0.9075, Adjusted R-squared: 0.837
## F-statistic: 12.88 on 16 and 21 DF, p-value: 1.883e-07summary(impacts(reg.SLX))## Impact measures (SlX, estimable, n-k):
## Direct Indirect Total
## X1 1.459065e-06 -7.246157e-06 -5.787092e-06
## X2 1.149198e-01 -3.082281e-02 8.409699e-02
## X3 9.172361e-02 2.523143e-01 3.440379e-01
## X4 -6.493476e-01 -1.210958e+00 -1.860305e+00
## X5 -3.516444e-01 7.302236e-01 3.785791e-01
## X6 1.559181e+01 1.953278e+01 3.512460e+01
## X7 3.817770e+01 -2.456914e+01 1.360856e+01
## X8 -2.717666e+00 -1.983995e+00 -4.701662e+00
## ========================================================
## Standard errors:
## Direct Indirect Total
## X1 5.539873e-06 1.483425e-05 1.534032e-05
## X2 2.009878e-01 2.932203e-01 2.412908e-01
## X3 1.576647e-01 3.113766e-01 3.623638e-01
## X4 3.463057e-01 9.400714e-01 1.103039e+00
## X5 5.927942e-01 1.189046e+00 1.121541e+00
## X6 6.142045e+00 1.198532e+01 1.508667e+01
## X7 1.115918e+01 2.071073e+01 1.885669e+01
## X8 1.097027e+00 1.444513e+00 1.736561e+00
## ========================================================
## Z-values:
## Direct Indirect Total
## X1 0.2633751 -0.4884748 -0.3772471
## X2 0.5717750 -0.1051182 0.3485297
## X3 0.5817636 0.8103187 0.9494267
## X4 -1.8750703 -1.2881549 -1.6865276
## X5 -0.5931982 0.6141256 0.3375525
## X6 2.5385376 1.6297253 2.3281874
## X7 3.4211923 -1.1863000 0.7216833
## X8 -2.4773013 -1.3734700 -2.7074559
##
## p-values:
## Direct Indirect Total
## X1 0.79226148 0.62521 0.7059900
## X2 0.56747445 0.91628 0.7274424
## X3 0.56072591 0.41776 0.3424037
## X4 0.06078305 0.19769 0.0916942
## X5 0.55304852 0.53913 0.7357004
## X6 0.01113168 0.10316 0.0199022
## X7 0.00062347 0.23550 0.4704892
## X8 0.01323801 0.16961 0.0067801#Perhitungan AIC
AIC(reg.SLX)## [1] -153.2672AIC(reg.OLS)## [1] -152.3906#Perbandingan Model dengan OLS
LR.Sarlm(reg.SLX,reg.OLS) # OLS vs SCR##
## Likelihood ratio for spatial linear models
##
## data:
## Likelihood ratio = 16.877, df = 8, p-value = 0.03142
## sample estimates:
## Log likelihood of reg.SLX Log likelihood of reg.OLS
## 94.63362 86.19530#Perbaikan Model SCR
## Mengambil Matriks Variabel Prediktor dari model OLS
X=model.matrix(reg.OLS)
#Membentuk matriks lag spasial
lagX = create_WX(X,queen.jatim,prefix = "lag")
spJatim2=cbind(spJatim,lagX)
#Membuat model SCR secara manual
reg.SLX2=lm(Y~X1+X2+X3+X4+X5+X6+X7+X8+lag.X4
+lag.X6,data=spJatim2)
summary(reg.SLX2)##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + lag.X4 +
## lag.X6, data = spJatim2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.038196 -0.017165 0.000368 0.018216 0.045651
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.699e+00 7.057e-01 3.825 0.000702 ***
## X1 -2.325e-06 4.004e-06 -0.581 0.566334
## X2 5.075e-02 1.494e-01 0.340 0.736697
## X3 4.178e-02 1.406e-01 0.297 0.768589
## X4 -6.389e-01 3.191e-01 -2.003 0.055365 .
## X5 1.861e-02 3.605e-01 0.052 0.959216
## X6 1.187e+01 4.722e+00 2.514 0.018217 *
## X7 3.610e+01 8.453e+00 4.270 0.000216 ***
## X8 -2.884e+00 7.174e-01 -4.020 0.000419 ***
## lag.X4 -8.769e-01 7.801e-01 -1.124 0.270821
## lag.X6 2.119e+01 7.934e+00 2.670 0.012674 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02562 on 27 degrees of freedom
## Multiple R-squared: 0.8927, Adjusted R-squared: 0.853
## F-statistic: 22.47 on 10 and 27 DF, p-value: 1.448e-10AIC(reg.SLX2)## [1] -159.6419