Library

Input

Input Data

setwd("C:/Unpad/Semester 5/Spatial")
Kemiskinan <- read.csv("Data Spatial Garis Kemiskinan.csv", sep = ";")

Input Map

Jateng<-st_read("RBI_50K_2023_Jawa Tengah/RBI_50K_2023_Jawa Tengah.shp")
## Reading layer `RBI_50K_2023_Jawa Tengah' from data source 
##   `C:\Unpad\Semester 5\Spatial\RBI_50K_2023_Jawa Tengah\RBI_50K_2023_Jawa Tengah.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 35 features and 25 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY, XYZ
## Bounding box:  xmin: 108.5558 ymin: -8.211806 xmax: 111.6914 ymax: -5.725371
## z_range:       zmin: 0 zmax: 0
## Geodetic CRS:  WGS 84 + EGM2008 height
Jateng<-st_zm(Jateng)
Jateng<-as_Spatial(Jateng)
plot(Jateng)

Jateng$NAMOBJ
##  [1] "Banjarnegara"    "Banyumas"        "Batang"          "Blora"          
##  [5] "Boyolali"        "Brebes"          "Cilacap"         "Demak"          
##  [9] "Grobogan"        "Jepara"          "Karanganyar"     "Kebumen"        
## [13] "Kendal"          "Klaten"          "Kota Magelang"   "Kota Pekalongan"
## [17] "Kota Salatiga"   "Kota Semarang"   "Kota Surakarta"  "Kota Tegal"     
## [21] "Kudus"           "Magelang"        "Pati"            "Pekalongan"     
## [25] "Pemalang"        "Purbalingga"     "Purworejo"       "Rembang"        
## [29] "Semarang"        "Sragen"          "Sukoharjo"       "Tegal"          
## [33] "Temanggung"      "Wonogiri"        "Wonosobo"

Data Preparation

Matriks Bobot

id<-c(1:35)
Jateng$id<-c(1:35)
Jateng_sf<-st_as_sf(Jateng)

Koordinat Centroid

CoordK<-coordinates(Jateng)

Rename

Jateng_merged<-Jateng_sf %>%
  left_join(Kemiskinan, by = "NAMOBJ") 

Mapping Data

ggplot() +
  geom_sf(data=Jateng_merged, aes(fill = Kemiskinan),color=NA) +
  theme_bw() +
  scale_fill_gradient(low = "yellow", high = "red") +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())+ 
  theme(legend.position = "right",
        axis.text.x = element_blank(),  # Remove x-axis labels
        axis.text.y = element_blank())+  # Remove y-axis labels   
  labs(title = "",
       fill = "Garis Kemiskinan") 

Mapping data migrasi masuk di 35 Kabupaten/Kota di Provinsi Jawa Tengah berdasarkan gradasi warna dari kuning hingga merah.

Mapping data by Interval

summary(Kemiskinan$Kemiskinan)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  378858  425809  470728  480044  512242  642456
breaks <- c(-Inf, 425809, 470728, 512242, Inf)
labels <- c("Very Low", "Low", "High", "Very High")
Jateng_merged$Kemiskinan_Discrete <- cut(Jateng_merged$Kemiskinan, breaks = breaks, labels = labels, right = TRUE)

ggplot() +
  geom_sf(data=Jateng_merged, aes(fill = Kemiskinan_Discrete),color=NA) +
  theme_bw() +
  scale_fill_manual(values = c("Very Low" = "yellow", 
                               "Low" = "orange",
                               "High" = "red",
                               "Very High" = "red3"))+
  labs(fill = "Garis Kemiskinan")+theme(legend.position = "right",
                                  axis.text.x = element_blank(),  # Remove x-axis labels
                                  axis.text.y = element_blank())+  # Remove y-axis labels   
  labs(title = "",
       fill = "Garis Kemiskinan")  

Mapping data migrasi masuk di 35 Kabupaten/Kota di Provinsi Jawa Tengah berdasarkan kelompok interval (Q1, Q2, dan Q3) dan membentuk gradasi warna dari kuning hingga merah.

K-Nearest Neighbors

W <- knn2nb(knearneigh(CoordK, k = 5), row.names = id) ; W
## Neighbour list object:
## Number of regions: 35 
## Number of nonzero links: 175 
## Percentage nonzero weights: 14.28571 
## Average number of links: 5 
## Non-symmetric neighbours list
WB <- nb2mat(W, style='B', zero.policy = TRUE) #menyajikan dalam bentuk matrix biner "B"
WB[1:5,1:5]
##   [,1] [,2] [,3] [,4] [,5]
## 1    0    0    1    0    0
## 2    0    0    0    0    0
## 3    1    0    0    0    0
## 4    0    0    0    0    0
## 5    0    0    0    0    0
WL<-nb2listw(W) ; WL # Moran's I (List Neighbours)
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 35 
## Number of nonzero links: 175 
## Percentage nonzero weights: 14.28571 
## Average number of links: 5 
## Non-symmetric neighbours list
## 
## Weights style: W 
## Weights constants summary:
##    n   nn S0    S1     S2
## W 35 1225 35 12.44 144.48
plot(Jateng, axes=T, col="gray90")
plot(W,coordinates(Jateng),col="red",add=TRUE)
text(CoordK[,1], CoordK[,2], row.names(Jateng), col="black", cex=0.8, pos=1.5)
points(CoordK[,1], CoordK[,2], pch=19, cex=0.7,col="blue") #hubungan antar lokasi

Dengan 35 Kabupaten/Kota di Provinsi Jawa Tengah yang dianalisis, terdapat hubungan signifikan antar Kabupaten/Kota, dengan 14.29% hubungan memiliki bobot nonzero dan rata-rata 5 hubungan per Kabupaten/Kota.

Data Generating

Mengkonversi Koordinat

# Mengkonversi ke format spasial
Jateng <- st_read("RBI_50K_2023_Jawa Tengah/RBI_50K_2023_Jawa Tengah.shp")
## Reading layer `RBI_50K_2023_Jawa Tengah' from data source 
##   `C:\Unpad\Semester 5\Spatial\RBI_50K_2023_Jawa Tengah\RBI_50K_2023_Jawa Tengah.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 35 features and 25 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY, XYZ
## Bounding box:  xmin: 108.5558 ymin: -8.211806 xmax: 111.6914 ymax: -5.725371
## z_range:       zmin: 0 zmax: 0
## Geodetic CRS:  WGS 84 + EGM2008 height
Jateng <- st_zm(Jateng)

# Mengambil koordinat centroid
centroids <- st_centroid(Jateng)
centroid_df <- as.data.frame(st_coordinates(centroids))
centroid_df$NAMOBJ <- Jateng$NAMOBJ

# Menggabungkan data dengan Kemiskinan
Jateng_merged <- Jateng %>%
  left_join(Kemiskinan, by = "NAMOBJ")

Derajat (Latitude longitude)

data_combined <- merge(Kemiskinan, centroid_df, by = "NAMOBJ")
data_spasial <- data.frame(
  y = data_combined[,2],
  x1 = data_combined[,3],
  x2 = data_combined[,4],
  x3 = data_combined[,5],
  x4 = data_combined[,6],
  longitude = data_combined[,"X"],
  latitude = data_combined[,"Y"]
)

coordinates(data_spasial) <- ~longitude+latitude

Regresi Liniear

OLS Model

model_ols <- lm(y ~ x1 + x2 + x3 + x4, data = data_spasial)

summary(model_ols)
## 
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4, data = data_spasial)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -72239 -39321  -2323  27999  91297 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)   
## (Intercept) -487217.4   175664.2  -2.774  0.00944 **
## x1            11473.0     5186.0   2.212  0.03470 * 
## x2            10820.9     4832.3   2.239  0.03270 * 
## x3             4274.3    22295.5   0.192  0.84926   
## x4              488.4     1799.0   0.271  0.78788   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 46140 on 30 degrees of freedom
## Multiple R-squared:  0.5961, Adjusted R-squared:  0.5423 
## F-statistic: 11.07 on 4 and 30 DF,  p-value: 1.233e-05
Dengan taraf signifikansi 5%, diketahui bahwa variabel tingkat pengangguran terbuka (x1) dan indeks pembangunan manusia (x2) terbukti signifikan, dengan model secara keseluruhan menjelaskan 54.23% variasi dalam variabel dependen.

Model Diagnostics

par(mfrow=c(2,2))
plot(model_ols)

Uji Asumsi Klasik

Multikolinearitas

vif(model_ols)
##       x1       x2       x3       x4 
## 1.226246 7.368807 6.617879 1.335688
Dengan taraf signifikansi 5%, diketahui bahwa tidak terdapat multikolinearitas karena 1 < VIF < 10.

Normalitas

resids <- residuals(model_ols)
shapiro.test(resids)
## 
##  Shapiro-Wilk normality test
## 
## data:  resids
## W = 0.96891, p-value = 0.414
Dengan taraf signifikansi 5%, diketahui bahwa data residual berdistribusi normal karena p-value > 0.05.

Non-Autokorelasi

durbinWatsonTest(model_ols)
##  lag Autocorrelation D-W Statistic p-value
##    1      -0.1148293      2.193399   0.706
##  Alternative hypothesis: rho != 0
Dengan taraf signifikansi 5%, diketahui bahwa tidak terdapat autokorelasi pada residual dari model regresi karena p-value > 0.05.

Uji Asumsi Efek Spasial

Heterogenitas Spasial

bptest(model_ols)
## 
##  studentized Breusch-Pagan test
## 
## data:  model_ols
## BP = 10.547, df = 4, p-value = 0.03215
Dengan taraf signifikansi 5%, diketahui bahwa terdapat heteroskedastisitas pada residual dari model regresi karena p-value < 0.05.

Autokorelasi Spasial

moran_test <- moran.test(resids, WL) ; moran_test
## 
##  Moran I test under randomisation
## 
## data:  resids  
## weights: WL    
## 
## Moran I statistic standard deviate = 2.3619, p-value = 0.009092
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.189432319      -0.029411765       0.008585353
moran_plot <- moran.plot(resids, WL) ; moran_plot

##              x          wx is_inf labels       dfb.1_        dfb.x       dffit
## 1  -48111.9429  -5691.7586  FALSE      1  0.034096286 -0.038405211  0.05135676
## 2   -6270.6629   5143.1163  FALSE      2  0.055848972 -0.008198961  0.05644759
## 3  -72239.1623   -662.4676  FALSE      3  0.119212399 -0.201615686  0.23422314
## 4   17947.5455  22490.9158  FALSE      4  0.159156309  0.066874278  0.17263516
## 5  -71505.8920 -13253.0365   TRUE      5  0.008879215 -0.014864372  0.01731445
## 6   67284.6057   9156.0711  FALSE      6 -0.024213327 -0.038141728 -0.04517828
## 7  -69260.1262  27637.7392  FALSE      7  0.378412060 -0.613590231  0.72089435
## 8   21475.4741  25479.4882  FALSE      8  0.178778458  0.089885189  0.20010268
## 9   27204.2809  -4675.5533  FALSE      9 -0.072659682 -0.046276546 -0.08614493
## 10  26286.9372  32063.5450  FALSE     10  0.227649049  0.140099316  0.26730490
## 11 -39817.2573 -21234.9963  FALSE     11 -0.105308057  0.098166460 -0.14396680
## 12  -2645.1966 -20107.2137  FALSE     12 -0.150417125  0.009315055 -0.15070528
## 13 -24638.7913  -4001.2559  FALSE     13  0.011201077 -0.006461143  0.01293099
## 14 -22443.4848 -34861.6076  FALSE     14 -0.244611679  0.128527975 -0.27632284
## 15  34941.0767 -20153.5584  FALSE     15 -0.214140327 -0.175172098 -0.27666106
## 16  50559.0398 -10726.5287  FALSE     16 -0.162905526 -0.192825658 -0.25242810
## 17 -45473.6342  -6321.0547  FALSE     17  0.024741289 -0.026339782  0.03613745
## 18  28793.5095 -12027.7310  FALSE     18 -0.134947317 -0.090968083 -0.16274511
## 19  21313.7946 -45756.6985   TRUE     19 -0.424739322 -0.211940246 -0.47468111
## 20  91296.6535  28903.9917   TRUE     20  0.113071354  0.241678289  0.26682115
## 21  25949.5224  32131.0280  FALSE     21  0.228715522  0.138948927  0.26761464
## 22 -42997.1589 -13013.8077  FALSE     22 -0.033534200  0.033756538 -0.04758200
## 23  56894.9382  31902.7723  FALSE     23  0.184862594  0.246236745  0.30790699
## 24  49544.2509  -6060.4905  FALSE     24 -0.121525590 -0.140958475 -0.18611223
## 25  41813.0019  35455.6088  FALSE     25  0.235721793  0.230749768  0.32986394
## 26  -2323.3892   5035.2272  FALSE     26  0.049047096 -0.002667875  0.04911960
## 27 -43034.7775 -11921.7035  FALSE     27 -0.024597434  0.024782183 -0.03491691
## 28  58597.6468  30856.6448  FALSE     28  0.173520300  0.238045872  0.29457619
## 29  19163.1910 -19248.4198  FALSE     29 -0.178368412 -0.080023209 -0.19549681
## 30 -56191.8094 -20326.0246  FALSE     30 -0.074417378  0.097898845 -0.12297207
## 31 -38825.0492 -24496.2545  FALSE     31 -0.133645528  0.121477721 -0.18060444
## 32 -11798.5111  38360.0418  FALSE     32  0.346561156 -0.095727623  0.35953917
## 33 -30664.8942  -2865.3957  FALSE     33  0.029493546 -0.021173792  0.03630701
## 34 -10028.4330 -27192.7612  FALSE     34 -0.197719526  0.046420822 -0.20309580
## 35   -795.2958 -25973.5394  FALSE     35 -0.201836617  0.003758018 -0.20187160
##        cov.r       cook.d        hat
## 1  1.1344902 1.358354e-03 0.06482060
## 2  1.0882294 1.637528e-03 0.02918720
## 3  1.1629234 2.790157e-02 0.11029319
## 4  1.0438263 1.496641e-02 0.03361575
## 5  1.1929144 1.545674e-04 0.10864256
## 6  1.1795795 1.051823e-03 0.09946776
## 7  0.9123710 2.349782e-01 0.10369200
## 8  1.0321984 1.997296e-02 0.03579377
## 9  1.0957934 3.805335e-03 0.04016099
## 10 0.9957039 3.493992e-02 0.03939255
## 11 1.0981076 1.056576e-02 0.05339904
## 12 1.0440830 1.143603e-02 0.02868100
## 13 1.1052830 8.620654e-05 0.03807816
## 14 0.9766609 3.703484e-02 0.03645954
## 15 1.0172402 3.766751e-02 0.04769042
## 16 1.0824856 3.199073e-02 0.06860183
## 17 1.1310848 6.729396e-04 0.06095402
## 18 1.0684032 1.340100e-02 0.04155462
## 19 0.7784210 9.760916e-02 0.03568543
## 20 1.2354607 3.628252e-02 0.15909894
## 21 0.9944234 3.500344e-02 0.03911654
## 22 1.1257731 1.166047e-03 0.05752298
## 23 1.0794503 4.725830e-02 0.07926345
## 24 1.1062644 1.759493e-02 0.06701103
## 25 1.0075720 5.306097e-02 0.05595027
## 26 1.0892756 1.240895e-03 0.02865596
## 27 1.1270391 6.282532e-04 0.05757366
## 28 1.0919112 4.343093e-02 0.08234299
## 29 1.0308314 1.906599e-02 0.03432222
## 30 1.1406920 7.754045e-03 0.07801825
## 31 1.0815952 1.651288e-02 0.05217710
## 32 0.8633653 5.912584e-02 0.03075138
## 33 1.1095850 6.790781e-04 0.04329711
## 34 1.0109485 2.042159e-02 0.03014635
## 35 1.0058116 2.014104e-02 0.02858133
locm <- localmoran(resids,WL) ; locm
##              Ii          E.Ii       Var.Ii       Z.Ii Pr(z != E(Ii))
## 1   0.150092760 -3.731533e-02 0.2209802223  0.3986682     0.69013768
## 2  -0.017676659 -6.338819e-04 0.0038968626 -0.2730127     0.78484347
## 3   0.026229964 -8.412535e-02 0.4739636122  0.1602953     0.87264848
## 4   0.221244548 -5.192680e-03 0.0317769824  1.2702567     0.20399320
## 5   0.519418723 -8.242617e-02 0.4652519841  0.8823491     0.37758804
## 6   0.337664375 -7.298152e-02 0.4161820787  0.6365401     0.52442443
## 7  -1.049170739 -7.733000e-02 0.4389110496 -1.4669214     0.14239744
## 8   0.299911996 -7.434760e-03 0.0453950124  1.4425290     0.14915319
## 9  -0.069715708 -1.193043e-02 0.0725146217 -0.2145875     0.83008893
## 10  0.461968008 -1.113939e-02 0.0677608132  1.8174833     0.06914314
## 11  0.463429748 -2.555784e-02 0.1532012490  1.2492998     0.21155543
## 12  0.029152125 -1.127970e-04 0.0006937939  1.1110463     0.26654844
## 13  0.054035118 -9.786339e-03 0.0596116676  0.2613973     0.79378615
## 14  0.428842755 -8.120112e-03 0.0495453866  1.9631023     0.04963429
## 15 -0.385965422 -1.968131e-02 0.1186870693 -1.0632034     0.28768970
## 16 -0.297247626 -4.120777e-02 0.2430444451 -0.5193555     0.60351286
## 17  0.157546802 -3.333503e-02 0.1982252023  0.4287313     0.66811881
## 18 -0.189818565 -1.336505e-02 0.0811165102 -0.6195490     0.53555474
## 19 -0.534534616 -7.323236e-03 0.0447190911 -2.4930919     0.01266361
## 20  1.446348864 -1.343666e-01 0.7154961541  1.8687450     0.06165830
## 21  0.456998070 -1.085526e-02 0.0660514173  1.8204077     0.06869695
## 22  0.306693464 -2.980307e-02 0.1778701103  0.7978644     0.42494916
## 23  0.994861056 -5.218296e-02 0.3042533226  1.8982232     0.05766669
## 24 -0.164574074 -3.957018e-02 0.2337845073 -0.2585328     0.79599572
## 25  0.812562307 -2.818410e-02 0.1684884969  2.0482355     0.04053693
## 26 -0.006412117 -8.702126e-05 0.0005352660 -0.2733898     0.78455363
## 27  0.281201903 -2.985524e-02 0.1781719051  0.7369204     0.46117077
## 28  0.991035597 -5.535308e-02 0.3216573063  1.8449998     0.06503757
## 29 -0.202173061 -5.919937e-03 0.0362010007 -1.0314693     0.30232082
## 30  0.626016586 -5.090114e-02 0.2971810026  1.2417248     0.21433811
## 31  0.521281217 -2.429996e-02 0.1458491532  1.4285894     0.15312229
## 32 -0.248065662 -2.244067e-03 0.0137734356 -2.0945893     0.03620751
## 33  0.048159997 -1.515879e-02 0.0918359482  0.2089422     0.83449335
## 34  0.149467497 -1.621241e-03 0.0099569206  1.5141523     0.12998721
## 35  0.011321922 -1.019621e-05 0.0000627215  1.4308789     0.15246492
## attr(,"call")
## localmoran(x = resids, listw = WL)
## attr(,"class")
## [1] "localmoran" "matrix"     "array"     
## attr(,"quadr")
##         mean    median     pysal
## 1    Low-Low   Low-Low   Low-Low
## 2   Low-High  Low-High  Low-High
## 3   Low-High  Low-High   Low-Low
## 4  High-High High-High High-High
## 5    Low-Low   Low-Low   Low-Low
## 6  High-High High-High High-High
## 7   Low-High  Low-High  Low-High
## 8  High-High High-High High-High
## 9   High-Low High-High  High-Low
## 10 High-High High-High High-High
## 11   Low-Low   Low-Low   Low-Low
## 12   Low-Low   Low-Low   Low-Low
## 13   Low-Low  Low-High   Low-Low
## 14   Low-Low   Low-Low   Low-Low
## 15  High-Low  High-Low  High-Low
## 16  High-Low  High-Low  High-Low
## 17   Low-Low   Low-Low   Low-Low
## 18  High-Low  High-Low  High-Low
## 19  High-Low  High-Low  High-Low
## 20 High-High High-High High-High
## 21 High-High High-High High-High
## 22   Low-Low   Low-Low   Low-Low
## 23 High-High High-High High-High
## 24  High-Low  High-Low  High-Low
## 25 High-High High-High High-High
## 26  Low-High  Low-High  Low-High
## 27   Low-Low   Low-Low   Low-Low
## 28 High-High High-High High-High
## 29  High-Low  High-Low  High-Low
## 30   Low-Low   Low-Low   Low-Low
## 31   Low-Low   Low-Low   Low-Low
## 32  Low-High  Low-High  Low-High
## 33   Low-Low  Low-High   Low-Low
## 34   Low-Low   Low-Low   Low-Low
## 35   Low-Low  High-Low   Low-Low
Dengan taraf signifikansi 5%, diketahui bahwa terdapat autokorelasi spasial pada residual dari model regresi karena p-value < 0.05.

Geographically Weighted Regression (GWR)

Adaptive Kernel Model

bandwidth_optimal <- gwr.sel(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude, data_spasial$latitude),
  adapt = TRUE)
## Adaptive q: 0.381966 CV score: 87338505781 
## Adaptive q: 0.618034 CV score: 86084845360 
## Adaptive q: 0.763932 CV score: 87105024261 
## Adaptive q: 0.5824378 CV score: 85841638787 
## Adaptive q: 0.5058644 CV score: 86872922685 
## Adaptive q: 0.5813597 CV score: 85848989510 
## Adaptive q: 0.5910577 CV score: 85785843627 
## Adaptive q: 0.6013617 CV score: 85760260700 
## Adaptive q: 0.6020974 CV score: 85774905743 
## Adaptive q: 0.5968218 CV score: 85751345850 
## Adaptive q: 0.5946201 CV score: 85764266759 
## Adaptive q: 0.5982467 CV score: 85743148932 
## Adaptive q: 0.5986954 CV score: 85740593823 
## Adaptive q: 0.5997138 CV score: 85734841728 
## Adaptive q: 0.6003433 CV score: 85740036612 
## Adaptive q: 0.5995394 CV score: 85735822504 
## Adaptive q: 0.5997895 CV score: 85734417179 
## Adaptive q: 0.600001 CV score: 85733257070 
## Adaptive q: 0.6001317 CV score: 85735845446 
## Adaptive q: 0.5999202 CV score: 85733684072 
## Adaptive q: 0.6000509 CV score: 85734245564 
## Adaptive q: 0.600001 CV score: 85733257070
model_gwr1 <- gwr(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude, data_spasial$latitude),
  adapt = bandwidth_optimal,
  hatmatrix = TRUE,
  se.fit = TRUE)

coef_gwr <- model_gwr1$SDF@data[, c("x1", "x2", "x3", "x4")]

summary(model_gwr1)
##           Length Class                  Mode     
## SDF         35   SpatialPointsDataFrame S4       
## lhat      1225   -none-                 numeric  
## lm          11   -none-                 list     
## results     14   -none-                 list     
## bandwidth   35   -none-                 numeric  
## adapt        1   -none-                 numeric  
## hatmatrix    1   -none-                 logical  
## gweight      1   -none-                 character
## gTSS         1   -none-                 numeric  
## this.call    7   -none-                 call     
## fp.given     1   -none-                 logical  
## timings     12   -none-                 numeric
model_gwr1
## Call:
## gwr(formula = y ~ x1 + x2 + x3 + x4, data = data_spasial, coords = cbind(data_spasial$longitude, 
##     data_spasial$latitude), adapt = bandwidth_optimal, hatmatrix = TRUE, 
##     se.fit = TRUE)
## Kernel function: gwr.Gauss 
## Adaptive quantile: 0.600001 (about 21 of 35 data points)
## Summary of GWR coefficient estimates at data points:
##                     Min.     1st Qu.      Median     3rd Qu.        Max.
## X.Intercept. -642415.427 -596738.749 -510504.220 -379737.827 -293035.809
## x1              8689.713   10185.570   12523.116   13372.485   13994.476
## x2              7547.074    9467.097   10636.669   11720.864   12581.935
## x3             -1124.840    1015.517    9311.947   14487.159   23680.076
## x4             -1534.070    -821.946     -68.455    1285.765    1469.481
##                  Global
## X.Intercept. -487217.37
## x1             11472.98
## x2             10820.90
## x3              4274.32
## x4               488.41
## Number of data points: 35 
## Effective number of parameters (residual: 2traceS - traceS'S): 8.655994 
## Effective degrees of freedom (residual: 2traceS - traceS'S): 26.34401 
## Sigma (residual: 2traceS - traceS'S): 44870.84 
## Effective number of parameters (model: traceS): 7.169543 
## Effective degrees of freedom (model: traceS): 27.83046 
## Sigma (model: traceS): 43656.1 
## Sigma (ML): 38928.8 
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 861.3293 
## AIC (GWR p. 96, eq. 4.22): 846.3595 
## Residual sum of squares: 53040811982 
## Quasi-global R2: 0.6645492
Dengan adaptive quantile 0.600001 dan CV score 85733257070, model GWR efektif dalam menjelaskan 21 titik data (AIC 846.3595 Quasi-global R² 0.6645492).

Fixed Kernel Model

bw1 <- gwr.sel(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude,data_spasial$latitude))
## Bandwidth: 1.114881 CV score: 88816064524 
## Bandwidth: 1.802115 CV score: 88874898868 
## Bandwidth: 0.6901475 CV score: 93355239861 
## Bandwidth: 1.45408 CV score: 88713224675 
## Bandwidth: 1.42018 CV score: 88701837029 
## Bandwidth: 1.356895 CV score: 88686932022 
## Bandwidth: 1.264454 CV score: 88689065437 
## Bandwidth: 1.317623 CV score: 88683504214 
## Bandwidth: 1.316234 CV score: 88683486429 
## Bandwidth: 1.313435 CV score: 88683474301 
## Bandwidth: 1.313762 CV score: 88683474071 
## Bandwidth: 1.313803 CV score: 88683474073 
## Bandwidth: 1.313722 CV score: 88683474075 
## Bandwidth: 1.313762 CV score: 88683474071
model_gwr2 <- gwr(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude, data_spasial$latitude),
  bandwidth = bw1,
  hatmatrix = TRUE,
  se.fit = TRUE)

coef_gwr2 <- model_gwr2$SDF@data[, c("x1", "x2", "x3", "x4")]

summary(model_gwr2)
##           Length Class                  Mode     
## SDF         35   SpatialPointsDataFrame S4       
## lhat      1225   -none-                 numeric  
## lm          11   -none-                 list     
## results     14   -none-                 list     
## bandwidth    1   -none-                 numeric  
## adapt        0   -none-                 NULL     
## hatmatrix    1   -none-                 logical  
## gweight      1   -none-                 character
## gTSS         1   -none-                 numeric  
## this.call    7   -none-                 call     
## fp.given     1   -none-                 logical  
## timings     12   -none-                 numeric
model_gwr2
## Call:
## gwr(formula = y ~ x1 + x2 + x3 + x4, data = data_spasial, coords = cbind(data_spasial$longitude, 
##     data_spasial$latitude), bandwidth = bw1, hatmatrix = TRUE, 
##     se.fit = TRUE)
## Kernel function: gwr.Gauss 
## Fixed bandwidth: 1.313762 
## Summary of GWR coefficient estimates at data points:
##                    Min.    1st Qu.     Median    3rd Qu.       Max.     Global
## X.Intercept. -636210.82 -564750.85 -494367.35 -422970.94 -340589.99 -487217.37
## x1              8823.19   10765.38   12011.65   12910.80   13549.48   11472.98
## x2              9082.05   10044.11   10823.89   11602.38   12209.52   10820.90
## x3             -1184.56    2173.85    6393.37   10140.99   14155.17    4274.32
## x4              -873.30    -392.37     214.17     885.68    1797.39     488.41
## Number of data points: 35 
## Effective number of parameters (residual: 2traceS - traceS'S): 7.360732 
## Effective degrees of freedom (residual: 2traceS - traceS'S): 27.63927 
## Sigma (residual: 2traceS - traceS'S): 45293.43 
## Effective number of parameters (model: traceS): 6.334603 
## Effective degrees of freedom (model: traceS): 28.6654 
## Sigma (model: traceS): 44475.36 
## Sigma (ML): 40249.87 
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 860.7803 
## AIC (GWR p. 96, eq. 4.22): 847.8606 
## Residual sum of squares: 56701816456 
## Quasi-global R2: 0.6413955
Dengan fixed bandwidth 1.313762 dan CV score 88683474071, model GWR efektif dalam menjelaskan 35 titik data (AIC 847.8606 Quasi-global R² 0.6413955).

Bi-square

bw2 <- gwr.sel(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude,data_spasial$latitude),
  gweight = gwr.bisquare)
## Bandwidth: 1.114881 CV score: 110797595123 
## Bandwidth: 1.802115 CV score: 94429533715 
## Bandwidth: 2.226848 CV score: 90824090745 
## Bandwidth: 2.322375 CV score: 90402172510 
## Bandwidth: 2.548387 CV score: 89801240658 
## Bandwidth: 2.67855 CV score: 89647200259 
## Bandwidth: 2.756314 CV score: 89584021082 
## Bandwidth: 2.816576 CV score: 89534987367 
## Bandwidth: 2.85382 CV score: 89505603846 
## Bandwidth: 2.876838 CV score: 89487962388 
## Bandwidth: 2.891064 CV score: 89477286050 
## Bandwidth: 2.899856 CV score: 89470779724 
## Bandwidth: 2.90529 CV score: 89466794849 
## Bandwidth: 2.908648 CV score: 89464346147 
## Bandwidth: 2.910723 CV score: 89462838199 
## Bandwidth: 2.912006 CV score: 89461908323 
## Bandwidth: 2.912799 CV score: 89461334429 
## Bandwidth: 2.913289 CV score: 89460980049 
## Bandwidth: 2.913592 CV score: 89460761147 
## Bandwidth: 2.913779 CV score: 89460625903 
## Bandwidth: 2.913895 CV score: 89460542335 
## Bandwidth: 2.913966 CV score: 89460490693 
## Bandwidth: 2.91401 CV score: 89460458780 
## Bandwidth: 2.91401 CV score: 89460458780
model_gwr3 <- gwr(
  y ~ x1 + x2 + x3 + x4,
  data = data_spasial,
  coords = cbind(data_spasial$longitude, data_spasial$latitude),
  bandwidth = bw2,
  gweight = gwr.bisquare,
  hatmatrix = TRUE,
  se.fit = TRUE)

coef_gwr3 <- model_gwr3$SDF@data[, c("x1", "x2", "x3", "x4")]

summary(model_gwr3)
##           Length Class                  Mode     
## SDF         35   SpatialPointsDataFrame S4       
## lhat      1225   -none-                 numeric  
## lm          11   -none-                 list     
## results     14   -none-                 list     
## bandwidth    1   -none-                 numeric  
## adapt        0   -none-                 NULL     
## hatmatrix    1   -none-                 logical  
## gweight      1   -none-                 character
## gTSS         1   -none-                 numeric  
## this.call    8   -none-                 call     
## fp.given     1   -none-                 logical  
## timings     12   -none-                 numeric
model_gwr3
## Call:
## gwr(formula = y ~ x1 + x2 + x3 + x4, data = data_spasial, coords = cbind(data_spasial$longitude, 
##     data_spasial$latitude), bandwidth = bw2, gweight = gwr.bisquare, 
##     hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.bisquare 
## Fixed bandwidth: 2.91401 
## Summary of GWR coefficient estimates at data points:
##                    Min.    1st Qu.     Median    3rd Qu.       Max.     Global
## X.Intercept. -664578.52 -562084.32 -493575.57 -424303.72 -324093.64 -487217.37
## x1              8393.64   10641.63   11960.47   12983.60   14240.27   11472.98
## x2              8925.83   10166.15   10841.47   11593.32   12305.73   10820.90
## x3             -1192.41    2422.28    6170.38    9682.00   15624.55    4274.32
## x4             -1151.86    -401.69     226.99     841.74    1915.51     488.41
## Number of data points: 35 
## Effective number of parameters (residual: 2traceS - traceS'S): 7.206494 
## Effective degrees of freedom (residual: 2traceS - traceS'S): 27.79351 
## Sigma (residual: 2traceS - traceS'S): 45438.54 
## Effective number of parameters (model: traceS): 6.266229 
## Effective degrees of freedom (model: traceS): 28.73377 
## Sigma (model: traceS): 44688.91 
## Sigma (ML): 40491.33 
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 860.9707 
## AIC (GWR p. 96, eq. 4.22): 848.2109 
## Residual sum of squares: 57384175114 
## Quasi-global R2: 0.63708
Dengan fixed bandwidth 2.91401 dan CV score 89460458780, model GWR efektif dalam menjelaskan 35 titik data (AIC 848.2109 Quasi-global R² 0.63708).

Comparing GWR Model

data.frame("MODEL" = c("Adaptive Gaussian", "Fixed Gaussian", "Bi-Square"), 
           "AIC" = c(model_gwr1$results$AICh, model_gwr2$results$AICh, model_gwr3$results$AICh),
           "R2" = c(0.6645492,0.6413955,0.63708))%>% arrange(AIC)
##               MODEL      AIC        R2
## 1 Adaptive Gaussian 846.3595 0.6645492
## 2    Fixed Gaussian 847.8606 0.6413955
## 3         Bi-Square 848.2109 0.6370800
Dapat disimpulkan bahwa pemodelan GWR dengan kernel bandwidth terbaik adalah Fixed Gaussian karena memiliki nilai AIC terendah dan R2 tertinggi.

GWR Model

# Menggunakan GWR dengan AIC terkecil dan R2 terbesar
model_gwr = model_gwr1
results <-as.data.frame(model_gwr$SDF)
head(results)
##      sum.w X.Intercept.        x1        x2         x3         x4
## 1 21.11540    -642415.4 13994.476 12107.310   772.8846  1430.3832
## 2 23.20282    -609590.0 13331.286 11651.287  1164.2501  1443.5054
## 3 22.81590    -590773.7 13413.683 12300.019   866.7842   769.0229
## 4 21.28866    -340066.3  9044.685  8921.756 15039.6979  -937.1239
## 5 21.91916    -383349.6 10834.337  8933.491 18243.0763 -1073.6904
## 6 23.64115    -595379.3 13104.700 11790.441   346.0111  1325.1237
##   X.Intercept._se    x1_se    x2_se    x3_se    x4_se     gwr.e     pred
## 1        183378.5 5502.015 4934.367 22389.06 1900.742 -38984.19 419030.2
## 2        176406.4 5237.872 4816.832 21984.57 1828.188  -5010.27 484037.3
## 3        174109.8 5349.301 4809.355 21894.04 1768.548 -71397.74 450255.7
## 4        180949.9 5583.321 5192.311 23643.19 1889.079  10322.88 414812.1
## 5        180081.8 5881.205 5250.351 23804.19 1885.640 -57043.16 477382.2
## 6        173920.9 5163.243 4781.393 21838.95 1790.301  69638.34 443700.7
##    pred.se   localR2 X.Intercept._se_EDF x1_se_EDF x2_se_EDF x3_se_EDF
## 1 16693.30 0.6874768            188481.0  5655.109  5071.666  23012.04
## 2 11937.24 0.6566978            181315.0  5383.616  4950.861  22596.29
## 3 14498.34 0.6998012            178954.5  5498.146  4943.176  22503.25
## 4 19790.69 0.7034268            185984.9  5738.678  5336.788  24301.06
## 5 20436.55 0.7468688            185092.7  6044.851  5396.443  24466.55
## 6 23253.65 0.6503113            178760.3  5306.912  4914.436  22446.62
##   x4_se_EDF pred.se.1 longitude  latitude
## 1  1953.631  17157.80  109.6571 -7.351319
## 2  1879.058  12269.40  109.1758 -7.455506
## 3  1817.758  14901.76  109.8616 -7.020996
## 4  1941.643  20341.37  111.3874 -7.074421
## 5  1938.108  21005.20  110.6525 -7.416504
## 6  1840.116  23900.69  108.9295 -7.062734

Mapping the Local R2

gwr.map <- cbind(Jateng, as.matrix(results)) 
qtm(gwr.map, fill = "localR2") 

Plotting residuals model GWR

gwr.map_sp <- as(gwr.map, "Spatial")
spplot(gwr.map_sp, "localR2")

Comparing OLS and GWR

data.frame("MODEL" = c("GWR","OLS"), 
           "AIC" = c(model_gwr$results$AICh,AIC(model_ols)),
           "R2" = c(0.6645492,0.5423))%>% arrange(AIC)
##   MODEL      AIC        R2
## 1   GWR 846.3595 0.6645492
## 2   OLS 857.6854 0.5423000
#Attach coefficients to original dataframe 
data_spasial$coefx1<-results$x1 
data_spasial$coefx2<-results$x2 
data_spasial$coefx3<-results$x3 
data_spasial$coefx4<-results$x4 
Dapat disimpulkan bahwa pemodelan terbaik adalah GWR Fixed Gaussian karena memiliki nilai AIC terendah dan R² tertinggi, yaitu 68.3% variasi dalam model dijelaskan oleh variabel independennya.

Global Significance Test

LMZ.F3GWR.test(model_gwr)
## 
## Leung et al. (2000) F(3) test
## 
##             F statistic Numerator d.f. Denominator d.f.    Pr(>)   
## (Intercept)     3.73018       11.94610           28.698 0.001889 **
## x1              0.68844       15.83113           28.698 0.780481   
## x2              0.61915        8.79171           28.698 0.767626   
## x3              0.92534        9.19236           28.698 0.519485   
## x4              2.38428        7.48818           28.698 0.044373 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Karena p-value untuk variabel x4 < 0.05, maka terdapat perbedaan signifikan pada x4 antara lokasi satu dengan lainnya.

Local Significance Test

t1 = model_gwr$SDF$x1 / model_gwr$SDF$x1_se
t2 = model_gwr$SDF$x2 / model_gwr$SDF$x2_se
t3 = model_gwr$SDF$x3 / model_gwr$SDF$x3_se
t4 = model_gwr$SDF$x4 / model_gwr$SDF$x4_se

pvalue1 = 2*pt(q=abs(t1), df=30, lower.tail=FALSE) 
pvalue2 = 2*pt(q=abs(t2), df=30, lower.tail=FALSE) 
pvalue3 = 2*pt(q=abs(t3), df=30, lower.tail=FALSE) 
pvalue4 = 2*pt(q=abs(t4), df=30, lower.tail=FALSE) 

Variabel x1

df <- data.frame(id = id, t1 = t1, pvalue1 = pvalue1) %>%
  mutate(ket = ifelse(pvalue1 < 0.05, "Reject", "Not Reject"))
df
##    id       t1    pvalue1        ket
## 1   1 2.543518 0.01636589     Reject
## 2   2 2.545172 0.01630257     Reject
## 3   3 2.507558 0.01779990     Reject
## 4   4 1.619947 0.11570961 Not Reject
## 5   5 1.842197 0.07535078 Not Reject
## 6   6 2.538075 0.01657588     Reject
## 7   7 2.531973 0.01681423     Reject
## 8   8 1.719448 0.09583610 Not Reject
## 9   9 1.346576 0.18820280 Not Reject
## 10 10 1.881068 0.06970018 Not Reject
## 11 11 1.747473 0.09078583 Not Reject
## 12 12 2.572452 0.01529028     Reject
## 13 13 2.291448 0.02913213     Reject
## 14 14 2.104573 0.04381250     Reject
## 15 15 2.363630 0.02477095     Reject
## 16 16 2.541532 0.01644223     Reject
## 17 17 2.022388 0.05213087 Not Reject
## 18 18 2.021700 0.05220601 Not Reject
## 19 19 1.792552 0.08313481 Not Reject
## 20 20 2.546782 0.01624114     Reject
## 21 21 1.703405 0.09883203 Not Reject
## 22 22 2.334702 0.02644214     Reject
## 23 23 1.722633 0.09525050 Not Reject
## 24 24 2.550486 0.01610066     Reject
## 25 25 2.561516 0.01568893     Reject
## 26 26 2.557777 0.01582738     Reject
## 27 27 2.557109 0.01585226     Reject
## 28 28 1.715003 0.09665847 Not Reject
## 29 29 2.002043 0.05439314 Not Reject
## 30 30 1.520749 0.13879379 Not Reject
## 31 31 1.836160 0.07626276 Not Reject
## 32 32 2.550397 0.01610405     Reject
## 33 33 2.379923 0.02387263     Reject
## 34 34 1.946300 0.06103776 Not Reject
## 35 35 2.552506 0.01602453     Reject
Variabel x1 berpengaruh signifikan di Banjarnegara, Banyumas, Batang, Brebes, Cilacap, Kebumen, Kendal, Klaten, Kota Magelang, Kota Pekalongan, Kota Tegal, Magelang, Pekalongan, Pemalang, Purbalingga, Purworejo di Provinsi Jawa Tengah.

Variabel x2

df2 <- data.frame(id = id, t2 = t2, pvalue2 = pvalue2) %>%
  mutate(ket = ifelse(pvalue2 < 0.05, "Reject", "Not Reject"))
df2
##    id       t2    pvalue2        ket
## 1   1 2.453671 0.02016625     Reject
## 2   2 2.418869 0.02184419     Reject
## 3   3 2.557520 0.01583696     Reject
## 4   4 1.718263 0.09605475 Not Reject
## 5   5 1.701503 0.09919225 Not Reject
## 6   6 2.465901 0.01960517     Reject
## 7   7 2.416127 0.02198169     Reject
## 8   8 1.888879 0.06860981 Not Reject
## 9   9 1.318752 0.19723104 Not Reject
## 10 10 2.128294 0.04164126     Reject
## 11 11 1.672271 0.10487031 Not Reject
## 12 12 2.366917 0.02458727     Reject
## 13 13 2.424099 0.02158417     Reject
## 14 14 1.951054 0.06044464 Not Reject
## 15 15 2.197375 0.03585341     Reject
## 16 16 2.586885 0.01477855     Reject
## 17 17 1.934780 0.06249577 Not Reject
## 18 18 2.165001 0.03846977     Reject
## 19 19 1.668690 0.10558413 Not Reject
## 20 20 2.502428 0.01801366     Reject
## 21 21 1.890885 0.06833215 Not Reject
## 22 22 2.149231 0.03980493     Reject
## 23 23 1.905786 0.06630016 Not Reject
## 24 24 2.547085 0.01622962     Reject
## 25 25 2.513092 0.01757196     Reject
## 26 26 2.451401 0.02027198     Reject
## 27 27 2.283678 0.02964071     Reject
## 28 28 1.872949 0.07084925 Not Reject
## 29 29 1.953260 0.06017115 Not Reject
## 30 30 1.431189 0.16271173 Not Reject
## 31 31 1.698996 0.09966904 Not Reject
## 32 32 2.486285 0.01870178     Reject
## 33 33 2.360322 0.02495708     Reject
## 34 34 1.849684 0.07423280 Not Reject
## 35 35 2.416214 0.02197732     Reject
Variabel x2 berpengaruh signifikan di Banjarnegara, Banyumas, Batang, Brebes, Cilacap, Jepara, Kebumen, Kendal, Kota Magelang, Kota Pekalongan, Kota Semarang, Kota Tegal, Magelang, Pekalongan, Pemalang, Purbalingga, Purworejo di Provinsi Jawa Tengah.

Variabel x3

df3 <- data.frame(id = id, t3 = t3, pvalue3 = pvalue3) %>%
  mutate(ket = ifelse(pvalue3 < 0.05, "Reject", "Not Reject"))
df3
##    id           t3   pvalue3        ket
## 1   1  0.034520630 0.9726907 Not Reject
## 2   2  0.052957603 0.9581168 Not Reject
## 3   3  0.039589961 0.9686823 Not Reject
## 4   4  0.636111281 0.5295274 Not Reject
## 5   5  0.766380843 0.4494385 Not Reject
## 6   6  0.015843761 0.9874639 Not Reject
## 7   7  0.062659476 0.9504533 Not Reject
## 8   8  0.626273645 0.5358720 Not Reject
## 9   9  0.915407797 0.3672775 Not Reject
## 10 10  0.418783622 0.6783575 Not Reject
## 11 11  0.741920611 0.4639047 Not Reject
## 12 12  0.139439086 0.8900351 Not Reject
## 13 13  0.239605319 0.8122658 Not Reject
## 14 14  0.603826499 0.5504980 Not Reject
## 15 15  0.414131569 0.6817238 Not Reject
## 16 16 -0.050897261 0.9597448 Not Reject
## 17 17  0.628671467 0.5343219 Not Reject
## 18 18  0.467125425 0.6437852 Not Reject
## 19 19  0.767939186 0.4485261 Not Reject
## 20 20 -0.012819955 0.9898563 Not Reject
## 21 21  0.575414080 0.5693025 Not Reject
## 22 22  0.448144884 0.6572685 Not Reject
## 23 23  0.543150422 0.5910395 Not Reject
## 24 24 -0.042096459 0.9667007 Not Reject
## 25 25 -0.026448999 0.9790744 Not Reject
## 26 26  0.022587270 0.9821291 Not Reject
## 27 27  0.285162224 0.7774794 Not Reject
## 28 28  0.533200869 0.5978229 Not Reject
## 29 29  0.617517764 0.5415526 Not Reject
## 30 30  0.873550296 0.3893023 Not Reject
## 31 31  0.748016766 0.4602741 Not Reject
## 32 32 -0.005397868 0.9957289 Not Reject
## 33 33  0.287024860 0.7760665 Not Reject
## 34 34  0.641743032 0.5259135 Not Reject
## 35 35  0.179949755 0.8584020 Not Reject
Variabel x3 tidak berpengaruh signifikan di setiap Kabupaten/Kota di Provinsi Jawa Tengah.

Variabel x4

df4 <- data.frame(id = id, t4 = t4, pvalue4 = pvalue4) %>%
  mutate(ket = ifelse(pvalue3 < 0.05, "Reject", "Not Reject"))
df4
##    id           t4   pvalue4        ket
## 1   1  0.752539241 0.4575915 Not Reject
## 2   2  0.789582394 0.4359683 Not Reject
## 3   3  0.434832885 0.6667958 Not Reject
## 4   4 -0.496074399 0.6234570 Not Reject
## 5   5 -0.569403857 0.5733213 Not Reject
## 6   6  0.740168216 0.4649515 Not Reject
## 7   7  0.760063910 0.4531485 Not Reject
## 8   8 -0.530515366 0.5996601 Not Reject
## 9   9 -0.757203736 0.4548343 Not Reject
## 10 10 -0.302799365 0.7641325 Not Reject
## 11 11 -0.546086201 0.5890451 Not Reject
## 12 12  0.704510975 0.4865482 Not Reject
## 13 13  0.002837640 0.9977547 Not Reject
## 14 14 -0.305454433 0.7621295 Not Reject
## 15 15  0.002585021 0.9979546 Not Reject
## 16 16  0.599288833 0.5534796 Not Reject
## 17 17 -0.420826615 0.6768812 Not Reject
## 18 18 -0.347864501 0.7303710 Not Reject
## 19 19 -0.561804204 0.5784229 Not Reject
## 20 20  0.716292981 0.4793497 Not Reject
## 21 21 -0.457401234 0.6506780 Not Reject
## 22 22 -0.038273066 0.9697235 Not Reject
## 23 23 -0.414393191 0.6815343 Not Reject
## 24 24  0.702416654 0.4878342 Not Reject
## 25 25  0.745249936 0.4619198 Not Reject
## 26 26  0.793066386 0.4339668 Not Reject
## 27 27  0.401669587 0.6907744 Not Reject
## 28 28 -0.390438327 0.6989712 Not Reject
## 29 29 -0.432300147 0.6686149 Not Reject
## 30 30 -0.708781497 0.4839319 Not Reject
## 31 31 -0.518569513 0.6078651 Not Reject
## 32 32  0.749073613 0.4596464 Not Reject
## 33 33  0.084182063 0.9334709 Not Reject
## 34 34 -0.375352338 0.7100393 Not Reject
## 35 35  0.434358456 0.6671364 Not Reject
Variabel x4 tidak berpengaruh signifikan di setiap Kabupaten/Kota di Provinsi Jawa Tengah.

Spatial Distribution of Each Variable

map1 <- tm_shape(gwr.map) + 
  tm_fill("x1", 
          n = 5, 
          style = "quantile", 
          title = "x1 Coefficient") + 
  tm_borders(col = "black", lwd = 0.5) + 
  tm_layout(frame = FALSE, 
            legend.text.size = 0.5, 
            legend.title.size = 0.6)

map2 <- tm_shape(gwr.map) + 
  tm_fill("x2", 
          n = 5, 
          style = "quantile", 
          title = "x2 Coefficient") + 
  tm_borders(col = "black", lwd = 0.5) + 
  tm_layout(frame = FALSE, 
            legend.text.size = 0.5, 
            legend.title.size = 0.6)

map3 <- tm_shape(gwr.map) + 
  tm_fill("x3", 
          n = 5, 
          style = "quantile", 
          title = "x3 Coefficient") + 
  tm_borders(col = "black", lwd = 0.5) + 
  tm_layout(frame = FALSE, 
            legend.text.size = 0.5, 
            legend.title.size = 0.6)

map4 <- tm_shape(gwr.map) + 
  tm_fill("x4", 
          n = 5, 
          style = "quantile", 
          title = "x4 Coefficient") + 
  tm_borders(col = "black", lwd = 0.5) + 
  tm_layout(frame = FALSE, 
            legend.text.size = 0.5, 
            legend.title.size = 0.6)

tmap_arrange(map1, map2, map3, map4, ncol = 2, nrow = 2)

Membuat plot -spplot dan -ggplot2

Ekstrak koefisien lokal untuk variabel x1, x2, x3

# Ekstrak koefisien lokal untuk variabel x1, x2, x3, x4
coef_x1 <- model_gwr$SDF$x1
coef_x2 <- model_gwr$SDF$x2
coef_x3 <- model_gwr$SDF$x3
coef_x4 <- model_gwr$SDF$x4

# Konversi SpatialPointsDataFrame menjadi data frame biasa
gwr_df <- as.data.frame(model_gwr$SDF)

Plot koefisien regresi lokal untuk variabel x1

# Koefisien Lokal x1 dengan ssplot
spplot(model_gwr$SDF, "x1", main = "Koefisien GWR untuk x1 (ssplot)", col.regions = terrain.colors(100))

# Koefisien Lokal x1 dengan ggplot2
ggplot(gwr_df, aes(x = data_spasial$longitude, y = data_spasial$latitude)) +
  geom_point(aes(color = x1), size = 2) +
  scale_color_viridis_c() +
  labs(title = "Koefisien GWR untuk x1 (ggplot2)", color = "Koefisien") +
  theme_minimal()

Plot koefisien regresi lokal untuk variabel x2

# Koefisien Lokal x2 dengan ssplot
spplot(model_gwr$SDF, "x2", main = "Koefisien GWR untuk x2 (ssplot)", col.regions = terrain.colors(100))

# Koefisien Lokal x2 dengan ggplot2
ggplot(gwr_df, aes(x = data_spasial$longitude, y = data_spasial$latitude)) +
  geom_point(aes(color = x2), size = 2) +
  scale_color_viridis_c() +
  labs(title = "Koefisien GWR untuk x2 (ggplot2)", color = "Koefisien") +
  theme_minimal()

Plot koefisien regresi lokal untuk variabel x3

# Koefisien Lokal x3 dengan ssplot
spplot(model_gwr$SDF, "x3", main = "Koefisien GWR untuk x3 (ssplot)", col.regions = terrain.colors(100))

# Koefisien Lokal x3 dengan ggplot2
ggplot(gwr_df, aes(x = data_spasial$longitude, y = data_spasial$latitude)) +
  geom_point(aes(color = x3), size = 2) +
  scale_color_viridis_c() +
  labs(title = "Koefisien GWR untuk x3 (ggplot)", color = "Koefisien") +
  theme_minimal()

Plot koefisien regresi lokal untuk variabel x4

# Koefisien Lokal x4 dengan ssplot
spplot(model_gwr$SDF, "x4", main = "Koefisien GWR untuk x4 (ssplot)", col.regions = terrain.colors(100))

# Koefisien Lokal x4 dengan ggplot2
ggplot(gwr_df, aes(x = data_spasial$longitude, y = data_spasial$latitude)) +
  geom_point(aes(color = x4), size = 2) +
  scale_color_viridis_c() +
  labs(title = "Koefisien GWR untuk x3 (ggplot)", color = "Koefisien") +
  theme_minimal()