library(spdep)
## Loading required package: spData
## Warning: package 'spData' was built under R version 4.4.1
## 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.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE
library(raster)
## Loading required package: sp
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.ZeroPolicyOption
library(sf)
library(readr)
library(spgwr)
## Warning: package 'spgwr' was built under R version 4.4.1
## NOTE: This package does not constitute approval of GWR
## as a method of spatial analysis; see example(gwr)
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ purrr 1.0.2
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::expand() masks Matrix::expand()
## ✖ tidyr::extract() masks raster::extract()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tidyr::pack() masks Matrix::pack()
## ✖ dplyr::select() masks raster::select()
## ✖ tidyr::unpack() masks Matrix::unpack()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sp)
library(tmap)
## Breaking News: tmap 3.x is retiring. Please test v4, e.g. with
## remotes::install_github('r-tmap/tmap')
library(tmaptools)
library(grid)
library(gridExtra)
##
## Attaching package: 'gridExtra'
##
## The following object is masked from 'package:dplyr':
##
## combine
library(dplyr)
library(GWmodel)
## Warning: package 'GWmodel' was built under R version 4.4.1
## Loading required package: robustbase
## Loading required package: Rcpp
## Welcome to GWmodel version 2.4-2.
Indo_Kec<-readRDS('gadm36_IDN_2_sp.rds')
Jabar<-Indo_Kec[Indo_Kec$NAME_1 == "Jawa Barat",]
plot(Jabar)
text(coordinates(Jabar), labels = Jabar$NAME_2, cex = 0.7, col = "black")
setwd("C:\\Users\\lenovo\\OneDrive\\Documents\\Sem 5\\Spasial")
Pneumonia<-read.csv("Pneumonia 2022.csv", sep=";")
Pneumonia
## Kabupaten.Kota Pneumonia.pada.Balita Imunisasi Gizi.Buruk
## 1 KABUPATEN BANDUNG 5345 61605 221317
## 2 KABUPATEN BANDUNG BARAT 2409 27620 223235
## 3 KOTA BANJAR 1426 70329 213501
## 4 KABUPATEN BEKASI 3605 93549 367107
## 5 KABUPATEN BOGOR 3104 18319 51344
## 6 KABUPATEN CIAMIS 4176 30907 155701
## 7 KABUPATEN CIANJUR 5737 42509 143607
## 8 KOTA CIMAHI 5297 45081 171481
## 9 KABUPATEN CIREBON 4652 24583 107748
## 10 KOTA DEPOK 4951 41966 136465
## 11 KABUPATEN GARUT 1529 18426 64789
## 12 KABUPATEN INDRAMAYU 2590 19641 73783
## 13 KABUPATEN KARAWANG 200 5975 20932
## 14 KOTA BANDUNG 3557 16391 63811
## 15 KOTA BEKASI 3205 29978 88722
## 16 KOTA BOGOR 4882 39056 169577
## 17 KOTA CIREBON 2867 18781 64757
## 18 KOTA SUKABUMI 1772 29105 105675
## 19 KOTA TASIKMALAYA 3238 36030 93176
## 20 KABUPATEN KUNINGAN 498 2617 10177
## 21 KABUPATEN MAJALENGKA 1489 44376 141352
## 22 KABUPATEN PURWAKARTA 1833 16188 65944
## 23 KABUPATEN SUBANG 947 9463 28370
## 24 KABUPATEN SUKABUMI 1543 4997 16815
## 25 KABUPATEN SUMEDANG 824 40175 99205
## 26 KABUPATEN TASIKMALAYA 1042 5569 19583
## 27 KABUPATEN PANGANDARAN 1118 9550 41040
Jabar$Pneumonia.pada.Balita <- Pneumonia$Pneumonia.pada.Balita
Jabar$Imunisasi <- Pneumonia$Imunisasi
Jabar$Gizi.Buruk <- Pneumonia$Gizi.Buruk
y <- Jabar$Pneumonia.pada.Balita
x1 <- Jabar$Imunisasi
x2 <- Jabar$Gizi.Buruk
jabarlm <- lm(y ~ (x1+x2))
summary(jabarlm)
##
## Call:
## lm(formula = y ~ (x1 + x2))
##
## Residuals:
## Min 1Q Median 3Q Max
## -2391.6 -1068.8 -222.7 1139.4 2647.8
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.590e+03 4.949e+02 3.213 0.00372 **
## x1 -1.670e-04 3.529e-02 -0.005 0.99626
## x2 1.049e-02 9.231e-03 1.136 0.26711
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1479 on 24 degrees of freedom
## Multiple R-squared: 0.266, Adjusted R-squared: 0.2048
## F-statistic: 4.349 on 2 and 24 DF, p-value: 0.02445
Model diagnostics
par(mfrow=c(2,2))
plot(jabarlm)
Multicollinearity
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
## The following object is masked from 'package:purrr':
##
## some
vif(jabarlm)
## x1 x2
## 6.790208 6.790208
Normality
resids <- residuals(jabarlm)
shapiro.test(resids)
##
## Shapiro-Wilk normality test
##
## data: resids
## W = 0.95014, p-value = 0.2161
Heteroskedasticity
library(lmtest)
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
bptest(jabarlm)
##
## studentized Breusch-Pagan test
##
## data: jabarlm
## BP = 9.3666, df = 2, p-value = 0.009248
Non-autocorrelation
durbinWatsonTest(jabarlm)
## lag Autocorrelation D-W Statistic p-value
## 1 0.4284018 1.087989 0.012
## Alternative hypothesis: rho != 0
Jabar$id<-c(1:27)
Pneumonia$id<-c(1:27)
Pneumonia$x<-coordinates(Jabar)[,1]
Pneumonia$y<-coordinates(Jabar)[,2]
Jabar_sf<-st_as_sf(Jabar)
Jabar_merged <- Jabar_sf %>%
left_join(Pneumonia, by = "id")
#Plot
ggplot() +
geom_sf(data = Jabar_merged, aes(fill = Pneumonia.pada.Balita.y), color = NA) +
theme_bw() +
scale_fill_gradient(low = "yellow", high = "red") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "right",
axis.text.x = element_blank(), # Menghilangkan label sumbu x
axis.text.y = element_blank()) + # Menghilangkan label sumbu y
labs(title = "",
fill = "Pneumonia") +
geom_sf_text(data = Jabar_merged, aes(label = Jabar$NAME_2), size = 3, color = "black")
## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Mapping by Interval
summary(Pneumonia$Pneumonia.pada.Balita)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 200 1458 2590 2735 3890 5737
breaks <- c(-Inf, 1458, 2735, 3890,Inf)
# Define labels for each interval
labels <- c("Very Low", "Low", "High",
"Very High")
# Create a new column with discretized Pneumonia
Jabar_merged$Pneumonia_Discrete <-
cut(Jabar_merged$Pneumonia.pada.Balita.y, breaks =
breaks, labels = labels, right = TRUE)
ggplot() +
geom_sf(data = Jabar_merged, aes(fill = Pneumonia_Discrete), color = NA) +
theme_bw() +
scale_fill_manual(values = c("Very Low" = "yellow",
"Low" = "orange",
"High" = "red",
"Very High" = "red3")) +
labs(fill = "Pneumonia") +
theme(legend.position = "right",
axis.text.x = element_blank(), # Menghapus label sumbu x
axis.text.y = element_blank()) + # Menghapus label sumbu y
geom_sf_text(data = Jabar_merged, aes(label = NAME_2), size = 3, color = "black") +
labs(title = "", fill = "Pneumonia")
## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
ID <- c(1:27)
CoordK<-coordinates(Jabar)
row.names(Jabar) = as.character(1:nrow(Jabar))
W <- poly2nb(Jabar, row.names=ID,
queen=FALSE) ; W
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 108
## Percentage nonzero weights: 14.81481
## Average number of links: 4
WB <- nb2mat(W, style='B', zero.policy = TRUE) #menyajikan dalam bentuk matrix biner "B"
WL_Rook <- nb2listw(W); WL_Rook #List neighbours
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 108
## Percentage nonzero weights: 14.81481
## Average number of links: 4
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 27 729 27 16.57685 119.6581
plot(Jabar, axes=T, col="white")
text(CoordK[, 1], CoordK[, 2],
label = Jabar$NAME_2, col = "black", cex = 0.8, pos = 1.5)
points(CoordK[,1], CoordK[,2], pch=19,
cex=0.5,col="blue")
plot(W, coordinates(Jabar), add=TRUE, col = "blue", lwd=1)
#global moran Rook
MI_Rook <- moran.test(Pneumonia$Pneumonia.pada.Balita, WL_Rook)
moran.plot(Pneumonia$Pneumonia.pada.Balita, WL_Rook)
#local moran Rook
locm<-localmoran(Pneumonia$Pneumonia.pada.Balita, WL_Rook)
head(locm)
## Ii E.Ii Var.Ii Z.Ii Pr(z != E(Ii))
## 1 0.16299652 -0.098978426 0.261429776 0.5123683 0.60839328
## 2 -0.08243156 -0.001540619 0.004509257 -1.2046126 0.22835288
## 3 -0.71238608 -0.024877497 0.654982400 -0.8494999 0.39560321
## 4 -0.18571886 -0.011003221 0.090104204 -0.5820486 0.56053393
## 5 0.07111622 -0.001981461 0.004805411 1.0544801 0.29166318
## 6 -0.65109112 -0.030177098 0.132752573 -1.7041585 0.08835147
W <- poly2nb(Jabar, row.names=ID,
queen=TRUE) #Mendapatkan W
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 1 0 0 0
## 2 1 0 0 0 0
## 3 0 0 0 0 0
## 4 0 0 0 0 1
## 5 0 0 0 1 0
WL_Queen<-nb2listw(W); WL_Queen
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 108
## Percentage nonzero weights: 14.81481
## Average number of links: 4
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 27 729 27 16.57685 119.6581
plot(Jabar, axes=T, col="white")
text(CoordK[, 1], CoordK[, 2],
label = Jabar$NAME_2, col = "black", cex = 0.8, pos = 1.5)
points(CoordK[,1], CoordK[,2], pch=19,
cex=0.5,col="blue")
plot(W, coordinates(Jabar), add=TRUE, col = "blue", lwd=1)
#global moran Queen
MI_Queen = moran.test(Pneumonia$Pneumonia.pada.Balita, WL_Queen)
moran.plot(Pneumonia$Pneumonia.pada.Balita, WL_Queen)
#local moran Queen
locm<-localmoran(Pneumonia$Pneumonia.pada.Balita, WL_Queen)
head(locm)
## Ii E.Ii Var.Ii Z.Ii Pr(z != E(Ii))
## 1 0.16299652 -0.098978426 0.261429776 0.5123683 0.60839328
## 2 -0.08243156 -0.001540619 0.004509257 -1.2046126 0.22835288
## 3 -0.71238608 -0.024877497 0.654982400 -0.8494999 0.39560321
## 4 -0.18571886 -0.011003221 0.090104204 -0.5820486 0.56053393
## 5 0.07111622 -0.001981461 0.004805411 1.0544801 0.29166318
## 6 -0.65109112 -0.030177098 0.132752573 -1.7041585 0.08835147
#3. Berdasarkan Jarak
#K-nearesr neigbhors k = 3
W <- knn2nb(knearneigh(CoordK, k =
3), row.names = ID); W
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 81
## Percentage nonzero weights: 11.11111
## Average number of links: 3
## 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 1 0 0 0
## 2 0 0 0 0 0
## 3 0 0 0 0 0
## 4 0 0 0 0 0
## 5 0 0 0 0 0
WL_k3 <- nb2listw(W); WL_k3 #untuk Morans's I
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 81
## Percentage nonzero weights: 11.11111
## Average number of links: 3
## Non-symmetric neighbours list
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 27 729 27 15.66667 113.1111
plot(Jabar, axes=T, col="white")
text(CoordK[, 1], CoordK[, 2],
label = Jabar$NAME_2, col = "black", cex = 0.8, pos = 1.5)
points(CoordK[,1], CoordK[,2], pch=19,
cex=0.5,col="blue")
plot(W, coordinates(Jabar), add=TRUE, col = "blue", lwd=1)
#global moran K-nearesr neigbhors k = 3
MI_k3 <- moran.test(Pneumonia$Pneumonia.pada.Balita, WL_k3); MI_k3
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_k3
##
## Moran I statistic standard deviate = 0.48824, p-value = 0.3127
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02913092 -0.03846154 0.01916568
moran.plot(Pneumonia$Pneumonia.pada.Balita, WL_k3)
#local moran k = 3
locm<-localmoran(Pneumonia$Pneumonia.pada.Balita, WL_k3)
head(locm)
## Ii E.Ii Var.Ii Z.Ii Pr(z != E(Ii))
## 1 1.00525478 -0.098978426 0.73842445 1.28501373 0.19878749
## 2 -0.07248636 -0.001540619 0.01273667 -0.62863466 0.52958826
## 3 0.04810748 -0.024877497 0.20086127 0.16284910 0.87063725
## 4 0.01665444 -0.011003221 0.09010420 0.09213886 0.92658771
## 5 0.22476353 -0.001981461 0.01637399 1.77198658 0.07639679
## 6 -0.45327061 -0.030177098 0.24232613 -0.85948094 0.39007523
#K-nearesr neigbhors k = 5
W <- knn2nb(knearneigh(CoordK, k =
5), row.names = ID); W
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 135
## Percentage nonzero weights: 18.51852
## 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 1 0 0 0
## 2 1 0 0 0 0
## 3 0 0 0 0 0
## 4 0 0 0 0 1
## 5 0 0 0 1 0
WL_k5 <- nb2listw(W); WL_k5 #untuk Morans's I
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 27
## Number of nonzero links: 135
## Percentage nonzero weights: 18.51852
## Average number of links: 5
## Non-symmetric neighbours list
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 27 729 27 10.04 109.76
plot(Jabar, axes=T, col="white")
text(CoordK[, 1], CoordK[, 2],
label = Jabar$NAME_2, col = "black", cex = 0.8, pos = 1.5)
points(CoordK[,1], CoordK[,2], pch=19,
cex=0.5,col="blue")
plot(W, coordinates(Jabar), add=TRUE, col = "blue", lwd=1)
#global moran K-nearesr neigbhors k = 5
MI_k5 <- moran.test(Pneumonia$Pneumonia.pada.Balita, WL_k5); MI_k5
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_k5
##
## Moran I statistic standard deviate = 0.38552, p-value = 0.3499
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.002517668 -0.038461538 0.011298650
moran.plot(Pneumonia$Pneumonia.pada.Balita, WL_k5)
#local moran k = 5
locm<-localmoran(Pneumonia$Pneumonia.pada.Balita, WL_k5)
head(locm)
## Ii E.Ii Var.Ii Z.Ii Pr(z != E(Ii))
## 1 0.957407410 -0.098978426 0.404528180 1.66091795 0.09672993
## 2 -0.085524061 -0.001540619 0.006977481 -1.00541262 0.31469821
## 3 0.006655966 -0.024877497 0.110037043 0.09506096 0.92426640
## 4 -0.025003566 -0.011003221 0.049361433 -0.06301513 0.94975445
## 5 0.132282220 -0.001981461 0.008970100 1.41762028 0.15630166
## 6 -0.651091119 -0.030177098 0.132752573 -1.70415853 0.08835147
Dipilih nilai Indeks Morans yang paling kuat
MI_Rook
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_Rook
##
## Moran I statistic standard deviate = 0.47161, p-value = 0.3186
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02841137 -0.03846154 0.02010622
MI_Queen
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_Queen
##
## Moran I statistic standard deviate = 0.47161, p-value = 0.3186
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02841137 -0.03846154 0.02010622
MI_k3
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_k3
##
## Moran I statistic standard deviate = 0.48824, p-value = 0.3127
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.02913092 -0.03846154 0.01916568
MI_k5
##
## Moran I test under randomisation
##
## data: Pneumonia$Pneumonia.pada.Balita
## weights: WL_k5
##
## Moran I statistic standard deviate = 0.38552, p-value = 0.3499
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.002517668 -0.038461538 0.011298650
moran_index<-data.frame(
"Weight matrix"=c("Rook Contiguity", "Queen Contiguity", "KNN (k=3)", "KNN (k=5)"),
"Moran's I"=c(-0.07660499, -0.07660499, -0.10633730, -0.02836378),
"P-value"=c(MI_Rook$p.value, MI_Queen$p.value, MI_k3$p.value, MI_k5$p.value))
moran_index
## Weight.matrix Moran.s.I P.value
## 1 Rook Contiguity -0.07660499 0.3186018
## 2 Queen Contiguity -0.07660499 0.3186018
## 3 KNN (k=3) -0.10633730 0.3126888
## 4 KNN (k=5) -0.02836378 0.3499250
Pilih W dari Moran’s I yang signifikan dan terbesar. Dalam kasus ini adalah matriks knn k=5. Sehingga dapat disimpulkan bahwa Pneumonia memiliki autokorelasi spasial paling besar adalah matriks bobot KNN k=5.
library(spgwr)
adapt_gaussian <- gwr.sel(y ~ x1+x2, data = Jabar,
coords = CoordK, adapt = T)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = CoordK, adapt = T): data
## is Spatial* object, ignoring coords argument
## Adaptive q: 0.381966 CV score: 90036390
## Adaptive q: 0.618034 CV score: 84509261
## Adaptive q: 0.763932 CV score: 81872703
## Adaptive q: 0.854102 CV score: 82456138
## Adaptive q: 0.7778974 CV score: 81671777
## Adaptive q: 0.7971971 CV score: 81762187
## Adaptive q: 0.7834619 CV score: 81698226
## Adaptive q: 0.7725631 CV score: 81740894
## Adaptive q: 0.7792171 CV score: 81678079
## Adaptive q: 0.7758598 CV score: 81696028
## Adaptive q: 0.7780767 CV score: 81672635
## Adaptive q: 0.7771191 CV score: 81679627
## Adaptive q: 0.7776001 CV score: 81673467
## Adaptive q: 0.7778567 CV score: 81671583
## Adaptive q: 0.7777587 CV score: 81671449
## Adaptive q: 0.777718 CV score: 81671966
## Adaptive q: 0.7777994 CV score: 81671309
## Adaptive q: 0.7777994 CV score: 81671309
adapt_gaussian
## [1] 0.7777994
gwr.model1 = gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
adapt=adapt_gaussian,
hatmatrix=TRUE,
se.fit=TRUE) ; gwr.model1
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, adapt =
## adapt_gaussian, : data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, adapt = adapt_gaussian,
## hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.Gauss
## Adaptive quantile: 0.7777994 (about 21 of 27 data points)
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.4132e+03 1.4830e+03 1.5249e+03 1.5483e+03 1.5678e+03
## x1 -1.2603e-02 -8.7491e-03 1.1985e-02 1.8522e-02 2.2247e-02
## x2 5.6255e-03 6.2869e-03 8.5215e-03 1.3035e-02 1.3822e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 4.356949
## Effective degrees of freedom (residual: 2traceS - traceS'S): 22.64305
## Sigma (residual: 2traceS - traceS'S): 1440.936
## Effective number of parameters (model: traceS): 3.749144
## Effective degrees of freedom (model: traceS): 23.25086
## Sigma (model: traceS): 1421.978
## Sigma (ML): 1319.564
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 476.6837
## AIC (GWR p. 96, eq. 4.22): 468.3649
## Residual sum of squares: 47013708
## Quasi-global R2: 0.3423684
results1 <-as.data.frame(gwr.model1$SDF)
head(results1)
## sum.w X.Intercept. x1 x2 X.Intercept._se x1_se
## 1 20.48688 1439.216 0.01063526 0.008987312 483.0679 0.03441223
## 2 20.04687 1442.880 0.02207656 0.006156950 488.6686 0.03522475
## 3 20.86152 1541.792 -0.01248086 0.013821812 488.0247 0.03431295
## 4 21.44953 1564.995 0.01622163 0.006527719 486.6023 0.03472520
## 5 21.58849 1561.628 0.01835760 0.006188642 486.4319 0.03472294
## 6 20.75807 1534.936 -0.01203003 0.013772669 488.6962 0.03431517
## x2_se gwr.e pred pred.se localR2 X.Intercept._se_EDF
## 1 0.009002881 1261.5542 4083.446 509.2915 0.3694669 489.5084
## 2 0.009136122 -1018.0819 3427.082 1142.4565 0.3740187 495.1838
## 3 0.008974231 -2188.9964 3614.996 706.9299 0.3377665 494.5313
## 4 0.009022125 -1873.8832 5478.883 940.5491 0.3530078 493.0900
## 5 0.009023925 888.3295 2215.671 356.9170 0.3562021 492.9173
## 6 0.008977960 868.4579 3307.542 468.7027 0.3393769 495.2118
## x1_se_EDF x2_se_EDF pred.se.1
## 1 0.03487104 0.009122912 516.0817
## 2 0.03569439 0.009257930 1157.6884
## 3 0.03477043 0.009093881 716.3551
## 4 0.03518817 0.009142413 953.0890
## 5 0.03518589 0.009144237 361.6756
## 6 0.03477268 0.009097659 474.9517
adapt_bisquare <- gwr.sel(y ~ x1+x2, data = Jabar,
coords = coords, gweight = gwr.bisquare, adapt = T)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.bisquare, : data is Spatial* object, ignoring coords argument
## Adaptive q: 0.381966 CV score: 158693409
## Adaptive q: 0.618034 CV score: 148534921
## Adaptive q: 0.763932 CV score: 105142732
## Adaptive q: 0.854102 CV score: 93339097
## Adaptive q: 0.9018111 CV score: 92517907
## Adaptive q: 0.8883936 CV score: 92586356
## Adaptive q: 0.9023247 CV score: 92515422
## Adaptive q: 0.9396334 CV score: 91100239
## Adaptive q: 0.9626914 CV score: 89340850
## Adaptive q: 0.976942 CV score: 88718434
## Adaptive q: 0.9947889 CV score: 88031542
## Adaptive q: 0.987972 CV score: 88282865
## Adaptive q: 0.9921851 CV score: 88125958
## Adaptive q: 0.9967794 CV score: 87960664
## Adaptive q: 0.9980095 CV score: 87917416
## Adaptive q: 0.9987698 CV score: 87890898
## Adaptive q: 0.9992397 CV score: 87874590
## Adaptive q: 0.9995301 CV score: 87864541
## Adaptive q: 0.9997096 CV score: 87858342
## Adaptive q: 0.9998205 CV score: 87854516
## Adaptive q: 0.9998891 CV score: 87852152
## Adaptive q: 0.9999314 CV score: 87850692
## Adaptive q: 0.9999314 CV score: 87850692
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.bisquare, : Bandwidth converged to upper bound:1
adapt_bisquare
## [1] 0.9999314
gwr.model2 = gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
adapt=adapt_bisquare,
hatmatrix=TRUE,
se.fit=TRUE) ; gwr.model2
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, adapt =
## adapt_bisquare, : data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, adapt = adapt_bisquare,
## hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.Gauss
## Adaptive quantile: 0.9999314 (about 26 of 27 data points)
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.4577e+03 1.5212e+03 1.5468e+03 1.5600e+03 1.5722e+03
## x1 -9.4817e-03 -6.2722e-03 9.8115e-03 1.2813e-02 1.4338e-02
## x2 7.3242e-03 7.6152e-03 8.8865e-03 1.2302e-02 1.2949e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 3.945491
## Effective degrees of freedom (residual: 2traceS - traceS'S): 23.05451
## Sigma (residual: 2traceS - traceS'S): 1453.574
## Effective number of parameters (model: traceS): 3.505711
## Effective degrees of freedom (model: traceS): 23.49429
## Sigma (model: traceS): 1439.905
## Sigma (ML): 1343.177
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 476.8932
## AIC (GWR p. 96, eq. 4.22): 469.0792
## Residual sum of squares: 48711362
## Quasi-global R2: 0.3186214
results2 <-as.data.frame(gwr.model2$SDF)
head(results2)
## sum.w X.Intercept. x1 x2 X.Intercept._se x1_se
## 1 21.15189 1457.685 0.009326040 0.009158540 487.4370 0.03474310
## 2 22.71018 1506.380 0.012957891 0.007876703 486.2054 0.03483482
## 3 22.32573 1553.186 -0.009061118 0.012903278 488.3289 0.03456885
## 4 22.98640 1571.271 0.011353632 0.007716401 487.2282 0.03476033
## 5 22.94806 1567.875 0.013244486 0.007378749 487.5422 0.03478506
## 6 22.36949 1549.304 -0.008485803 0.012799824 488.2627 0.03455903
## x2_se gwr.e pred pred.se localR2 X.Intercept._se_EDF
## 1 0.009089769 1285.8439 4059.156 512.7769 0.3411872 492.0642
## 2 0.009084970 -1213.6332 3622.633 1135.8623 0.3369778 490.8208
## 3 0.009040909 -2244.7890 3670.789 708.0647 0.3155279 492.9645
## 4 0.009062990 -1861.1369 5466.137 942.4283 0.3268258 491.8534
## 5 0.009068448 914.6449 2189.355 357.9916 0.3292861 492.1703
## 6 0.009039847 896.0213 3279.979 471.3448 0.3167066 492.8977
## x1_se_EDF x2_se_EDF pred.se.1
## 1 0.03507291 0.009176056 517.6446
## 2 0.03516550 0.009171212 1146.6448
## 3 0.03489700 0.009126732 714.7862
## 4 0.03509030 0.009149023 951.3746
## 5 0.03511527 0.009154533 361.3899
## 6 0.03488709 0.009125660 475.8191
adapt_tricube <- gwr.sel(y ~ x1+x2, data = Jabar,
coords = coords, gweight = gwr.tricube, adapt = T)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.tricube, : data is Spatial* object, ignoring coords argument
## Adaptive q: 0.381966 CV score: 166515037
## Adaptive q: 0.618034 CV score: 145081681
## Adaptive q: 0.763932 CV score: 108968393
## Adaptive q: 0.854102 CV score: 92684803
## Adaptive q: 0.9098301 CV score: 91919398
## Adaptive q: 0.8879709 CV score: 92013096
## Adaptive q: 0.9065831 CV score: 91931604
## Adaptive q: 0.9442719 CV score: 89977636
## Adaptive q: 0.9655581 CV score: 88329347
## Adaptive q: 0.9787138 CV score: 87715084
## Adaptive q: 0.9868444 CV score: 87355731
## Adaptive q: 0.9918694 CV score: 87141717
## Adaptive q: 0.994975 CV score: 87012606
## Adaptive q: 0.9968944 CV score: 86934030
## Adaptive q: 0.9980806 CV score: 86885936
## Adaptive q: 0.9988138 CV score: 86856392
## Adaptive q: 0.9992669 CV score: 86838201
## Adaptive q: 0.9995469 CV score: 86826985
## Adaptive q: 0.99972 CV score: 86820063
## Adaptive q: 0.9998269 CV score: 86815788
## Adaptive q: 0.999893 CV score: 86813148
## Adaptive q: 0.9999339 CV score: 86811517
## Adaptive q: 0.9999339 CV score: 86811517
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.tricube, : Bandwidth converged to upper bound:1
adapt_tricube
## [1] 0.9999339
gwr.model3 = gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
adapt=adapt_tricube,
hatmatrix=TRUE,
se.fit=TRUE) ; gwr.model3
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, adapt =
## adapt_tricube, : data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, adapt = adapt_tricube,
## hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.Gauss
## Adaptive quantile: 0.9999339 (about 26 of 27 data points)
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.4577e+03 1.5212e+03 1.5468e+03 1.5600e+03 1.5722e+03
## x1 -9.4817e-03 -6.2722e-03 9.8114e-03 1.2813e-02 1.4338e-02
## x2 7.3242e-03 7.6152e-03 8.8865e-03 1.2302e-02 1.2949e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 3.945488
## Effective degrees of freedom (residual: 2traceS - traceS'S): 23.05451
## Sigma (residual: 2traceS - traceS'S): 1453.574
## Effective number of parameters (model: traceS): 3.505709
## Effective degrees of freedom (model: traceS): 23.49429
## Sigma (model: traceS): 1439.906
## Sigma (ML): 1343.177
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 476.8932
## AIC (GWR p. 96, eq. 4.22): 469.0792
## Residual sum of squares: 48711373
## Quasi-global R2: 0.3186213
results3 <-as.data.frame(gwr.model3$SDF)
head(results3)
## sum.w X.Intercept. x1 x2 X.Intercept._se x1_se
## 1 21.15191 1457.685 0.009326012 0.009158544 487.4370 0.03474310
## 2 22.71020 1506.381 0.012957805 0.007876720 486.2053 0.03483481
## 3 22.32574 1553.186 -0.009061099 0.012903273 488.3289 0.03456885
## 4 22.98641 1571.271 0.011353602 0.007716408 487.2282 0.03476033
## 5 22.94808 1567.875 0.013244432 0.007378762 487.5422 0.03478506
## 6 22.36951 1549.304 -0.008485781 0.012799818 488.2627 0.03455903
## x2_se gwr.e pred pred.se localR2 X.Intercept._se_EDF
## 1 0.009089769 1285.8444 4059.156 512.7769 0.3411870 492.0642
## 2 0.009084970 -1213.6351 3622.635 1135.8623 0.3369775 490.8208
## 3 0.009040910 -2244.7893 3670.789 708.0647 0.3155277 492.9645
## 4 0.009062990 -1861.1368 5466.137 942.4283 0.3268256 491.8534
## 5 0.009068448 914.6452 2189.355 357.9916 0.3292859 492.1703
## 6 0.009039848 896.0215 3279.979 471.3448 0.3167064 492.8976
## x1_se_EDF x2_se_EDF pred.se.1
## 1 0.03507291 0.009176056 517.6446
## 2 0.03516549 0.009171211 1146.6447
## 3 0.03489700 0.009126733 714.7862
## 4 0.03509030 0.009149023 951.3746
## 5 0.03511526 0.009154533 361.3899
## 6 0.03488709 0.009125661 475.8191
library(GWmodel)
adapt_exponential <- bw.gwr(y ~ x1+x2,
data=Jabar,
approach="CV",
kernel="exponential",
adaptive=T)
## Adaptive bandwidth: 24 CV score: 84014504
## Adaptive bandwidth: 23 CV score: 84134897
## Adaptive bandwidth: 26 CV score: 83590612
## Adaptive bandwidth: 26 CV score: 83590612
adapt_exponential
## [1] 26
gwr.model4 <- gwr.basic(y ~ x1+x2,
data=Jabar,
bw=adapt_exponential,
kernel="exponential") ; gwr.model4
## ***********************************************************************
## * Package GWmodel *
## ***********************************************************************
## Program starts at: 2024-11-07 01:01:04.46428
## Call:
## gwr.basic(formula = y ~ x1 + x2, data = Jabar, bw = adapt_exponential,
## kernel = "exponential")
##
## Dependent (y) variable: y
## Independent variables: x1 x2
## Number of data points: 27
## ***********************************************************************
## * Results of Global Regression *
## ***********************************************************************
##
## Call:
## lm(formula = formula, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2391.6 -1068.8 -222.7 1139.4 2647.8
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.590e+03 4.949e+02 3.213 0.00372 **
## x1 -1.670e-04 3.529e-02 -0.005 0.99626
## x2 1.049e-02 9.231e-03 1.136 0.26711
##
## ---Significance stars
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Residual standard error: 1479 on 24 degrees of freedom
## Multiple R-squared: 0.266
## Adjusted R-squared: 0.2048
## F-statistic: 4.349 on 2 and 24 DF, p-value: 0.02445
## ***Extra Diagnostic information
## Residual sum of squares: 52472584
## Sigma(hat): 1448.759
## AIC: 475.5817
## AICc: 477.3999
## BIC: 466.9484
## ***********************************************************************
## * Results of Geographically Weighted Regression *
## ***********************************************************************
##
## *********************Model calibration information*********************
## Kernel function: exponential
## Fixed bandwidth: 26
## Regression points: the same locations as observations are used.
## Distance metric: Euclidean distance metric is used.
##
## ****************Summary of GWR coefficient estimates:******************
## Min. 1st Qu. Median 3rd Qu. Max.
## Intercept 1.5825e+03 1.5840e+03 1.5860e+03 1.5888e+03 1591.6241
## x1 -1.8739e-03 -9.9438e-04 5.4461e-04 1.2875e-03 0.0015
## x2 1.0064e-02 1.0149e-02 1.0365e-02 1.0722e-02 0.0109
## ************************Diagnostic information*************************
## Number of data points: 27
## Effective number of parameters (2trace(S) - trace(S'S)): 3.170263
## Effective degrees of freedom (n-2trace(S) + trace(S'S)): 23.82974
## AICc (GWR book, Fotheringham, et al. 2002, p. 61, eq 2.33): 477.3419
## AIC (GWR book, Fotheringham, et al. 2002,GWR p. 96, eq. 4.22): 470.3598
## BIC (GWR book, Fotheringham, et al. 2002,GWR p. 61, eq. 2.34): 450.4441
## Residual sum of squares: 51878099
## R-square value: 0.2743249
## Adjusted R-square value: 0.1735536
##
## ***********************************************************************
## Program stops at: 2024-11-07 01:01:04.471529
results4 <-as.data.frame(gwr.model4$SDF)
head(results4)
## Intercept x1 x2 y yhat residual CV_Score
## 1 1582.520 0.0005446119 0.01037599 5345 3912.454 1432.5456 0
## 2 1583.136 0.0012034519 0.01019552 2409 3892.372 -1483.3715 0
## 3 1588.875 -0.0018738502 0.01090756 1426 3785.865 -2359.8645 0
## 4 1591.596 0.0012056519 0.01011219 3605 5416.638 -1811.6383 0
## 5 1590.462 0.0015490743 0.01007815 3104 2136.292 967.7077 0
## 6 1587.869 -0.0016715070 0.01087000 4176 3228.677 947.3231 0
## Stud_residual Intercept_SE x1_SE x2_SE Intercept_TV x1_TV
## 1 1.0403057 493.8682 0.03521382 0.009211814 3.204336 0.01546586
## 2 -1.6259412 493.8769 0.03521724 0.009212614 3.205527 0.03417223
## 3 -1.8540428 493.9714 0.03521914 0.009212743 3.216531 -0.05320545
## 4 -1.6343297 493.9632 0.03521725 0.009212363 3.222095 0.03423470
## 5 0.6786129 493.9602 0.03521882 0.009212751 3.219818 0.04398427
## 6 0.6841664 493.9615 0.03521828 0.009212705 3.214560 -0.04746134
## x2_TV Local_R2
## 1 1.126379 0.2758868
## 2 1.106691 0.2756193
## 3 1.183965 0.2728984
## 4 1.097676 0.2749558
## 5 1.093935 0.2755768
## 6 1.179892 0.2733910
bw5 <- gwr.sel(y ~ x1+x2, data = Jabar, coords = coords)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords): data is
## Spatial* object, ignoring coords argument
## Bandwidth: 95.40776 CV score: 88939359
## Bandwidth: 154.2189 CV score: 84158233
## Bandwidth: 190.5661 CV score: 83432138
## Bandwidth: 187.8928 CV score: 83466944
## Bandwidth: 213.03 CV score: 83206136
## Bandwidth: 226.9134 CV score: 83109381
## Bandwidth: 235.4938 CV score: 83060855
## Bandwidth: 240.7968 CV score: 83034334
## Bandwidth: 244.0743 CV score: 83019109
## Bandwidth: 246.0998 CV score: 83010111
## Bandwidth: 247.3517 CV score: 83004701
## Bandwidth: 248.1254 CV score: 83001412
## Bandwidth: 248.6036 CV score: 82999401
## Bandwidth: 248.8991 CV score: 82998166
## Bandwidth: 249.0817 CV score: 82997406
## Bandwidth: 249.1946 CV score: 82996937
## Bandwidth: 249.2644 CV score: 82996648
## Bandwidth: 249.3075 CV score: 82996469
## Bandwidth: 249.3341 CV score: 82996359
## Bandwidth: 249.3506 CV score: 82996291
## Bandwidth: 249.3608 CV score: 82996248
## Bandwidth: 249.3671 CV score: 82996222
## Bandwidth: 249.371 CV score: 82996206
## Bandwidth: 249.3734 CV score: 82996196
## Bandwidth: 249.3749 CV score: 82996190
## Bandwidth: 249.3758 CV score: 82996186
## Bandwidth: 249.3763 CV score: 82996184
## Bandwidth: 249.3767 CV score: 82996183
## Bandwidth: 249.3769 CV score: 82996182
## Bandwidth: 249.377 CV score: 82996181
## Bandwidth: 249.3771 CV score: 82996181
## Bandwidth: 249.3772 CV score: 82996181
## Bandwidth: 249.3772 CV score: 82996181
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords): Bandwidth
## converged to upper bound:249.377265284254
bw5
## [1] 249.3772
gwr.model5 = gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
bandwidth = bw5,
hatmatrix=TRUE,
se.fit=TRUE) ; gwr.model5
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw5, :
## data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw5,
## hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.Gauss
## Fixed bandwidth: 249.3772
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.5561e+03 1.5619e+03 1.5655e+03 1.5710e+03 1.5766e+03
## x1 -7.2214e-03 -3.4590e-03 2.0447e-03 6.8521e-03 1.0248e-02
## x2 8.0659e-03 8.8629e-03 1.0035e-02 1.1472e-02 1.2407e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 3.499937
## Effective degrees of freedom (residual: 2traceS - traceS'S): 23.50006
## Sigma (residual: 2traceS - traceS'S): 1467.354
## Effective number of parameters (model: traceS): 3.260544
## Effective degrees of freedom (model: traceS): 23.73946
## Sigma (model: traceS): 1459.937
## Sigma (ML): 1368.95
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 477.1829
## AIC (GWR p. 96, eq. 4.22): 469.8604
## Residual sum of squares: 50598653
## Quasi-global R2: 0.2922219
results5 <-as.data.frame(gwr.model5$SDF)
head(results5)
## sum.w X.Intercept. x1 x2 X.Intercept._se x1_se
## 1 25.51818 1560.827 0.001949978 0.010177940 488.8990 0.03486209
## 2 25.55345 1564.191 0.004050783 0.009632970 489.0701 0.03489155
## 3 23.18272 1559.977 -0.007221423 0.012406739 492.7142 0.03497462
## 4 24.15999 1576.573 0.007805744 0.008575165 491.1945 0.03503445
## 5 23.78761 1572.197 0.010247596 0.008074540 492.0794 0.03510030
## 6 23.40545 1558.593 -0.006409531 0.012226360 492.2353 0.03495539
## x2_se gwr.e pred pred.se localR2 X.Intercept._se_EDF
## 1 0.009119991 1411.4932 3933.507 510.9066 0.2968131 491.3829
## 2 0.009124675 -1417.4892 3826.489 1140.5664 0.2977337 491.5548
## 3 0.009147465 -2274.9530 3700.953 714.6311 0.2894648 495.2175
## 4 0.009150576 -1849.7955 5454.795 950.2156 0.2980293 493.6901
## 5 0.009163239 929.4986 2174.501 361.4721 0.3004378 494.5794
## 6 0.009143413 911.8502 3264.150 476.5236 0.2904420 494.7361
## x1_se_EDF x2_se_EDF pred.se.1
## 1 0.03503921 0.009166326 513.5023
## 2 0.03506882 0.009171033 1146.3611
## 3 0.03515231 0.009193939 718.2618
## 4 0.03521245 0.009197066 955.0432
## 5 0.03527863 0.009209793 363.3086
## 6 0.03513298 0.009189867 478.9446
bw6 <- gwr.sel(y ~ x1+x2, data = Jabar,
coords = coords, gweight = gwr.bisquare)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.bisquare): data is Spatial* object, ignoring coords argument
## Bandwidth: 95.40776 CV score: 168283364
## Bandwidth: 154.2189 CV score: 96723830
## Bandwidth: 190.5661 CV score: 92978823
## Bandwidth: 176.7942 CV score: 93455965
## Bandwidth: 189.3986 CV score: 92990896
## Bandwidth: 192.6633 CV score: 92966339
## Bandwidth: 214.3261 CV score: 92694248
## Bandwidth: 227.7145 CV score: 91252705
## Bandwidth: 235.9889 CV score: 90279510
## Bandwidth: 241.1028 CV score: 89706894
## Bandwidth: 244.2634 CV score: 89369439
## Bandwidth: 246.2167 CV score: 89167829
## Bandwidth: 247.4239 CV score: 89045975
## Bandwidth: 248.17 CV score: 88971730
## Bandwidth: 248.6312 CV score: 88926254
## Bandwidth: 248.9161 CV score: 88898304
## Bandwidth: 249.0923 CV score: 88881091
## Bandwidth: 249.2011 CV score: 88870475
## Bandwidth: 249.2684 CV score: 88863923
## Bandwidth: 249.31 CV score: 88859877
## Bandwidth: 249.3357 CV score: 88857377
## Bandwidth: 249.3516 CV score: 88855833
## Bandwidth: 249.3614 CV score: 88854879
## Bandwidth: 249.3674 CV score: 88854289
## Bandwidth: 249.3712 CV score: 88853925
## Bandwidth: 249.3735 CV score: 88853700
## Bandwidth: 249.3749 CV score: 88853561
## Bandwidth: 249.3758 CV score: 88853475
## Bandwidth: 249.3764 CV score: 88853422
## Bandwidth: 249.3767 CV score: 88853389
## Bandwidth: 249.3769 CV score: 88853368
## Bandwidth: 249.3771 CV score: 88853356
## Bandwidth: 249.3771 CV score: 88853348
## Bandwidth: 249.3772 CV score: 88853343
## Bandwidth: 249.3772 CV score: 88853343
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.bisquare): Bandwidth converged to upper bound:249.377265284254
bw6
## [1] 249.3772
gwr.model6<-gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
bandwidth = bw6,
gweight = gwr.bisquare,
se.fit=T, hatmatrix=T) ; gwr.model6
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw6, :
## data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw6,
## gweight = gwr.bisquare, hatmatrix = T, se.fit = T)
## Kernel function: gwr.bisquare
## Fixed bandwidth: 249.3772
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.4210e+03 1.4566e+03 1.4791e+03 1.5013e+03 1.5365e+03
## x1 -3.4154e-02 -1.4704e-02 9.3211e-03 3.1951e-02 4.9516e-02
## x2 -7.8729e-04 3.1745e-03 8.5980e-03 1.4963e-02 1.9931e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 5.135199
## Effective degrees of freedom (residual: 2traceS - traceS'S): 21.8648
## Sigma (residual: 2traceS - traceS'S): 1428.776
## Effective number of parameters (model: traceS): 4.29701
## Effective degrees of freedom (model: traceS): 22.70299
## Sigma (model: traceS): 1402.153
## Sigma (ML): 1285.746
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 477.03
## AIC (GWR p. 96, eq. 4.22): 467.5108
## Residual sum of squares: 44634821
## Quasi-global R2: 0.3756444
results6 <-as.data.frame(gwr.model6$SDF)
results6
## sum.w X.Intercept. x1 x2 X.Intercept._se x1_se
## 1 21.33041 1454.418 0.0091436603 0.0092138356 474.6135 0.03382529
## 2 21.50522 1471.695 0.0179798927 0.0069110670 477.3841 0.03432761
## 3 14.31363 1449.782 -0.0341535656 0.0199310819 578.1429 0.03593563
## 4 17.04557 1530.718 0.0364959676 0.0017009824 517.5919 0.03731685
## 5 16.14198 1522.013 0.0484813808 -0.0007297033 533.5566 0.03853612
## 6 14.86810 1439.622 -0.0296073278 0.0188697352 565.3089 0.03551080
## 7 19.28499 1446.111 0.0274061889 0.0050538374 487.5175 0.03491664
## 8 21.85593 1472.930 0.0130723675 0.0080448023 474.6162 0.03404919
## 9 15.95426 1490.892 -0.0278915619 0.0177654640 525.6531 0.03483443
## 10 15.89979 1536.494 0.0476900075 -0.0007872878 535.9572 0.03876108
## 11 19.49235 1429.242 0.0005842356 0.0115475905 484.9344 0.03375319
## 12 18.41874 1502.300 -0.0083542593 0.0126997872 487.6914 0.03390849
## 13 18.55175 1518.761 0.0262264772 0.0042250902 499.4795 0.03585963
## 14 21.90012 1470.623 0.0093211220 0.0089605360 473.6755 0.03387726
## 15 16.54502 1534.516 0.0419440718 0.0004588219 526.2779 0.03803217
## 16 16.56555 1515.530 0.0468066597 -0.0002740121 528.7059 0.03817810
## 17 15.92164 1491.170 -0.0279904908 0.0177867620 525.8492 0.03484117
## 18 17.98587 1470.093 0.0388747448 0.0021238858 506.6036 0.03639891
## 19 17.48628 1437.258 -0.0170768285 0.0157364662 518.2718 0.03438866
## 20 15.91036 1474.656 -0.0299722997 0.0185271540 540.1060 0.03514676
## 21 19.08291 1481.544 -0.0143088298 0.0144911156 494.7173 0.03400550
## 22 21.13265 1495.025 0.0195704991 0.0061947282 482.1764 0.03469492
## 23 20.76940 1500.264 0.0089313160 0.0085979570 478.6559 0.03413030
## 24 15.24962 1458.743 0.0495161415 -0.0001014527 530.0915 0.03775726
## 25 21.06225 1479.097 -0.0030927503 0.0117996931 478.8931 0.03372431
## 26 16.82733 1421.005 -0.0150985364 0.0154351234 521.9320 0.03437880
## 27 20.97399 1485.302 0.0237110718 0.0053833324 484.2423 0.03490825
## x2_se gwr.e pred pred.se localR2 X.Intercept._se_EDF
## 1 0.008849605 1288.1088 4056.891 499.2223 0.3918788 483.6252
## 2 0.008923388 -1102.0921 3511.092 1115.7393 0.3918011 486.4483
## 3 0.009626890 -1877.1020 3303.102 807.3070 0.3604172 589.1202
## 4 0.009436168 -1964.3220 5569.322 989.4484 0.3821471 527.4195
## 5 0.009667038 731.3222 2372.678 390.7342 0.3887056 543.6874
## 6 0.009483524 713.4149 3462.585 506.4762 0.3672118 576.0426
## 7 0.009031996 2400.1127 3336.887 344.9554 0.3999919 496.7741
## 8 0.008876635 1855.2242 3441.776 345.5276 0.3894007 483.6279
## 9 0.009148027 1932.5731 2719.427 339.1492 0.3519106 535.6338
## 10 0.009706888 1520.5842 3430.416 407.5592 0.3852700 546.1335
## 11 0.008874372 -659.1641 2188.164 317.9807 0.3926344 494.1420
## 12 0.008867321 314.7576 2275.242 312.6039 0.3650927 496.9514
## 13 0.009177723 -1563.9039 1763.904 434.7796 0.3811860 508.9633
## 14 0.008850167 1361.8141 2195.186 322.8778 0.3881431 482.6693
## 15 0.009568271 372.3766 2832.623 375.0250 0.3835794 536.2705
## 16 0.009600665 1584.8552 3297.145 405.7106 0.3897361 538.7446
## 17 0.009149598 749.7018 2117.298 346.0668 0.3517024 535.8336
## 18 0.009284687 -1053.9843 2825.984 292.5302 0.3975594 516.2226
## 19 0.009085073 949.7591 2288.241 466.4545 0.3783240 528.1123
## 20 0.009279668 -1086.7691 1584.769 507.1502 0.3558894 550.3611
## 21 0.008918146 -1405.9234 2894.923 389.5229 0.3681914 504.1106
## 22 0.008979920 -387.3369 2220.337 331.9667 0.3862807 491.3316
## 23 0.008884582 -881.7051 1828.705 391.6137 0.3788333 487.7443
## 24 0.009530571 -161.4695 1704.469 474.0418 0.4001178 540.1565
## 25 0.008835338 -1701.4345 2525.434 519.5928 0.3780640 487.9859
## 26 0.009101125 -597.1877 1639.188 457.0686 0.3837248 531.8421
## 27 0.009019082 -814.6744 1932.674 387.5975 0.3906047 493.4367
## x1_se_EDF x2_se_EDF pred.se.1
## 1 0.03446754 0.009017635 508.7011
## 2 0.03497939 0.009092819 1136.9242
## 3 0.03661795 0.009809678 822.6355
## 4 0.03802539 0.009615335 1008.2354
## 5 0.03926782 0.009850589 398.1532
## 6 0.03618505 0.009663590 516.0928
## 7 0.03557961 0.009203489 351.5051
## 8 0.03469569 0.009045178 352.0883
## 9 0.03549585 0.009321723 345.5888
## 10 0.03949705 0.009891195 415.2977
## 11 0.03439407 0.009042872 324.0183
## 12 0.03455232 0.009035687 318.5394
## 13 0.03654050 0.009351983 443.0349
## 14 0.03452050 0.009018208 329.0084
## 15 0.03875429 0.009749947 382.1457
## 16 0.03890300 0.009782955 413.4139
## 17 0.03550271 0.009323324 352.6376
## 18 0.03709003 0.009460978 298.0846
## 19 0.03504161 0.009257574 475.3112
## 20 0.03581410 0.009455863 516.7796
## 21 0.03465117 0.009087477 396.9189
## 22 0.03535368 0.009150424 338.2698
## 23 0.03477834 0.009053276 399.0493
## 24 0.03847417 0.009711531 483.0426
## 25 0.03436464 0.009003097 529.4584
## 26 0.03503156 0.009273931 465.7470
## 27 0.03557106 0.009190330 394.9569
bw7 <- gwr.sel(y ~ x1+x2, data = Jabar,
coords = coords, gweight = gwr.tricube)
## Warning in gwr.sel(y ~ x1 + x2, data = Jabar, coords = coords, gweight =
## gwr.tricube): data is Spatial* object, ignoring coords argument
## Bandwidth: 95.40776 CV score: 185396875
## Bandwidth: 154.2189 CV score: 95297119
## Bandwidth: 190.5661 CV score: 90436663
## Bandwidth: 176.9426 CV score: 90907544
## Bandwidth: 187.7148 CV score: 90455920
## Bandwidth: 190.4485 CV score: 90436757
## Bandwidth: 190.691 CV score: 90436626
## Bandwidth: 190.6999 CV score: 90436625
## Bandwidth: 190.7001 CV score: 90436625
## Bandwidth: 190.7002 CV score: 90436625
## Bandwidth: 190.7001 CV score: 90436625
## Bandwidth: 190.7001 CV score: 90436625
bw7
## [1] 190.7001
gwr.model7<-gwr(y ~ x1+x2,
data = Jabar,
coords=coords,
bandwidth = bw7,
gweight = gwr.tricube,
se.fit=T, hatmatrix=T) ; gwr.model7
## Warning in gwr(y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw7, :
## data is Spatial* object, ignoring coords argument
## Call:
## gwr(formula = y ~ x1 + x2, data = Jabar, coords = coords, bandwidth = bw7,
## gweight = gwr.tricube, hatmatrix = T, se.fit = T)
## Kernel function: gwr.tricube
## Fixed bandwidth: 190.7001
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max.
## X.Intercept. 1.2785e+03 1.3913e+03 1.4413e+03 1.4970e+03 1.5810e+03
## x1 -6.3080e-02 -2.2999e-02 1.6294e-02 4.9696e-02 7.5654e-02
## x2 -5.8397e-03 -5.9949e-04 7.0038e-03 1.7728e-02 2.7627e-02
## Global
## X.Intercept. 1590.1581
## x1 -0.0002
## x2 0.0105
## Number of data points: 27
## Effective number of parameters (residual: 2traceS - traceS'S): 5.844947
## Effective degrees of freedom (residual: 2traceS - traceS'S): 21.15505
## Sigma (residual: 2traceS - traceS'S): 1385.869
## Effective number of parameters (model: traceS): 4.946715
## Effective degrees of freedom (model: traceS): 22.05329
## Sigma (model: traceS): 1357.352
## Sigma (ML): 1226.726
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 476.6897
## AIC (GWR p. 96, eq. 4.22): 465.623
## Residual sum of squares: 40631098
## Quasi-global R2: 0.4316488
results7 <-as.data.frame(gwr.model7$SDF)
bw8 <- bw.gwr(y ~ x1+x2,
data=Jabar,
approach="CV",
kernel="exponential",
adaptive=F)
## Fixed bandwidth: 1.23965 CV score: 86035383
## Fixed bandwidth: 0.766299 CV score: 90521019
## Fixed bandwidth: 1.532197 CV score: 84986673
## Fixed bandwidth: 1.713001 CV score: 84575952
## Fixed bandwidth: 1.824744 CV score: 84377441
## Fixed bandwidth: 1.893805 CV score: 84270876
## Fixed bandwidth: 1.936487 CV score: 84210279
## Fixed bandwidth: 1.962866 CV score: 84174661
## Fixed bandwidth: 1.979169 CV score: 84153310
## Fixed bandwidth: 1.989245 CV score: 84140358
## Fixed bandwidth: 1.995472 CV score: 84132444
## Fixed bandwidth: 1.999321 CV score: 84127588
## Fixed bandwidth: 2.001699 CV score: 84124600
## Fixed bandwidth: 2.003169 CV score: 84122758
## Fixed bandwidth: 2.004078 CV score: 84121621
## Fixed bandwidth: 2.004639 CV score: 84120920
## Fixed bandwidth: 2.004986 CV score: 84120486
## Fixed bandwidth: 2.005201 CV score: 84120219
## Fixed bandwidth: 2.005333 CV score: 84120053
## Fixed bandwidth: 2.005415 CV score: 84119951
## Fixed bandwidth: 2.005466 CV score: 84119888
## Fixed bandwidth: 2.005497 CV score: 84119849
bw8
## [1] 2.005497
gwr.model8 <- gwr.basic(y ~ x1+x2,
data=Jabar,
bw=bw8,
kernel="exponential") ; gwr.model8
## ***********************************************************************
## * Package GWmodel *
## ***********************************************************************
## Program starts at: 2024-11-07 01:01:05.070294
## Call:
## gwr.basic(formula = y ~ x1 + x2, data = Jabar, bw = bw8, kernel = "exponential")
##
## Dependent (y) variable: y
## Independent variables: x1 x2
## Number of data points: 27
## ***********************************************************************
## * Results of Global Regression *
## ***********************************************************************
##
## Call:
## lm(formula = formula, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2391.6 -1068.8 -222.7 1139.4 2647.8
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.590e+03 4.949e+02 3.213 0.00372 **
## x1 -1.670e-04 3.529e-02 -0.005 0.99626
## x2 1.049e-02 9.231e-03 1.136 0.26711
##
## ---Significance stars
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Residual standard error: 1479 on 24 degrees of freedom
## Multiple R-squared: 0.266
## Adjusted R-squared: 0.2048
## F-statistic: 4.349 on 2 and 24 DF, p-value: 0.02445
## ***Extra Diagnostic information
## Residual sum of squares: 52472584
## Sigma(hat): 1448.759
## AIC: 475.5817
## AICc: 477.3999
## BIC: 466.9484
## ***********************************************************************
## * Results of Geographically Weighted Regression *
## ***********************************************************************
##
## *********************Model calibration information*********************
## Kernel function: exponential
## Fixed bandwidth: 2.005497
## Regression points: the same locations as observations are used.
## Distance metric: Euclidean distance metric is used.
##
## ****************Summary of GWR coefficient estimates:******************
## Min. 1st Qu. Median 3rd Qu. Max.
## Intercept 1.4883e+03 1.5099e+03 1.5391e+03 1.5888e+03 1624.6939
## x1 -2.2692e-02 -1.0260e-02 9.1261e-03 1.7909e-02 0.0209
## x2 5.2296e-03 6.4350e-03 9.0109e-03 1.3421e-02 0.0159
## ************************Diagnostic information*************************
## Number of data points: 27
## Effective number of parameters (2trace(S) - trace(S'S)): 5.270634
## Effective degrees of freedom (n-2trace(S) + trace(S'S)): 21.72937
## AICc (GWR book, Fotheringham, et al. 2002, p. 61, eq 2.33): 476.9794
## AIC (GWR book, Fotheringham, et al. 2002,GWR p. 96, eq. 4.22): 467.6144
## BIC (GWR book, Fotheringham, et al. 2002,GWR p. 61, eq. 2.34): 450.3233
## Residual sum of squares: 44919669
## R-square value: 0.37166
## Adjusted R-square value: 0.2118987
##
## ***********************************************************************
## Program stops at: 2024-11-07 01:01:05.075995
results8 <-as.data.frame(gwr.model8$SDF)
head(results8)
## Intercept x1 x2 y yhat residual CV_Score
## 1 1489.517 0.009126127 0.009095625 5345 4064.748 1280.2518 0
## 2 1499.529 0.016668727 0.007006212 2409 3523.951 -1114.9507 0
## 3 1589.825 -0.022691931 0.015925497 1426 3394.034 -1968.0339 0
## 4 1620.130 0.016188821 0.005850927 3605 5282.494 -1677.4945 0
## 5 1610.748 0.020607390 0.005522323 3104 2271.793 832.2068 0
## 6 1574.003 -0.019836053 0.015413632 4176 3360.848 815.1519 0
## Stud_residual Intercept_SE x1_SE x2_SE Intercept_TV x1_TV
## 1 1.0092292 485.7765 0.03464532 0.009058239 3.066259 0.2634159
## 2 -1.4364465 486.6285 0.03509008 0.009162545 3.081465 0.4750268
## 3 -1.8834691 503.0382 0.03555728 0.009248407 3.160446 -0.6381795
## 4 -1.8579963 498.3776 0.03514325 0.009131485 3.250809 0.4606524
## 5 0.6217145 498.1327 0.03537477 0.009195799 3.233572 0.5825448
## 6 0.6681317 501.3655 0.03543335 0.009251699 3.139432 -0.5598131
## x2_TV Local_R2
## 1 1.0041273 0.3888639
## 2 0.7646579 0.3837277
## 3 1.7219718 0.3560830
## 4 0.6407421 0.3717502
## 5 0.6005267 0.3789652
## 6 1.6660326 0.3630271
Perbandingan_AIC_GWR = data.frame("Model"=c("Model 1", "Model 2", "Model 3", "Model 4", "Model 5", "Model 6", "Model 7", "Model 8"),
"AIC"=c(gwr.model1$results$AICh,gwr.model2$results$AICh,gwr.model3$results$AICh,gwr.model4$GW.diagnostic$AIC,gwr.model5$results$AICh,gwr.model6$results$AICh,gwr.model7$results$AICh,gwr.model8$GW.diagnostic$AIC))
Perbandingan_AIC_GWR
## Model AIC
## 1 Model 1 468.3649
## 2 Model 2 469.0792
## 3 Model 3 469.0792
## 4 Model 4 470.3598
## 5 Model 5 469.8604
## 6 Model 6 467.5108
## 7 Model 7 465.6230
## 8 Model 8 467.6144
resid_gwr7 <- results7$gwr.e
Jabar$resid_gwr7 <- resid_gwr7
spplot(Jabar, zcol="resid_gwr7", main="Peta Sebaran Residual Model Tri-cube")
#Residuals diagnostics
shapiro.test(resid_gwr7) #Normality
##
## Shapiro-Wilk normality test
##
## data: resid_gwr7
## W = 0.94936, p-value = 0.2068
library(lmtest)
bptest(gwr.model7$lm, weights = gwr.model7$gweight) #Heteroskedasticity
##
## studentized Breusch-Pagan test
##
## data: gwr.model7$lm
## BP = 9.3666, df = 2, p-value = 0.009248
gwr.morantest(gwr.model7, WL_k5,zero.policy = TRUE) #Non-autocorrelation
##
## Leung et al. 2000 three moment approximation for Moran's I
##
## data: GWR residuals
## statistic = 15.71, df = 17.417, p-value = 0.5734
## sample estimates:
## I
## -0.1034682
#Comparing OLS and GWR models under an inferential framework
BFC02.gwr.test(gwr.model7)
##
## Brunsdon, Fotheringham & Charlton (2002, pp. 91-2) ANOVA
##
## data: gwr.model7
## F = 1.2914, df1 = 24.000, df2 = 21.155, p-value = 0.2781
## alternative hypothesis: greater
## sample estimates:
## SS OLS residuals SS GWR residuals
## 52472584 40631098
#Parameter Significance Test
LMZ.F3GWR.test(gwr.model7)
##
## Leung et al. (2000) F(3) test
##
## F statistic Numerator d.f. Denominator d.f. Pr(>)
## (Intercept) 6.3352e+10 9.0959e+00 22.088 < 2.2e-16 ***
## x1 5.2749e+12 5.9963e+00 22.088 < 2.2e-16 ***
## x2 6.2586e+12 6.7852e+00 22.088 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
data.frame("MODEL" = c("GWR","OLS"),
"AIC" = c(gwr.model7$results$AICh,AIC(jabarlm)))%>% arrange(AIC)
## MODEL AIC
## 1 GWR 465.6230
## 2 OLS 475.5817
Koefisien pada tiap lokasi
Pneumonia_koefisien <-data.frame(
Nama_Daerah = Pneumonia$Kabupaten.Kota,
Koefisien_X1 = results7$x1, # Koefisien untuk variabel X1
Koefisien_x2 = results7$x2 # Koefisien untuk variabel x2
)
Pneumonia_koefisien
## Nama_Daerah Koefisien_X1 Koefisien_x2
## 1 KABUPATEN BANDUNG 0.0133903571 0.0088996016
## 2 KABUPATEN BANDUNG BARAT 0.0283851191 0.0049467099
## 3 KOTA BANJAR -0.0630797491 0.0276269682
## 4 KABUPATEN BEKASI 0.0518548589 -0.0020361813
## 5 KABUPATEN BOGOR 0.0673966316 -0.0054863577
## 6 KABUPATEN CIAMIS -0.0526607274 0.0250142124
## 7 KABUPATEN CIANJUR 0.0475362402 0.0013405903
## 8 KOTA CIMAHI 0.0200785659 0.0068245099
## 9 KABUPATEN CIREBON -0.0521736745 0.0247104302
## 10 KOTA DEPOK 0.0652957794 -0.0052651472
## 11 KABUPATEN GARUT 0.0007943289 0.0126669732
## 12 KABUPATEN INDRAMAYU -0.0135804272 0.0142248980
## 13 KABUPATEN KARAWANG 0.0405483891 0.0008371993
## 14 KOTA BANDUNG 0.0140089691 0.0083002672
## 15 KOTA BEKASI 0.0584359003 -0.0036009305
## 16 KOTA BOGOR 0.0659645510 -0.0050060729
## 17 KOTA CIREBON -0.0524965389 0.0247914484
## 18 KOTA SUKABUMI 0.0615954081 -0.0028292154
## 19 KOTA TASIKMALAYA -0.0268584803 0.0189965033
## 20 KABUPATEN KUNINGAN -0.0537624429 0.0251855156
## 21 KABUPATEN MAJALENGKA -0.0217496938 0.0170090010
## 22 KABUPATEN PURWAKARTA 0.0322606479 0.0034949815
## 23 KABUPATEN SUBANG 0.0162937588 0.0070038484
## 24 KABUPATEN SUKABUMI 0.0756541182 -0.0058396463
## 25 KABUPATEN SUMEDANG -0.0034372349 0.0123562084
## 26 KABUPATEN TASIKMALAYA -0.0242487000 0.0184475806
## 27 KABUPATEN PANGANDARAN 0.0384760762 0.0022989508
#Map the results
gwr.map <- cbind(Jabar, as.matrix(results7))
gwr.map2 <- st_as_sf(gwr.map)
qtm(gwr.map2, fill = "localR2") +
tm_text("NAME_2", size = 0.8, col = "black")
class(gwr.map2)
## [1] "sf" "data.frame"
names(gwr.map2)
## [1] "GID_0" "NAME_0" "GID_1"
## [4] "NAME_1" "NL_NAME_1" "GID_2"
## [7] "NAME_2" "VARNAME_2" "NL_NAME_2"
## [10] "TYPE_2" "ENGTYPE_2" "CC_2"
## [13] "HASC_2" "Pneumonia.pada.Balita" "Imunisasi"
## [16] "Gizi.Buruk" "id" "resid_gwr7"
## [19] "sum.w" "X.Intercept." "x1"
## [22] "x2" "X.Intercept._se" "x1_se"
## [25] "x2_se" "gwr.e" "pred"
## [28] "pred.se" "localR2" "X.Intercept._se_EDF"
## [31] "x1_se_EDF" "x2_se_EDF" "pred.se.1"
## [34] "geometry"
#Attach coefficients to original dataframe
Jabar$coefx1<-results7$x1
Jabar$coefx2<-results7$x2
#Spatial distribution of x1
map1 <- tm_shape(gwr.map2) +
tm_fill("x1",
n = 5,
style = "quantile",
title = "x1 Coefficient") +
tm_borders(col = "black", lwd = 0.5) +
tm_text("NAME_2",
size = 0.6,
col = "black",
shadow = TRUE,
remove.overlap = TRUE) +
tm_layout(frame = FALSE,
legend.text.size = 0.5,
legend.title.size = 0.6)
map1
## Variable(s) "x1" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.
#Spatial distribution of x2
map2 <- tm_shape(gwr.map2) +
tm_fill("x2",
n = 5,
style = "quantile",
title = "x2 Coefficient") +
tm_borders(col = "black", lwd = 0.5)+
tm_text("NAME_2",
size = 0.6,
col = "black",
shadow = TRUE,
remove.overlap = TRUE) +
tm_layout(frame = FALSE,
legend.text.size = 0.5,
legend.title.size = 0.6)
map2
## Variable(s) "x2" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.
#Local significance test
t1 = gwr.model7$SDF$x1 / gwr.model7$SDF$x1_se
t2 = gwr.model7$SDF$x2 / gwr.model7$SDF$x2_se
pvalue1 = 2*pt(q=abs(t1), df=26, lower.tail=F)
pvalue2 = 2*pt(q=abs(t2), df=26, lower.tail=F)
library(dplyr)
data.frame(Kabupaten.Kota = Pneumonia$Kabupaten.Kota, t1, pvalue1) %>%
mutate(ket=ifelse(pvalue1 < 0.05, "Reject", "Not Reject"))
## Kabupaten.Kota t1 pvalue1 ket
## 1 KABUPATEN BANDUNG 0.40384488 0.68962927 Not Reject
## 2 KABUPATEN BANDUNG BARAT 0.82259260 0.41822044 Not Reject
## 3 KOTA BANJAR -1.60240861 0.12114641 Not Reject
## 4 KABUPATEN BEKASI 1.31365071 0.20043811 Not Reject
## 5 KABUPATEN BOGOR 1.64694806 0.11160402 Not Reject
## 6 KABUPATEN CIAMIS -1.40614037 0.17152191 Not Reject
## 7 KABUPATEN CIANJUR 1.30409577 0.20362873 Not Reject
## 8 KOTA CIMAHI 0.59560335 0.55658981 Not Reject
## 9 KABUPATEN CIREBON -1.44170777 0.16132237 Not Reject
## 10 KOTA DEPOK 1.59281364 0.12328824 Not Reject
## 11 KABUPATEN GARUT 0.02396477 0.98106369 Not Reject
## 12 KABUPATEN INDRAMAYU -0.40283182 0.69036519 Not Reject
## 13 KABUPATEN KARAWANG 1.06426604 0.29699288 Not Reject
## 14 KOTA BANDUNG 0.42124458 0.67703859 Not Reject
## 15 KOTA BEKASI 1.45449591 0.15777582 Not Reject
## 16 KOTA BOGOR 1.62499523 0.11622599 Not Reject
## 17 KOTA CIREBON -1.44917045 0.15924507 Not Reject
## 18 KOTA SUKABUMI 1.57669812 0.12695595 Not Reject
## 19 KOTA TASIKMALAYA -0.78525422 0.43940193 Not Reject
## 20 KABUPATEN KUNINGAN -1.46402533 0.15517387 Not Reject
## 21 KABUPATEN MAJALENGKA -0.64647263 0.52363949 Not Reject
## 22 KABUPATEN PURWAKARTA 0.90356359 0.37452283 Not Reject
## 23 KABUPATEN SUBANG 0.47463941 0.63900653 Not Reject
## 24 KABUPATEN SUKABUMI 1.84518191 0.07643089 Not Reject
## 25 KABUPATEN SUMEDANG -0.10423692 0.91778127 Not Reject
## 26 KABUPATEN TASIKMALAYA -0.70875150 0.48478562 Not Reject
## 27 KABUPATEN PANGANDARAN 1.06440208 0.29693243 Not Reject
data.frame(Kabupaten.Kota = Pneumonia$Kabupaten.Kota, t2, pvalue2) %>%
mutate(ket=ifelse(pvalue2<0.05,"Reject", "Not Reject"))
## Kabupaten.Kota t2 pvalue2 ket
## 1 KABUPATEN BANDUNG 1.02629115 0.31420842 Not Reject
## 2 KABUPATEN BANDUNG BARAT 0.55895967 0.58097228 Not Reject
## 3 KOTA BANJAR 2.47900272 0.01998319 Reject
## 4 KABUPATEN BEKASI -0.20860382 0.83638316 Not Reject
## 5 KABUPATEN BOGOR -0.54397161 0.59109572 Not Reject
## 6 KABUPATEN CIAMIS 2.38691084 0.02455384 Reject
## 7 KABUPATEN CIANJUR 0.14581791 0.88518978 Not Reject
## 8 KOTA CIMAHI 0.78232961 0.44108816 Not Reject
## 9 KABUPATEN CIREBON 2.50000880 0.01905709 Reject
## 10 KOTA DEPOK -0.52171919 0.60628225 Not Reject
## 11 KABUPATEN GARUT 1.43095893 0.16435241 Not Reject
## 12 KABUPATEN INDRAMAYU 1.61210952 0.11901235 Not Reject
## 13 KABUPATEN KARAWANG 0.08824923 0.93035482 Not Reject
## 14 KOTA BANDUNG 0.95844420 0.34666723 Not Reject
## 15 KOTA BEKASI -0.36332281 0.71930117 Not Reject
## 16 KOTA BOGOR -0.49981629 0.62140826 Not Reject
## 17 KOTA CIREBON 2.50522014 0.01883357 Reject
## 18 KOTA SUKABUMI -0.29188943 0.77268753 Not Reject
## 19 KOTA TASIKMALAYA 2.02435418 0.05332028 Not Reject
## 20 KABUPATEN KUNINGAN 2.47385067 0.02021659 Reject
## 21 KABUPATEN MAJALENGKA 1.89986795 0.06859580 Not Reject
## 22 KABUPATEN PURWAKARTA 0.38599946 0.70263763 Not Reject
## 23 KABUPATEN SUBANG 0.79287991 0.43502372 Not Reject
## 24 KABUPATEN SUKABUMI -0.57689318 0.56897337 Not Reject
## 25 KABUPATEN SUMEDANG 1.42371711 0.16641930 Not Reject
## 26 KABUPATEN TASIKMALAYA 1.95918102 0.06089662 Not Reject
## 27 KABUPATEN PANGANDARAN 0.25183315 0.80315003 Not Reject