Persebaran kasus hipertensi
Install packages
library(spdep)## Loading required package: spData## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`## Loading required package: sf## Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUElibrary(sp)
library(raster)
library(spatialreg)## Loading required package: Matrix##
## Attaching package: 'spatialreg'## The following objects are masked from 'package:spdep':
##
## get.ClusterOption, get.coresOption, get.mcOption,
## get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOptionlibrary(sf)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:raster':
##
## intersect, select, union## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(maps)
library(ggmap)## Loading required package: ggplot2## ℹ Google's Terms of Service: <https://mapsplatform.google.com>
## Stadia Maps' Terms of Service: <https://stadiamaps.com/terms-of-service/>
## OpenStreetMap's Tile Usage Policy: <https://operations.osmfoundation.org/policies/tiles/>
## ℹ Please cite ggmap if you use it! Use `citation("ggmap")` for details.library(gstat)
library(ggplot2)
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ lubridate 1.9.3 ✔ tibble 3.2.1
## ✔ purrr 1.0.2 ✔ tidyr 1.3.1
## ✔ readr 2.1.5
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::expand() masks Matrix::expand()
## ✖ tidyr::extract() masks raster::extract()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ purrr::map() masks maps::map()
## ✖ 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 errorslibrary(spgwr)## NOTE: This package does not constitute approval of GWR
## as a method of spatial analysis; see example(gwr)library(spdep)
library(car)## Loading required package: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:purrr':
##
## some
##
## The following object is masked from 'package:dplyr':
##
## recodelibrary(lmtest)## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(viridis)## Loading required package: viridisLite
##
## Attaching package: 'viridis'
##
## The following object is masked from 'package:maps':
##
## unemplibrary(readxl)Import data
Peta Indonesia berdasarkan provinsi
setwd("D:/DIAN/unpad/studi/SEMESTER 5/Spatial Statistics/Tugas/paper")
Indo<-st_read("gadm41_IDN.gpkg", layer = "ADM_ADM_1") #sf package, reading .gpkg file## Reading layer `ADM_ADM_1' from data source
## `D:\DIAN\unpad\studi\SEMESTER 5\Spatial Statistics\Tugas\paper\gadm41_IDN.gpkg'
## using driver `GPKG'
## Simple feature collection with 34 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 95.00971 ymin: -11.00761 xmax: 141.0194 ymax: 6.076941
## Geodetic CRS: WGS 84Indo$NAME_1## [1] "Aceh" "Bali" "Bangka Belitung"
## [4] "Banten" "Bengkulu" "Gorontalo"
## [7] "Jakarta Raya" "Jambi" "Jawa Barat"
## [10] "Jawa Tengah" "Jawa Timur" "Kalimantan Barat"
## [13] "Kalimantan Selatan" "Kalimantan Tengah" "Kalimantan Timur"
## [16] "Kalimantan Utara" "Kepulauan Riau" "Lampung"
## [19] "Maluku" "Maluku Utara" "Nusa Tenggara Barat"
## [22] "Nusa Tenggara Timur" "Papua" "Papua Barat"
## [25] "Riau" "Sulawesi Barat" "Sulawesi Selatan"
## [28] "Sulawesi Tengah" "Sulawesi Tenggara" "Sulawesi Utara"
## [31] "Sumatera Barat" "Sumatera Selatan" "Sumatera Utara"
## [34] "Yogyakarta"Indo$id = c(1:34)
indo_sf<-st_as_sf(Indo)
head(indo_sf)## Simple feature collection with 6 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 95.00971 ymin: -8.849262 xmax: 123.5523 ymax: 6.076941
## Geodetic CRS: WGS 84
## GID_1 GID_0 COUNTRY NAME_1 VARNAME_1 NL_NAME_1 TYPE_1
## 1 IDN.1_1 IDN Indonesia Aceh NA NA Propinisi
## 2 IDN.2_1 IDN Indonesia Bali NA NA Propinisi
## 3 IDN.3_1 IDN Indonesia Bangka Belitung NA NA Propinisi
## 4 IDN.4_1 IDN Indonesia Banten NA NA Propinisi
## 5 IDN.5_1 IDN Indonesia Bengkulu NA NA Propinisi
## 6 IDN.6_1 IDN Indonesia Gorontalo NA NA Propinisi
## ENGTYPE_1 CC_1 HASC_1 ISO_1 geom id
## 1 Province 11 ID.AC ID-AC MULTIPOLYGON (((98.14903 2.... 1
## 2 Province 51 ID.BA ID-BA MULTIPOLYGON (((115.5257 -8... 2
## 3 Province 19 ID.BB NA MULTIPOLYGON (((107.9614 -3... 3
## 4 Province 36 ID.BT ID-BT MULTIPOLYGON (((106.3868 -6... 4
## 5 Province 17 ID.BE ID-BE MULTIPOLYGON (((103.573 -4.... 5
## 6 Province 75 ID.GO ID-GO MULTIPOLYGON (((123.5491 0.... 6Data Hipertensi
data = read_xlsx("data hipertensi dan faktor.xlsx")
str(data)## tibble [34 × 10] (S3: tbl_df/tbl/data.frame)
## $ Provinsi: chr [1:34] "ACEH" "BALI" "BANGKA BELITUNG" "BANTEN" ...
## $ Y : num [1:34] 21.4 21.7 26 26.8 23.5 26.6 29.5 22.3 32.6 31.3 ...
## $ X1 : num [1:34] 72.6 90.5 88.2 82.8 81.4 ...
## $ X2 : num [1:34] 47.5 32.4 46.9 40 35.1 42.7 55.7 50.8 33.7 30.4 ...
## $ X3 : num [1:34] 18 6.2 12.9 15.8 21.6 30.7 10.9 17.9 15.4 5.1 ...
## $ X4 : num [1:34] 0.2 9.3 1.1 1 0.5 4.2 1.5 0.2 1.6 1.4 ...
## $ X5 : num [1:34] 20 17.7 12.3 41.2 22.5 37.4 35.7 18.4 51.1 54.2 ...
## $ X6 : num [1:34] 14.45 4.25 4.52 6.17 14.04 ...
## $ X7 : num [1:34] 9.89 9.74 8.66 9.48 9.35 ...
## $ X8 : num [1:34] 6.03 2.69 4.56 7.52 3.42 3.06 6.53 4.53 7.44 5.13 ...attach(data)
Ket_variabel = read_xlsx("data hipertensi dan faktor.xlsx",sheet = "vars")
Ket_variabel## # A tibble: 9 × 2
## Variabel Keterangan
## <chr> <chr>
## 1 Y Persentase dari populasi yang menderita hipertensi (tekanan darah ti…
## 2 X1 Persentase dari populasi rumah tangga yang memiliki akses terhadap f…
## 3 X2 Persentase dari populasi yang kurang beraktivitas fisik
## 4 X3 Persentase dari populasi yang tidak konsumsi buah atau sayur segar
## 5 X4 Persentase dari populasi yang memiliki kebiasaan mengonsumsi minuman…
## 6 X5 Persentase dari populasi yang memiliki kebiasaan mengonsumsi makanan…
## 7 X6 Persentase penduduk miskin
## 8 X7 Rata-rata lama sekolah (tahun)
## 9 X8 Persentase tingkat pengangguran terbukasummary(data)## Provinsi Y X1 X2
## Length:34 Min. :20.80 Min. :31.78 Min. :27.80
## Class :character 1st Qu.:24.50 1st Qu.:75.30 1st Qu.:36.10
## Mode :character Median :26.55 Median :79.13 Median :42.65
## Mean :27.09 Mean :77.58 Mean :42.57
## 3rd Qu.:29.18 3rd Qu.:84.12 3rd Qu.:49.25
## Max. :38.70 Max. :90.54 Max. :55.70
## X3 X4 X5 X6
## Min. : 5.00 Min. : 0.200 Min. :12.20 Min. : 4.250
## 1st Qu.:10.35 1st Qu.: 1.100 1st Qu.:21.75 1st Qu.: 6.240
## Median :14.10 Median : 2.000 Median :27.40 Median : 8.425
## Mean :14.28 Mean : 3.191 Mean :28.52 Mean :10.089
## 3rd Qu.:17.98 3rd Qu.: 4.150 3rd Qu.:33.20 3rd Qu.:12.252
## Max. :30.70 Max. :15.200 Max. :54.20 Max. :26.030
## X7 X8
## Min. : 7.340 Min. :2.270
## 1st Qu.: 8.675 1st Qu.:3.487
## Median : 9.285 Median :4.320
## Mean : 9.303 Mean :4.614
## 3rd Qu.: 9.852 3rd Qu.:5.763
## Max. :11.420 Max. :7.520sd = data.frame(Variabel = c("Y","X1","X2","X3","X4","X5","X6","X7","X8"),
"Standar Deviasi"=c(sd(Y),
sd(X1),
sd(X2),
sd(X3),
sd(X4),
sd(X5),
sd(X6),
sd(X7),
sd(X8)))
sd## Variabel Standar.Deviasi
## 1 Y 3.9082173
## 2 X1 11.6799752
## 3 X2 7.8635608
## 4 X3 6.1920799
## 5 X4 3.3365168
## 6 X5 10.1894322
## 7 X6 5.1835090
## 8 X7 0.8187008
## 9 X8 1.4190650data$id = c(1:34)Mapping kasus
indo_merged <- indo_sf %>% left_join(data, by = "id")
head(indo_merged)## Simple feature collection with 6 features and 22 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 95.00971 ymin: -8.849262 xmax: 123.5523 ymax: 6.076941
## Geodetic CRS: WGS 84
## GID_1 GID_0 COUNTRY NAME_1 VARNAME_1 NL_NAME_1 TYPE_1
## 1 IDN.1_1 IDN Indonesia Aceh NA NA Propinisi
## 2 IDN.2_1 IDN Indonesia Bali NA NA Propinisi
## 3 IDN.3_1 IDN Indonesia Bangka Belitung NA NA Propinisi
## 4 IDN.4_1 IDN Indonesia Banten NA NA Propinisi
## 5 IDN.5_1 IDN Indonesia Bengkulu NA NA Propinisi
## 6 IDN.6_1 IDN Indonesia Gorontalo NA NA Propinisi
## ENGTYPE_1 CC_1 HASC_1 ISO_1 id Provinsi Y X1 X2 X3 X4 X5
## 1 Province 11 ID.AC ID-AC 1 ACEH 21.4 72.59 47.5 18.0 0.2 20.0
## 2 Province 51 ID.BA ID-BA 2 BALI 21.7 90.54 32.4 6.2 9.3 17.7
## 3 Province 19 ID.BB NA 3 BANGKA BELITUNG 26.0 88.17 46.9 12.9 1.1 12.3
## 4 Province 36 ID.BT ID-BT 4 BANTEN 26.8 82.79 40.0 15.8 1.0 41.2
## 5 Province 17 ID.BE ID-BE 5 BENGKULU 23.5 81.41 35.1 21.6 0.5 22.5
## 6 Province 75 ID.GO ID-GO 6 GORONTALO 26.6 82.44 42.7 30.7 4.2 37.4
## X6 X7 X8 geom
## 1 14.45 9.89 6.03 MULTIPOLYGON (((98.14903 2....
## 2 4.25 9.74 2.69 MULTIPOLYGON (((115.5257 -8...
## 3 4.52 8.66 4.56 MULTIPOLYGON (((107.9614 -3...
## 4 6.17 9.48 7.52 MULTIPOLYGON (((106.3868 -6...
## 5 14.04 9.35 3.42 MULTIPOLYGON (((103.573 -4....
## 6 15.15 8.48 3.06 MULTIPOLYGON (((123.5491 0....ggplot() +
geom_sf(data=indo_merged, aes(fill = Y),color="grey",size=0.1) +
theme_bw() +
scale_fill_gradient(low = "yellow", high = "red") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank()) +
theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.key.width = unit(1,"cm"),
legend.key.height = unit(0.3,"cm"),
legend.text = element_text(size = 6),
plot.title = element_text(hjust = 0.5, size = 10),
legend.title = element_text(size=8),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank()) +
labs(title = "Persebaran Kasus Hipertensi di Indonesia 2023",fill = "% Prevalensi") +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, 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
id = c(1:34)
Provinsi = data.frame(id, Indo$NAME_1)
Provinsi## id Indo.NAME_1
## 1 1 Aceh
## 2 2 Bali
## 3 3 Bangka Belitung
## 4 4 Banten
## 5 5 Bengkulu
## 6 6 Gorontalo
## 7 7 Jakarta Raya
## 8 8 Jambi
## 9 9 Jawa Barat
## 10 10 Jawa Tengah
## 11 11 Jawa Timur
## 12 12 Kalimantan Barat
## 13 13 Kalimantan Selatan
## 14 14 Kalimantan Tengah
## 15 15 Kalimantan Timur
## 16 16 Kalimantan Utara
## 17 17 Kepulauan Riau
## 18 18 Lampung
## 19 19 Maluku
## 20 20 Maluku Utara
## 21 21 Nusa Tenggara Barat
## 22 22 Nusa Tenggara Timur
## 23 23 Papua
## 24 24 Papua Barat
## 25 25 Riau
## 26 26 Sulawesi Barat
## 27 27 Sulawesi Selatan
## 28 28 Sulawesi Tengah
## 29 29 Sulawesi Tenggara
## 30 30 Sulawesi Utara
## 31 31 Sumatera Barat
## 32 32 Sumatera Selatan
## 33 33 Sumatera Utara
## 34 34 YogyakartaUji Multikolinieritas
vif=diag(solve(cor(data[,c("X1","X2","X3","X4","X5","X6","X7","X8")])))
vif## X1 X2 X3 X4 X5 X6 X7 X8
## 3.346301 1.848012 1.315307 1.278346 1.573357 2.858372 1.966635 1.738952Regresi Metode OLS
model_ols = lm(Y~X1+X2+X3+X4+X5+X6+X7+X8, data = data)
summary(model_ols)##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.1968 -1.6255 -0.0305 1.5204 9.4329
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 38.53506 11.31253 3.406 0.00223 **
## X1 -0.07110 0.10088 -0.705 0.48744
## X2 0.07185 0.11136 0.645 0.52464
## X3 -0.05491 0.11930 -0.460 0.64931
## X4 0.12852 0.21828 0.589 0.56129
## X5 0.21263 0.07929 2.681 0.01279 *
## X6 -0.31325 0.21009 -1.491 0.14847
## X7 -1.17020 1.10335 -1.061 0.29902
## X8 -0.13623 0.59858 -0.228 0.82182
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.7 on 25 degrees of freedom
## Multiple R-squared: 0.3209, Adjusted R-squared: 0.1036
## F-statistic: 1.477 on 8 and 25 DF, p-value: 0.2155Asumsi Regresi OLS
Asumsi Homoskedastisitas
bptest(model_ols)##
## studentized Breusch-Pagan test
##
## data: model_ols
## BP = 17.077, df = 8, p-value = 0.02932Asumsi Independensi Spasial Residual
# membuat koordinat centroid dari peta
centroid = st_centroid(Indo) #sf package## Warning: st_centroid assumes attributes are constant over geometriesCoordK <- st_coordinates(centroid) #sf package
# matriks bobot
nb = knn2nb(knearneigh(CoordK, k = 7))
listw = nb2listw(nb, style = "W")
# moran test
moran.test(model_ols$residuals, listw)##
## Moran I test under randomisation
##
## data: model_ols$residuals
## weights: listw
##
## Moran I statistic standard deviate = 0.95336, p-value = 0.1702
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.037758514 -0.030303030 0.005096665Asumsi Normalitas
shapiro.test(model_ols$residuals)##
## Shapiro-Wilk normality test
##
## data: model_ols$residuals
## W = 0.95748, p-value = 0.2051Metode GWR
# koordinat centroid dari peta
head(CoordK)## X Y
## [1,] 96.91174 4.2240685
## [2,] 115.13202 -8.3697239
## [3,] 106.55223 -2.4486222
## [4,] 106.10941 -6.4564129
## [5,] 102.33629 -3.5327631
## [6,] 122.37616 0.6868824Mencari Optimal Bandwidth
Adaptive Bisquare
adapt_bisq = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.bisquare,
adapt = TRUE
)## Adaptive q: 0.381966 CV score: 1346.684
## Adaptive q: 0.618034 CV score: 950.7472
## Adaptive q: 0.763932 CV score: 840.8457
## Adaptive q: 0.8466892 CV score: 810.8515
## Adaptive q: 0.9113298 CV score: 753.611
## Adaptive q: 0.9451988 CV score: 722.3746
## Adaptive q: 0.966131 CV score: 680.2376
## Adaptive q: 0.9790678 CV score: 669.2915
## Adaptive q: 0.9848785 CV score: 666.0956
## Adaptive q: 0.9906544 CV score: 663.0093
## Adaptive q: 0.9942241 CV score: 661.2316
## Adaptive q: 0.9964303 CV score: 660.1961
## Adaptive q: 0.9977938 CV score: 659.5819
## Adaptive q: 0.9986365 CV score: 659.2123
## Adaptive q: 0.9991573 CV score: 658.9878
## Adaptive q: 0.9994792 CV score: 658.8506
## Adaptive q: 0.9996781 CV score: 658.7663
## Adaptive q: 0.9998011 CV score: 658.7145
## Adaptive q: 0.9998771 CV score: 658.6825
## Adaptive q: 0.999924 CV score: 658.6628
## Adaptive q: 0.999924 CV score: 658.6628## Warning in gwr.sel(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8, data =
## indo_merged, : Bandwidth converged to upper bound:1Fixed Bisquare
fix_bisq = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.bisquare,
adapt = F
)## Bandwidth: 16.79511 CV score: NA## Warning in optimize(gwr.cv.f, lower = beta1, upper = beta2, maximum = FALSE, :
## NA/Inf replaced by maximum positive value## Bandwidth: 27.14792 CV score: 673.8624
## Bandwidth: 27.14788 CV score: 673.8624
## Bandwidth: 27.14784 CV score: 673.8624
## Bandwidth: 23.19345 CV score: 792.4236
## Bandwidth: 27.11167 CV score: 673.8592
## Bandwidth: 27.11221 CV score: 673.8592
## Bandwidth: 27.11247 CV score: 673.8592
## Bandwidth: 27.11243 CV score: 673.8592
## Bandwidth: 27.11251 CV score: 673.8592
## Bandwidth: 27.11247 CV score: 673.8592Adaptive Gaussian
adapt_gauss = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.Gauss,
adapt = TRUE
)## Adaptive q: 0.381966 CV score: 593.919
## Adaptive q: 0.618034 CV score: 605.1165
## Adaptive q: 0.236068 CV score: 671.0322
## Adaptive q: 0.484272 CV score: 605.5408
## Adaptive q: 0.3262379 CV score: 601.7829
## Adaptive q: 0.4210434 CV score: 599.3685
## Adaptive q: 0.3779433 CV score: 594.6459
## Adaptive q: 0.392118 CV score: 594.9578
## Adaptive q: 0.3858437 CV score: 594.2222
## Adaptive q: 0.3827117 CV score: 593.8879
## Adaptive q: 0.3828828 CV score: 593.9054
## Adaptive q: 0.3824716 CV score: 593.8636
## Adaptive q: 0.3822785 CV score: 593.8646
## Adaptive q: 0.3823852 CV score: 593.8549
## Adaptive q: 0.3823446 CV score: 593.8531
## Adaptive q: 0.3823446 CV score: 593.8531Fixed Gaussian
fix_gauss = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.Gauss
)## Bandwidth: 16.79511 CV score: 630.0657
## Bandwidth: 27.14792 CV score: 642.2421
## Bandwidth: 10.39671 CV score: 635.6551
## Bandwidth: 17.1655 CV score: 630.7356
## Bandwidth: 14.69821 CV score: 626.2836
## Bandwidth: 13.05519 CV score: 624.5621
## Bandwidth: 11.2851 CV score: 628.3901
## Bandwidth: 13.31975 CV score: 624.6502
## Bandwidth: 12.86653 CV score: 624.5688
## Bandwidth: 12.98265 CV score: 624.5573
## Bandwidth: 12.98095 CV score: 624.5573
## Bandwidth: 12.97986 CV score: 624.5573
## Bandwidth: 12.9799 CV score: 624.5573
## Bandwidth: 12.97994 CV score: 624.5573
## Bandwidth: 12.9799 CV score: 624.5573Adapative Tricube
adapt_tric = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.tricube,
adapt = TRUE
)## Adaptive q: 0.381966 CV score: 1434.656
## Adaptive q: 0.618034 CV score: 978.0577
## Adaptive q: 0.763932 CV score: 858.7345
## Adaptive q: 0.8309029 CV score: 842.9178
## Adaptive q: 0.8406317 CV score: 837.7254
## Adaptive q: 0.9015049 CV score: 792.7834
## Adaptive q: 0.9391267 CV score: 762.8071
## Adaptive q: 0.9623782 CV score: 718.6995
## Adaptive q: 0.9767485 CV score: 700.4097
## Adaptive q: 0.9856298 CV score: 693.991
## Adaptive q: 0.9964649 CV score: 685.5883
## Adaptive q: 0.9923262 CV score: 688.7375
## Adaptive q: 0.9948841 CV score: 686.7739
## Adaptive q: 0.9978152 CV score: 684.596
## Adaptive q: 0.9986497 CV score: 683.9931
## Adaptive q: 0.9991655 CV score: 683.6246
## Adaptive q: 0.9994842 CV score: 683.3985
## Adaptive q: 0.9996812 CV score: 683.2594
## Adaptive q: 0.999803 CV score: 683.1737
## Adaptive q: 0.9998782 CV score: 683.1208
## Adaptive q: 0.9999248 CV score: 683.0882
## Adaptive q: 0.9999248 CV score: 683.0882## Warning in gwr.sel(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8, data =
## indo_merged, : Bandwidth converged to upper bound:1Fixed Tricube
fix_tric = gwr.sel(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
gweight = gwr.tricube
)## Bandwidth: 16.79511 CV score: NA## Warning in optimize(gwr.cv.f, lower = beta1, upper = beta2, maximum = FALSE, :
## NA/Inf replaced by maximum positive value## Bandwidth: 27.14792 CV score: 690.9788
## Bandwidth: 27.14788 CV score: 690.9793
## Bandwidth: 27.14796 CV score: 690.9784
## Bandwidth: 33.54634 CV score: 672.2615
## Bandwidth: 31.59586 CV score: 679.0926
## Bandwidth: 37.50075 CV score: 664.6169
## Bandwidth: 39.16107 CV score: 662.9698
## Bandwidth: 40.97085 CV score: 661.3451
## Bandwidth: 42.08935 CV score: 660.4533
## Bandwidth: 42.78063 CV score: 659.9525
## Bandwidth: 43.20786 CV score: 659.6624
## Bandwidth: 43.4719 CV score: 659.4889
## Bandwidth: 43.63509 CV score: 659.3834
## Bandwidth: 43.73594 CV score: 659.3189
## Bandwidth: 43.79827 CV score: 659.2793
## Bandwidth: 43.8368 CV score: 659.2549
## Bandwidth: 43.86061 CV score: 659.2398
## Bandwidth: 43.87532 CV score: 659.2305
## Bandwidth: 43.88442 CV score: 659.2248
## Bandwidth: 43.89004 CV score: 659.2213
## Bandwidth: 43.89351 CV score: 659.2191
## Bandwidth: 43.89566 CV score: 659.2177
## Bandwidth: 43.89698 CV score: 659.2169
## Bandwidth: 43.8978 CV score: 659.2164
## Bandwidth: 43.89831 CV score: 659.216
## Bandwidth: 43.89862 CV score: 659.2158
## Bandwidth: 43.89882 CV score: 659.2157
## Bandwidth: 43.89894 CV score: 659.2156
## Bandwidth: 43.89901 CV score: 659.2156
## Bandwidth: 43.89906 CV score: 659.2156
## Bandwidth: 43.89906 CV score: 659.2156## Warning in gwr.sel(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8, data =
## indo_merged, : Bandwidth converged to upper bound:43.8991300774903bandwidth optimal untuk tiap fungsi kernel
bandwidth = data.frame("Fungsi Kernel" = c("Adaptive Bisquare",
"Fixed Bisquare",
"Adaptive Gaussian",
"Fixed Gaussian",
"Adaptive Tricube",
"Fixed Tricube"),
Bandwidth = c(adapt_bisq[1],
fix_bisq[1],
adapt_gauss[1],
fix_gauss[1],
adapt_tric[1],
fix_tric[1]),
"CV Score" = c(658.6628,
673.8592,
593.8531,
624.5573,
683.0882,
659.2156))
bandwidth## Fungsi.Kernel Bandwidth CV.Score
## 1 Adaptive Bisquare 0.9999240 658.6628
## 2 Fixed Bisquare 27.1124694 673.8592
## 3 Adaptive Gaussian 0.3823446 593.8531
## 4 Fixed Gaussian 12.9798999 624.5573
## 5 Adaptive Tricube 0.9999248 683.0882
## 6 Fixed Tricube 43.8990561 659.2156Model GWR untuk setiap fungsi Kernel
model_gwr1 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
adapt = adapt_bisq,
gweight = gwr.bisquare,
hatmatrix = TRUE,
se.fit = TRUE
)model_gwr2 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
bandwidth = fix_bisq,
gweight = gwr.bisquare,
hatmatrix = TRUE,
se.fit = TRUE
)model_gwr3 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
adapt = adapt_gauss,
gweight = gwr.Gauss,
hatmatrix = TRUE,
se.fit = TRUE
)model_gwr4 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
bandwidth = fix_gauss,
gweight = gwr.Gauss,
hatmatrix = TRUE,
se.fit = TRUE
)model_gwr5 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
adapt = adapt_tric,
gweight = gwr.tricube,
hatmatrix = TRUE,
se.fit = TRUE
)model_gwr6 = gwr(
Y~X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8,
data = indo_merged,
coords = CoordK,
bandwidth = fix_tric,
gweight = gwr.tricube,
hatmatrix = TRUE,
se.fit = TRUE
)Memilih model terbaik berdasarkan nilai AIC dari tiap model dengan fungsi Kernel berbeda
aic_results = data.frame("Fungsi Kernel atau Bobot" =
c("Adaptive Bisquare",
"Fixed Bisquare",
"Adaptive Gaussian",
"Fixed Gaussian",
"Adaptive Tricube",
"Fixed Tricube"),
Bandwidth = c(adapt_bisq[1],
fix_bisq[1],
adapt_gauss[1],
fix_gauss[1],
adapt_tric[1],
fix_tric[1]),
"CV Score" = c(658.6628,
673.8592,
593.8531,
624.5573,
683.0882,
659.2156),
AIC = c(model_gwr1$results$AICh,
model_gwr2$results$AICh,
model_gwr3$results$AICh,
model_gwr4$results$AICh,
model_gwr5$results$AICh,
model_gwr6$results$AICh))
aic_results## Fungsi.Kernel.atau.Bobot Bandwidth CV.Score AIC
## 1 Adaptive Bisquare 0.9999240 658.6628 172.9896
## 2 Fixed Bisquare 27.1124694 673.8592 173.9039
## 3 Adaptive Gaussian 0.3823446 593.8531 163.8957
## 4 Fixed Gaussian 12.9798999 624.5573 174.4516
## 5 Adaptive Tricube 0.9999248 683.0882 175.9473
## 6 Fixed Tricube 43.8990561 659.2156 182.4531# Model GWR dengan fungsi kernel Adaptive Gaussian memiliki AIC terkecil
model_gwr3## Call:
## gwr(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8, data = indo_merged,
## coords = CoordK, gweight = gwr.Gauss, adapt = adapt_gauss,
## hatmatrix = TRUE, se.fit = TRUE)
## Kernel function: gwr.Gauss
## Adaptive quantile: 0.3823446 (about 12 of 34 data points)
## Summary of GWR coefficient estimates at data points:
## Min. 1st Qu. Median 3rd Qu. Max. Global
## X.Intercept. 24.793771 34.093045 43.473428 50.287971 58.693952 38.5351
## X1 -0.253711 -0.172589 -0.091684 -0.040971 0.021344 -0.0711
## X2 -0.081401 -0.035608 0.051043 0.185563 0.247851 0.0719
## X3 -0.044186 0.021116 0.066345 0.087961 0.178826 -0.0549
## X4 -0.254108 -0.103550 0.057029 0.286263 0.414017 0.1285
## X5 0.026732 0.099021 0.283479 0.347328 0.377729 0.2126
## X6 -0.732479 -0.543263 -0.436375 -0.365251 -0.227090 -0.3133
## X7 -1.606910 -1.485579 -1.320172 -1.059701 -0.878300 -1.1702
## X8 -1.029366 -0.824807 -0.513693 0.445124 0.779286 -0.1362
## Number of data points: 34
## Effective number of parameters (residual: 2traceS - traceS'S): 17.32682
## Effective degrees of freedom (residual: 2traceS - traceS'S): 16.67318
## Sigma (residual: 2traceS - traceS'S): 3.107758
## Effective number of parameters (model: traceS): 14.52962
## Effective degrees of freedom (model: traceS): 19.47038
## Sigma (model: traceS): 2.87587
## Sigma (ML): 2.176291
## AICc (GWR p. 61, eq 2.33; p. 96, eq. 4.21): 209.8121
## AIC (GWR p. 96, eq. 4.22): 163.8957
## Residual sum of squares: 161.0322
## Quasi-global R2: 0.6805216print(paste("AIC Model OLS:",AIC(model_ols)))## [1] "AIC Model OLS: 195.005226570829"print(paste("Adj R-Squared Model OLS:",0.1036))## [1] "Adj R-Squared Model OLS: 0.1036"print(paste("AIC Model GWR (adaptive gaussian):", model_gwr3$results$AICh)) ## [1] "AIC Model GWR (adaptive gaussian): 163.895741796543"print(paste("Quasi Global R-Squared Model OLS:", 0.6805216))## [1] "Quasi Global R-Squared Model OLS: 0.6805216"Asumsi GWR
residual_gwr3 = model_gwr3$SDF$gwr.eNormalitas Residual
shapiro.test(residual_gwr3)##
## Shapiro-Wilk normality test
##
## data: residual_gwr3
## W = 0.95565, p-value = 0.1811Autokorelasi Spasial Residual
nb = knn2nb(knearneigh(CoordK, k = 10))
listw = nb2listw(nb, style = "W")
moran.test(residual_gwr3, listw)##
## Moran I test under randomisation
##
## data: residual_gwr3
## weights: listw
##
## Moran I statistic standard deviate = -0.51303, p-value = 0.696
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## -0.059375523 -0.030303030 0.003211329Homoskedastisitas Spasial
bptest(residual_gwr3~model_gwr3$SDF$X1+model_gwr3$SDF$X2+model_gwr3$SDF$X3+model_gwr3$SDF$X4+model_gwr3$SDF$X5+model_gwr3$SDF$X6+model_gwr3$SDF$X7+model_gwr3$SDF$X8)##
## studentized Breusch-Pagan test
##
## data: residual_gwr3 ~ model_gwr3$SDF$X1 + model_gwr3$SDF$X2 + model_gwr3$SDF$X3 + model_gwr3$SDF$X4 + model_gwr3$SDF$X5 + model_gwr3$SDF$X6 + model_gwr3$SDF$X7 + model_gwr3$SDF$X8
## BP = 7.881, df = 8, p-value = 0.4452Koefisien dan p-values dari model GWR terpilih (Adaptive Gaussian)
Ekstrak koefisien lokal untuk variabel xi
coef_x1 = model_gwr3$SDF$X1
coef_x2 = model_gwr3$SDF$X2
coef_x3 = model_gwr3$SDF$X3
coef_x4 = model_gwr3$SDF$X4
coef_x5 = model_gwr3$SDF$X5
coef_x6 = model_gwr3$SDF$X6
coef_x7 = model_gwr3$SDF$X7
coef_x8 = model_gwr3$SDF$X8Ekstrak standar error untuk variabel xi (untuk tiap provinsi)
se_x1 = model_gwr3$SDF$X1_se
se_x2 = model_gwr3$SDF$X2_se
se_x3 = model_gwr3$SDF$X3_se
se_x4 = model_gwr3$SDF$X4_se
se_x5 = model_gwr3$SDF$X5_se
se_x6 = model_gwr3$SDF$X6_se
se_x7 = model_gwr3$SDF$X7_se
se_x8 = model_gwr3$SDF$X8_set-value untuk variabel xi (untuk tiap provinsi)
t_x1 = coef_x1/se_x1
t_x2 = coef_x2/se_x2
t_x3 = coef_x3/se_x3
t_x4 = coef_x4/se_x4
t_x5 = coef_x5/se_x5
t_x6 = coef_x6/se_x6
t_x7 = coef_x7/se_x7
t_x8 = coef_x8/se_x8p-values untuk variabel xi (untuk tiap provinsi)
n = nrow(data)
k = length(coef(model_gwr3$lm))
df = n-k-1
# p-values untuk masing2 variabel (untuk setiap provinsi)
p_values1 <- 2 * (1 - pt(abs(t_x1), df))
p_values2 <- 2 * (1 - pt(abs(t_x2), df))
p_values3 <- 2 * (1 - pt(abs(t_x3), df))
p_values4 <- 2 * (1 - pt(abs(t_x4), df))
p_values5 <- 2 * (1 - pt(abs(t_x5), df))
p_values6 <- 2 * (1 - pt(abs(t_x6), df))
p_values7 <- 2 * (1 - pt(abs(t_x7), df))
p_values8 <- 2 * (1 - pt(abs(t_x8), df))Satukan p-value dan koefisien dalam dataframe
prov_id = rep(1:34, times=8)
var_id = rep(1:8, each = 34)
# koefisien
coef = data.frame(ID = prov_id,
Variabel_X_ke = var_id,
coef = c(coef_x1, coef_x2, coef_x3, coef_x4, coef_x5, coef_x6, coef_x7, coef_x8))
head(coef)## ID Variabel_X_ke coef
## 1 1 1 -0.069841138
## 2 2 1 -0.206204509
## 3 3 1 0.010916377
## 4 4 1 -0.069621781
## 5 5 1 -0.008996719
## 6 6 1 -0.098741250# p-values
p_values = data.frame(ID = prov_id,
Variabel_X_ke = var_id,
p_value = c(p_values1,p_values2,p_values3,p_values4,p_values5,p_values6,p_values7,p_values8))
head(p_values)## ID Variabel_X_ke p_value
## 1 1 1 0.49908732
## 2 2 1 0.03493772
## 3 3 1 0.93750017
## 4 4 1 0.55248516
## 5 5 1 0.94417621
## 6 6 1 0.28315712# Untuk mapping signifikansi
indo_pv = data.frame(indo_merged,CoordK,p_values1,p_values2,p_values3,p_values4,p_values5,p_values6,p_values7,p_values8)
indo_pv_sf =st_as_sf(indo_pv)
# p-values & koefisien
pvalues_coef = merge(p_values, coef, by = c("ID","Variabel_X_ke"))
pvalues_coef## ID Variabel_X_ke p_value coef
## 1 1 1 4.990873e-01 -0.069841138
## 2 1 2 1.315215e-01 0.155561370
## 3 1 3 8.853604e-01 -0.015373173
## 4 1 4 1.980986e-01 0.251577214
## 5 1 5 1.113882e-04 0.304777176
## 6 1 6 6.011133e-02 -0.421431141
## 7 1 7 1.228030e-01 -1.452877682
## 8 1 8 2.630824e-01 -0.590993400
## 9 10 1 8.218157e-02 -0.184387541
## 10 10 2 3.588940e-01 0.101988546
## 11 10 3 3.118213e-01 0.117396996
## 12 10 4 6.454652e-01 0.106460310
## 13 10 5 9.114211e-05 0.343764200
## 14 10 6 7.873489e-03 -0.694671760
## 15 10 7 1.226507e-01 -1.490373668
## 16 10 8 1.304081e-01 -0.888180560
## 17 11 1 3.576039e-02 -0.207797428
## 18 11 2 6.110131e-01 0.050885780
## 19 11 3 2.861997e-01 0.114553998
## 20 11 4 8.176520e-01 0.046323672
## 21 11 5 3.282334e-04 0.289722872
## 22 11 6 6.500685e-03 -0.642339861
## 23 11 7 1.725614e-01 -1.281028603
## 24 11 8 3.207422e-01 -0.531517772
## 25 12 1 1.356183e-01 -0.151764738
## 26 12 2 2.638074e-01 0.112161390
## 27 12 3 4.260717e-01 0.084792542
## 28 12 4 3.979136e-01 0.167879688
## 29 12 5 9.809270e-05 0.318681611
## 30 12 6 4.299468e-03 -0.698778818
## 31 12 7 1.043325e-01 -1.530000051
## 32 12 8 1.871984e-01 -0.707812058
## 33 13 1 1.615258e-02 -0.253711239
## 34 13 2 8.043060e-01 -0.025651124
## 35 13 3 1.148887e-01 0.178825980
## 36 13 4 6.989166e-01 -0.077076233
## 37 13 5 4.369543e-03 0.228223581
## 38 13 6 4.581402e-03 -0.732479045
## 39 13 7 2.492651e-01 -1.104416465
## 40 13 8 6.824642e-01 -0.223672239
## 41 14 1 4.011181e-02 -0.202984869
## 42 14 2 5.978015e-01 0.051200888
## 43 14 3 2.993400e-01 0.109352982
## 44 14 4 7.198230e-01 0.067733582
## 45 14 5 4.283673e-04 0.277235612
## 46 14 6 3.842032e-03 -0.699227212
## 47 14 7 1.459265e-01 -1.359315857
## 48 14 8 3.457553e-01 -0.495869217
## 49 15 1 5.825303e-02 -0.180205479
## 50 15 2 9.031963e-01 0.011731037
## 51 15 3 5.212233e-01 0.065834442
## 52 15 4 9.262699e-01 -0.017075097
## 53 15 5 8.582670e-03 0.189160854
## 54 15 6 9.183865e-03 -0.534149176
## 55 15 7 2.033504e-01 -1.198711448
## 56 15 8 8.953820e-01 -0.066365712
## 57 16 1 1.081812e-01 -0.142909181
## 58 16 2 6.558937e-01 0.041031999
## 59 16 3 8.934739e-01 0.013099540
## 60 16 4 8.086624e-01 0.043306607
## 61 16 5 4.598600e-03 0.198682424
## 62 16 6 1.415789e-02 -0.462760398
## 63 16 7 1.682816e-01 -1.268709702
## 64 16 8 7.865299e-01 -0.132322609
## 65 17 1 8.310181e-01 -0.026390002
## 66 17 2 6.787709e-02 0.212013775
## 67 17 3 7.341862e-01 0.039206728
## 68 17 4 1.692737e-01 0.332997158
## 69 17 5 4.633050e-05 0.346025969
## 70 17 6 5.187295e-02 -0.480960725
## 71 17 7 1.096996e-01 -1.552518394
## 72 17 8 1.551212e-01 -0.799751378
## 73 18 1 8.930578e-01 -0.017477090
## 74 18 2 7.245523e-02 0.221852347
## 75 18 3 5.694250e-01 0.071688798
## 76 18 4 2.282360e-01 0.329041245
## 77 18 5 3.700135e-05 0.365494537
## 78 18 6 8.722957e-02 -0.452966744
## 79 18 7 1.182704e-01 -1.526257260
## 80 18 8 1.117054e-01 -0.939784081
## 81 19 1 5.975435e-01 -0.046561737
## 82 19 2 6.828558e-01 -0.041252805
## 83 19 3 8.898017e-01 -0.014968297
## 84 19 4 5.605692e-01 -0.112374724
## 85 19 5 4.498003e-01 0.056355426
## 86 19 6 2.194532e-01 -0.230787444
## 87 19 7 3.231382e-01 -0.955897353
## 88 19 8 2.763076e-01 0.583786998
## 89 2 1 3.493772e-02 -0.206204509
## 90 2 2 8.969165e-01 0.012621476
## 91 2 3 2.864891e-01 0.111532693
## 92 2 4 8.716937e-01 -0.031130822
## 93 2 5 2.005407e-03 0.235909262
## 94 2 6 9.409844e-03 -0.566243334
## 95 2 7 2.446983e-01 -1.083148059
## 96 2 8 6.760746e-01 -0.215505287
## 97 20 1 6.443482e-01 -0.042119746
## 98 20 2 6.287561e-01 -0.050709338
## 99 20 3 9.049696e-01 -0.013309644
## 100 20 4 4.646242e-01 -0.146179808
## 101 20 5 7.299632e-01 0.026732493
## 102 20 6 2.340087e-01 -0.228306111
## 103 20 7 2.949617e-01 -1.057917124
## 104 20 8 2.322320e-01 0.662117617
## 105 21 1 6.288205e-02 -0.174006950
## 106 21 2 9.713355e-01 -0.003431708
## 107 21 3 4.374671e-01 0.078813331
## 108 21 4 6.995725e-01 -0.071812900
## 109 21 5 1.137697e-02 0.184373460
## 110 21 6 1.972728e-02 -0.460112286
## 111 21 7 2.870598e-01 -0.981312742
## 112 21 8 8.968252e-01 0.065383414
## 113 22 1 1.724158e-01 -0.122436099
## 114 22 2 6.930179e-01 -0.038927525
## 115 22 3 6.477523e-01 0.047383957
## 116 22 4 4.700083e-01 -0.137267807
## 117 22 5 1.427410e-01 0.107214877
## 118 22 6 6.987922e-02 -0.341211677
## 119 22 7 3.477746e-01 -0.878300348
## 120 22 8 3.902461e-01 0.447299270
## 121 23 1 5.801719e-01 -0.046273278
## 122 23 2 9.453463e-01 0.006378370
## 123 23 3 6.602136e-01 -0.044186167
## 124 23 4 9.946732e-01 -0.001202139
## 125 23 5 8.277394e-02 0.119629059
## 126 23 6 1.804899e-01 -0.240625124
## 127 23 7 2.797980e-01 -0.988056142
## 128 23 8 5.720028e-01 0.280524671
## 129 24 1 6.379276e-01 -0.040588196
## 130 24 2 8.479301e-01 -0.018673965
## 131 24 3 7.422014e-01 -0.034530421
## 132 24 4 7.526026e-01 -0.058915887
## 133 24 5 2.569103e-01 0.081212985
## 134 24 6 2.188631e-01 -0.227090225
## 135 24 7 3.004342e-01 -0.979934700
## 136 24 8 3.965770e-01 0.438597269
## 137 25 1 9.653513e-01 0.005790646
## 138 25 2 6.151215e-02 0.231351779
## 139 25 3 8.023426e-01 0.030526292
## 140 25 4 1.394664e-01 0.391093199
## 141 25 5 6.005346e-05 0.347761508
## 142 25 6 1.452680e-01 -0.381595083
## 143 25 7 1.375932e-01 -1.468244752
## 144 25 8 1.743043e-01 -0.801317276
## 145 26 1 7.857268e-02 -0.175201850
## 146 26 2 6.014312e-01 -0.055134905
## 147 26 3 3.814834e-01 0.096865450
## 148 26 4 3.412753e-01 -0.191306629
## 149 26 5 2.031545e-01 0.096289928
## 150 26 6 2.662818e-02 -0.450916169
## 151 26 7 3.068977e-01 -1.001543214
## 152 26 8 4.016682e-01 0.451136098
## 153 27 1 9.389763e-02 -0.168335589
## 154 27 2 4.964772e-01 -0.073592505
## 155 27 3 3.665261e-01 0.102009717
## 156 27 4 2.446721e-01 -0.240296615
## 157 27 5 3.725970e-01 0.069605276
## 158 27 6 3.742361e-02 -0.421834811
## 159 27 7 3.322886e-01 -0.961779086
## 160 27 8 2.674471e-01 0.612990722
## 161 28 1 1.826552e-01 -0.127402028
## 162 28 2 6.175450e-01 -0.052618428
## 163 28 3 6.302982e-01 0.052855172
## 164 28 4 3.620641e-01 -0.183547928
## 165 28 5 4.117020e-01 0.062155877
## 166 28 6 6.243797e-02 -0.361835603
## 167 28 7 2.839889e-01 -1.065052982
## 168 28 8 2.976965e-01 0.565721081
## 169 29 1 1.924911e-01 -0.127341140
## 170 29 2 4.559860e-01 -0.081400714
## 171 29 3 5.253614e-01 0.072011913
## 172 29 4 2.231218e-01 -0.254108485
## 173 29 5 6.821854e-01 0.032777704
## 174 29 6 7.928193e-02 -0.345424238
## 175 29 7 3.214955e-01 -1.000230961
## 176 29 8 1.693670e-01 0.779286429
## 177 3 1 9.375002e-01 0.010916377
## 178 3 2 5.353740e-02 0.246929047
## 179 3 3 4.898783e-01 0.088499572
## 180 3 4 2.151746e-01 0.379132405
## 181 3 5 3.396013e-05 0.377729285
## 182 3 6 7.566147e-02 -0.488121278
## 183 3 7 1.009759e-01 -1.606910040
## 184 3 8 9.231969e-02 -1.029365652
## 185 30 1 4.486804e-01 -0.068859298
## 186 30 2 6.857872e-01 -0.041881629
## 187 30 3 9.991204e-01 0.000120180
## 188 30 4 4.819976e-01 -0.139546873
## 189 30 5 5.435853e-01 0.045732090
## 190 30 6 1.520711e-01 -0.272329140
## 191 30 7 2.751148e-01 -1.086857489
## 192 30 8 2.917559e-01 0.571416103
## 193 31 1 8.679837e-01 -0.020496472
## 194 31 2 7.503315e-02 0.208353807
## 195 31 3 8.535674e-01 0.021741214
## 196 31 4 1.530113e-01 0.338701168
## 197 31 5 6.093939e-05 0.339151566
## 198 31 6 1.151767e-01 -0.396764172
## 199 31 7 1.335006e-01 -1.458706448
## 200 31 8 1.876321e-01 -0.752658014
## 201 32 1 8.793304e-01 0.021344136
## 202 32 2 5.617693e-02 0.247850512
## 203 32 3 6.066426e-01 0.066856341
## 204 32 4 1.713670e-01 0.411397300
## 205 32 5 4.319594e-05 0.367090391
## 206 32 6 1.481649e-01 -0.394828554
## 207 32 7 1.292485e-01 -1.501431829
## 208 32 8 1.196849e-01 -0.937808679
## 209 33 1 7.057457e-01 -0.043336624
## 210 33 2 9.480941e-02 0.185018825
## 211 33 3 9.883265e-01 0.001643619
## 212 33 4 1.662851e-01 0.296155639
## 213 33 5 7.814434e-05 0.323557802
## 214 33 6 8.739487e-02 -0.413154898
## 215 33 7 1.275509e-01 -1.464120147
## 216 33 8 2.206356e-01 -0.678041001
## 217 34 1 6.698654e-02 -0.189385517
## 218 34 2 3.943100e-01 0.092169017
## 219 34 3 3.208005e-01 0.113025903
## 220 34 4 6.472649e-01 0.102041333
## 221 34 5 1.068201e-04 0.333551915
## 222 34 6 7.776731e-03 -0.674576870
## 223 34 7 1.332413e-01 -1.442085577
## 224 34 8 1.576706e-01 -0.809386092
## 225 4 1 5.524852e-01 -0.069621781
## 226 4 2 1.109261e-01 0.185744356
## 227 4 3 5.363593e-01 0.074749979
## 228 4 4 3.022236e-01 0.256586415
## 229 4 5 3.980069e-05 0.358460493
## 230 4 6 4.263162e-02 -0.516945674
## 231 4 7 1.141219e-01 -1.528444491
## 232 4 8 1.170438e-01 -0.907482931
## 233 5 1 9.441762e-01 -0.008996719
## 234 5 2 6.768520e-02 0.221795831
## 235 5 3 7.174952e-01 0.044390273
## 236 5 4 1.654247e-01 0.358682347
## 237 5 5 4.813125e-05 0.352309055
## 238 5 6 1.228117e-01 -0.397642903
## 239 5 7 1.327762e-01 -1.468772772
## 240 5 8 1.525212e-01 -0.829946805
## 241 6 1 2.831571e-01 -0.098741250
## 242 6 2 7.015697e-01 -0.039506173
## 243 6 3 8.463727e-01 0.020907178
## 244 6 4 4.673133e-01 -0.143777063
## 245 6 5 3.982832e-01 0.062739320
## 246 6 6 9.227770e-02 -0.320327590
## 247 6 7 2.689904e-01 -1.091135838
## 248 6 8 3.365990e-01 0.514875939
## 249 7 1 4.574635e-01 -0.084626918
## 250 7 2 1.240151e-01 0.176109423
## 251 7 3 5.120234e-01 0.078044953
## 252 7 4 3.251398e-01 0.238733972
## 253 7 5 4.144532e-05 0.356658279
## 254 7 6 3.073482e-02 -0.546300331
## 255 7 7 1.094953e-01 -1.541287366
## 256 7 8 1.172394e-01 -0.902729187
## 257 8 1 9.013643e-01 0.017038609
## 258 8 2 5.711979e-02 0.242033678
## 259 8 3 7.064588e-01 0.047496830
## 260 8 4 1.451076e-01 0.414016748
## 261 8 5 5.235317e-05 0.357150779
## 262 8 6 1.591235e-01 -0.375498183
## 263 8 7 1.373729e-01 -1.471195040
## 264 8 8 1.490461e-01 -0.859987252
## 265 9 1 3.248298e-01 -0.109520475
## 266 9 2 1.622878e-01 0.159097926
## 267 9 3 4.669484e-01 0.086346309
## 268 9 4 3.955128e-01 0.204450301
## 269 9 5 4.752358e-05 0.355553881
## 270 9 6 2.232709e-02 -0.581996823
## 271 9 7 1.114089e-01 -1.535906734
## 272 9 8 1.161875e-01 -0.912891961Plot Signifikansi Lokal untuk Setiap Variabel
Plot signifikansi lokal untuk variabel x1
p_value_sign1 <- ifelse(indo_pv_sf$p_values1 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign1), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase dari populasi rumah tangga yang memiliki akses terhadap fasilitas kesehatan dasar") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot signifikansi lokal untuk variabel x2
p_value_sign2 <- ifelse(indo_pv_sf$p_values2 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign2), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase dari populasi yang kurang beraktivitas fisik") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 6)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot signifikansi lokal untuk variabel x3
p_value_sign3 <- ifelse(indo_pv_sf$p_values3 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign3), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase Populasi Makan Sayur dan Buah") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot signifikansi lokal untuk variabel x4
p_value_sign4 <- ifelse(indo_pv_sf$p_values4 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign4), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Koefisien Lokal untuk variabel Persentase dari populasi yang memiliki kebiasaan mengonsumsi minuman beralkohol") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot signifikansi lokal untuk variabel x5
p_value_sign5 <- ifelse(indo_pv_sf$p_values5 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign5), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase dari populasi yang memiliki kebiasaan mengonsumsi makanan berlemak/berkolesterol/gorengan lebih dari 1 kali per hari") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 6.5)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot koefisien signifikansi lokal untuk variabel x6
p_value_sign6 <- ifelse(indo_pv_sf$p_values6 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign6), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase Penduduk Miskin") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot koefisien signifikansi lokal untuk variabel x7
p_value_sign7 <- ifelse(indo_pv_sf$p_values7 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign7), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Rata-Rata Lama Sekolah (tahun)") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
### Plot signifikansi lokal untuk variabel x8
p_value_sign8 <- ifelse(indo_pv_sf$p_values8 < 0.05, "Signifikan", "Tidak Signifikan")
# Plot with manual breaks
ggplot(indo_pv_sf) +
geom_sf(aes(fill = p_value_sign8), color = "black") +
scale_fill_manual(values = c("Signifikan" = "blue3", "Tidak Signifikan" = "grey"),
name = "p-values") +
labs(title = "Peta Signifikansi Lokal untuk variabel Persentase Tingkat Pengagguran Terbuka") +
theme_minimal() +
theme(legend.position = "bottom",
legend.direction = "horizontal")+
theme(plot.title = element_text(size = 7)) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "white")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot koefisien regresi untuk variabel signifikan secara lokal
Plot koefisien regresi lokal untuk variabel x1
indo_coef = data.frame(indo_merged,CoordK,coef_x1,coef_x2,coef_x3,coef_x4,coef_x5,coef_x6,coef_x7,coef_x8)
indo_coef_sf =st_as_sf(indo_coef)
ggplot(indo_coef_sf) +
geom_sf(aes(fill = coef_x1), color = "grey") +
scale_fill_gradient(low = "yellow", high = "red",
limits = c(min(indo_coef$coef_x1), max(indo_coef$coef_x1)),
name = "Koefisien Lokal") +
labs(title = "Peta Koefisien Lokal untuk variabel Persentase dari populasi rumah tangga yang memiliki akses terhadap fasilitas kesehatan dasar ") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank()) +
theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.key.width = unit(1,"cm"),
legend.key.height = unit(0.3,"cm"),
legend.text = element_text(size = 6),
plot.title = element_text(hjust = 0.5, size = 10),
legend.title = element_text(size=8),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank()) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "black")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot koefisien regresi lokal untuk variabel x5
indo_coef = data.frame(indo_merged,CoordK,coef_x1,coef_x2,coef_x3,coef_x4,coef_x5,coef_x6,coef_x7,coef_x8)
indo_coef_sf =st_as_sf(indo_coef)
ggplot(indo_coef_sf) +
geom_sf(aes(fill = coef_x5), color = "grey") +
scale_fill_gradient(low = "yellow", high = "red",
limits = c(min(indo_coef$coef_x5), max(indo_coef$coef_x5)),
name = "Koefisien Lokal") +
labs(title = "Peta Koefisien Lokal untuk variabel Persentase dari populasi yang memiliki kebiasaan mengonsumsi makanan berlemak/berkolesterol/gorengan lebih dari 1 kali per hari") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank()) +
theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.key.width = unit(1,"cm"),
legend.key.height = unit(0.3,"cm"),
legend.text = element_text(size = 6),
plot.title = element_text(hjust = 0.5, size = 10),
legend.title = element_text(size=8),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank()) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "black")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data
Plot koefisien regresi lokal untuk variabel x6
indo_coef = data.frame(indo_merged,CoordK,coef_x1,coef_x2,coef_x3,coef_x4,coef_x5,coef_x6,coef_x7,coef_x8)
indo_coef_sf =st_as_sf(indo_coef)
ggplot(indo_coef_sf) +
geom_sf(aes(fill = coef_x6), color = "grey") +
scale_fill_gradient(low = "yellow", high = "red",
limits = c(min(indo_coef$coef_x6), max(indo_coef$coef_x6)),
name = "Koefisien Lokal") +
labs(title = "Peta Koefisien Lokal untuk variabel Persentase penduduk miskin
") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank()) +
theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.key.width = unit(1,"cm"),
legend.key.height = unit(0.3,"cm"),
legend.text = element_text(size = 6),
plot.title = element_text(hjust = 0.5, size = 10),
legend.title = element_text(size=8),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.ticks = element_blank()) +
geom_sf_text(data = indo_merged, aes(label = id), size = 2, color = "black")+
coord_sf()## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data