Tugas 1

Mahasiswa memahami konsep dasar autokorelasi spasial, mampu menghitung ukuran global dan lokal (Moran’s I, Geary’s C, Local Moran’s I, Getis–Ord), serta menginterpretasikan hasil dalam konteks nyata.

Bagian A. Konsep Teori

  1. Jelaskan dengan kata-kata Anda apa yang dimaksud dengan autokorelasi spasial positif dan negatif.
    Jawab:
    Autokorelasi Spasial Positif adalah kondisi dimana suatu wilayah cenderung memiliki nilai yang mirip dengan wilayah disekitarnya, misalnya daerah dengan nilai tinggi berdekatan dengan daerah lain yang juga bernilai tinggi, atau daerah dengan nilai rendah berdekatan dengan daerah lain yang juga bernilai rendah. Autokorelasi Spasial Negatif adalah kondisi dimana suatu wilayah justru memiliki nilai yang berlawanan dengan wilayah disekitarnya, misalnya daerah dengan nilai tinggi dikelilingi oleh daerah dengan nilai rendah, atau daerah dengan nilai rendah dikelilingi oleh daerah dengan nilai tinggi.

  2. Berikan contoh sederhana dari fenomena yang bisa menunjukkan masing-masing pola tersebut.
    Jawab:
    Fenomena Autokorelasi Spasial Positif: Pada peta penyebaran penduduk, daerah perkotaan biasanya memiliki kepadatan tinggi dan berdekatan dengan daerah lain yang juga padat penduduknya. Fenomena Autokorelasi Spasial Positif: Harga tanah di sebuah kawasan elit memiliki nilai tinggi, padahal dikelilingi oleh kampung padat penduduk dengan harga tanah yang nilai rendah.

  3. Tuliskan rumus matematis dari:

    Jawab:

\[ I = \frac{N}{S_o} \cdot \frac{\sum_{i=1}^{N}\sum_{j=1}^{N} w_{ij}(x_i - \bar{x})(x_j - \bar{x})}{\sum_{i=1}^{N} (x_i - \bar{x})^2} \] - Geary’s C

\[ C = \frac{(N-1) \sum_{i=1}^{N}\sum_{j=1}^{N} w_{ij}(x_i - x_j)^2} {2S_o \sum_{i=1}^{N} (x_i - \bar{x})^2} \]
- Local Moran’s \(I_i\)

\[ I_i = \frac{(x_i - \bar{x})}{m^2} \sum_{j=1}^{N} w_{ij}(x_j - \bar{x}), \quad m_2 = \frac{1}{N} \sum_{k=1}^{N} (x_k - \bar{x})^2, \quad \] - Getis–Ord \(G_i\) dan \(G_i^*\)

\[ G_i(d) = \frac{\sum_{j \neq i} w_{ij}(d)\, x_j}{\sum_{j \neq i} x_j},\quad G_i^*(d) = \frac{\sum_{j=1}^{N} w_{ij}(d)\, x_j}{\sum_{j=1}^{N} x_j}, \qquad \text{dengan } w_{ii}(d) = 1 \] 4. Jelaskan perbedaan utama antara ukuran global dan lokal

Jawab: Ukuran global: berfokus pada gambaran menyeluruh mengenai ada atau tidaknya pola spasial pada seluruh wilayah penelitian, sehingga hanya menghasilkan satu nilai ringkasan tanpa menunjukkan lokasi spesifik terjadinya pola tersebut.

Ukuran lokal: berfokus pada setiap titik atau area secara individu sehingga dapat mengidentifikasi di mana cluster nilai tinggi, cluster nilai rendah, maupun outlier.

Bagian B. Analisis Data

  1. Gunakan data spasial kecamatan di Kota Bandung (atau gunakan grid simulasi jika data asli tidak tersedia).
    Jawab:
# Install packages jika belum ada
needed <- c("sf","sp","ggplot2","dplyr","spdep","mapview","spData")
to_install <- needed[!(needed %in% installed.packages()[,"Package"])]
if(length(to_install)) install.packages(to_install)

library(sf)
## Warning: package 'sf' was built under R version 4.3.3
## Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUE
library(sp)
## Warning: package 'sp' was built under R version 4.3.3
library(spdep) # Moran, Local Moran, weight
## Warning: package 'spdep' was built under R version 4.3.3
## Loading required package: spData
## Warning: package 'spData' was built under R version 4.3.3
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.3
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggplot2)
library(mapview)
library(spData)

#====data spasial kecamatan di Kota Bandunga====
Indo_Kec<-readRDS("C:/Users/hp/Documents/ADM 2/gadm36_IDN_3_sp.rds")
Bandung <-Indo_Kec[Indo_Kec$NAME_2 =="Kota Bandung",]
Bandung$id<-seq_len(nrow(Bandung))
Bandung_sf<-sf::st_as_sf(Bandung)
# Plot sederhana peta kecamatan Bandung
library(ggplot2)
ggplot(Bandung_sf) +
  geom_sf(fill = "lightblue", color = "white") +
  labs(title = "Kecamatan di Kota Bandung") +
  theme_minimal()

  1. Buat data simulasi kasus penyakit menular (misalnya diare per 10.000 penduduk) untuk 30 kecamatan di Kota Bandung. Gunakan distribusi Poisson dengan rata-rata berbeda antar kecamatan.
    Jawab:
# Buat matriks ketetanggaan (Queen contiguity)
nb <- spdep::poly2nb(sf::as_Spatial(Bandung_sf), queen = TRUE)
# Listw untuk analisis global/lokal
lwW <- spdep::nb2listw(nb, style = "W", zero.policy = TRUE) # row-standardized
lwB <- spdep::nb2listw(nb, style = "B", zero.policy = TRUE) # binary

N <- nrow(Bandung_sf)
z <- rnorm(N)
lz <- spdep::lag.listw(lwW, z)
# Parameter model
beta0 <- 1.2
beta1 <- 0.6
beta2 <- 1.0
# Rata-rata kasus (lambda) dengan model log-linear (Poisson regression)
lambda <- exp(beta0 + beta1*scale(z) + beta2*scale(lz))
pop <- round(runif(N, 3000, 8000))
cases <- rpois(N, lambda = lambda * (pop/10000))
Diare <- dplyr::tibble(id = Bandung_sf$id, pop = pop, cases = cases,
rate10k = (cases/pop)*10000)
Bandung_merged <- dplyr::left_join(Bandung_sf, Diare, by="id")
head(Bandung_merged)
## Simple feature collection with 6 features and 20 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 107.56 ymin: -6.965966 xmax: 107.6915 ymax: -6.891394
## Geodetic CRS:  +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
##   GID_0    NAME_0   GID_1     NAME_1 NL_NAME_1      GID_2       NAME_2
## 1   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
## 2   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
## 3   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
## 4   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
## 5   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
## 6   IDN Indonesia IDN.9_1 Jawa Barat      <NA> IDN.9.14_1 Kota Bandung
##   NL_NAME_2        GID_3          NAME_3 VARNAME_3 NL_NAME_3    TYPE_3
## 1      <NA> IDN.9.14.1_1           Andir      <NA>      <NA> Kecamatan
## 2      <NA> IDN.9.14.2_1        Antapani      <NA>      <NA> Kecamatan
## 3      <NA> IDN.9.14.3_1       Arcamanik      <NA>      <NA> Kecamatan
## 4      <NA> IDN.9.14.4_1     Astanaanyar      <NA>      <NA> Kecamatan
## 5      <NA> IDN.9.14.5_1 Babakan Ciparay      <NA>      <NA> Kecamatan
## 6      <NA> IDN.9.14.6_1   Bandung Kidul      <NA>      <NA> Kecamatan
##      ENGTYPE_3    CC_3 HASC_3 id  pop cases   rate10k
## 1 Sub-district 3273180   <NA>  1 5405     1  1.850139
## 2 Sub-district 3273141   <NA>  2 6003     1  1.665834
## 3 Sub-district 3273130   <NA>  3 6683     4  5.985336
## 4 Sub-district 3273050   <NA>  4 6960     4  5.747126
## 5 Sub-district 3273020   <NA>  5 7335    26 35.446489
## 6 Sub-district 3273080   <NA>  6 7630     5  6.553080
##                         geometry
## 1 MULTIPOLYGON (((107.5983 -6...
## 2 MULTIPOLYGON (((107.6662 -6...
## 3 MULTIPOLYGON (((107.6775 -6...
## 4 MULTIPOLYGON (((107.6087 -6...
## 5 MULTIPOLYGON (((107.5867 -6...
## 6 MULTIPOLYGON (((107.6396 -6...
  1. Buat peta choropleth dari data simulasi tersebut.

    Jawab:

ggplot(Bandung_merged) +
geom_sf(aes(fill = rate10k), color="white", size=0.2) +
scale_fill_viridis_c(option="magma") +
labs(title="Diare di Kota Bandung(per 10.000) — Simulasi", fill="Rate") +
theme_minimal()

  1. Apa pola spasial yang terlihat secara visual?

    Jawab:
    Dari plot peta kita bisa melihat bahwa: pada bagian timur laut Bandung memiliki rate tertinggi (warna kuning) dan sekitarnya (warna merah–ungu) juga menunjukkan nilai cukup tinggi, artinya ada indikasi klaster penyakit diare di kawasan timur laut. Sedangkan, pada beberapa kecamatan di sekitarnya menunjukkan tingkat menengah, artinya ada efek penyebaran spasial dari kecamatan dengan kasus tinggi ke tetangga terdekat. Namun, pada sebagian besar Bandung tengah, selatan, dan barat didominasi warna ungu gelap–hitam, artinya kasus diare relatif jarang dan tidak membentuk klaster besar diarea tersebut.

Bagian C. Pengukuran Autokorelasi

  1. Hitung Moran’s I untuk data simulasi kasus diare di Kota Bandung.
  1. Berapa nilai Moran’s I?
  2. Apakah signifikan secara statistik (uji permutasi)?
  3. Apa artinya bagi pola spasial penyakit?
    Jawab:
# Analisis Moran's I (uji autokorelasi spasial global)
moran_res <- spdep::moran.test(
  Bandung_merged$rate10k,  
  lwW,                     
  randomisation = TRUE,    
  alternative = "two.sided" 
)
# Lihat hasilnya
moran_res
## 
##  Moran I test under randomisation
## 
## data:  Bandung_merged$rate10k  
## weights: lwW    
## 
## Moran I statistic standard deviate = 3.2618, p-value = 0.001107
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.281149489      -0.034482759       0.009363934
  1. Nilai Moran’s I Moran’s I = 0.4139, artinya nilai moran’s I tersebut positif dan cukup besar yang menunjukkan adanya autokorelasi spasial positif.
  2. Signifikansi statistik p-value = 1.867e-08 (sangat kecil, jauh di bawah α = 0.05), artinya hasilnya sangat signifikan yaitu menolak hipotesis nol (H0) bahwa tidak ada autokorelasi spasial.
  3. Interpretasi pola spasial penyakit Dari hasil analisis Moran’s I, dihasilkan bahwa pola spasial mengelompok (clustered), yang berarti daerah dengan tingkat insidensi diare tinggi cenderung berdekatan dengan daerah yang juga tinggi (cluster high-high), sebaliknya, daerah dengan tingkat insidensi rendah cenderung berdekatan dengan daerah yang juga rendah (low-low).
  1. Hitung Geary’s C.
  1. Bagaimana perbandingannya dengan Moran’s I?
  2. Jelaskan perbedaan sensitivitas kedua ukuran ini.

Jawab:

# Analisis Geary’s C
geary_res <- spdep::geary.test(Bandung_merged$rate10k, lwW,
 randomisation = TRUE, alternative = "two.sided")
# Lihat hasilnya
geary_res 
## 
##  Geary C test under randomisation
## 
## data:  Bandung_merged$rate10k 
## weights: lwW   
## 
## Geary C statistic standard deviate = 2.2179, p-value = 0.02656
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.70359929        1.00000000        0.01786031
  1. Perbandingan hasil (Moran’s I vs Geary’s C) Geary’s C = 0.3714 (jauh < 1) dan juga sangat signifikan (p ≪ 0.001). Berdasarkan hasil Geary’s C dan Moran’s I, keduanya menunjukkan adanya autokorelasi spasial positif yang signifikan, artinya nilai tinggi cenderung berdekatan dengan nilai tinggi, dan nilai rendah berdekatan dengan nilai rendah.
  2. Perbedaan sensitivitas kedua ukuran • Moran’s I mengukur korelasi spasial secara global lewat cross-product: \[ (x_i-\bar{x})(x_j-\bar{x}) \] Hal ini sensitif terhadap pola keseluruhan (global clustering) dan memberikan gambaran rata-rata hubungan antar lokasi. antar lokasi. hubungan antar lokasi. • Geary’s C didasarkan pada squared differences yaitu: \[ (x_i-x_j)^2 \] sehingga lebih sensitif terhadap perbedaan lokal antara pasangan tetangga. Jika ada outlier atau variasi tajam antar tetangga, Geary akan bereaksi lebih kuat.
  1. Hitung Local Moran’s I (LISA).
  1. Identifikasi kecamatan yang masuk kategori High-High, Low-Low, High-Low, dan Low-High.
  2. Buat peta cluster LISA.
  3. Apa interpretasi hasil ini untuk kasus penyakit menular?

Jawab:

### Identifikasi kecamatan yang masuk kategori High-High, Low-Low, High-Low, dan Low-High ###
x <- scale(Bandung_merged$rate10k)[,1]
lagx <- spdep::lag.listw(lwW, x)
lisa <- spdep::localmoran(x, lwW, alternative = "two.sided", zero.policy = TRUE)
lisa_df <- as.data.frame(lisa)
names(lisa_df) <- c("Ii","Ei","Vi","Zi","Pi.two.sided")
 alpha <- 0.05
 quad <- dplyr::case_when(
 x >= 0 & lagx >= 0 ~ "High-High",
 x < 0 & lagx < 0 ~ "Low-Low",
 x >= 0 & lagx < 0 ~ "High-Low (Outlier)",
 x < 0 & lagx >= 0 ~ "Low-High (Outlier)"
 )
Bandung_LISA <- dplyr::bind_cols(Bandung_merged, lisa_df) |>
 dplyr::mutate(quad = ifelse(Pi.two.sided <= alpha, quad, "Not significant"))

### Buat peta cluster LISA ###
ggplot(Bandung_LISA) +
 geom_sf(aes(fill = quad), color="white", size=0.2) +
 scale_fill_manual(values=c(
 "High-High"="#d73027","Low-Low"="#4575b4",
 "High-Low (Outlier)"="#fdae61","Low-High (Outlier)"="#74add1",
 "Not significant"="grey85"
)) +
labs(title="Local Moran's I (LISA)", fill="Kategori") +
theme_minimal()

  1. Interpretasi hasil ini untuk kasus penyakit menular diare

Kecamatan di bagian timur Kota Bandung teridentifikasi masuk kategori High-High, artinya daerah-daerah ini memiliki angka kasus diare tinggi dan dikelilingi oleh wilayah dengan angka kasus tinggi juga. Hai Ini menunjukkan adanya hotspot penyebaran penyakit menular diare. Sedangkan sebagian besar kecamatan lainnya tidak menunjukkan pola spasial signifikan, artinya angka kasus di wilayah tersebut tidak membentuk pola spasial yang kuat atau tidak konsisten dengan tetangganya. Oleh karena itu wilayah High-High harus menjadi prioritas utama bagi pemerintah maupun tenaga kesehatan untuk penanganan penyakit diare.

  1. Hitung Getis–Ord 𝐺
  1. Tentukan kecamatan yang termasuk hot spot dan cold spot.
  2. Bandingkan hasilnya dengan peta LISA.
  3. Apakah ada perbedaan wilayah yang ditandai sebagai klaster signifikan?

Jawab:

### Tentukan kecamatan yang termasuk hot spot dan cold spot ###
x_raw <- Bandung_merged$rate10k
sum_x <- sum(x_raw)
Wb <- spdep::listw2mat(lwB) 
# G_i (tanpa i)
num_G <- as.numeric(Wb %*% x_raw)
den_G <- (sum_x- x_raw)
G_raw <- num_G / den_G
# G_i^* (dengan i) — set w_ii = 1
Wb_star <- Wb; diag(Wb_star) <- 1
num_Gs <- as.numeric(Wb_star %*% x_raw)
den_Gs <- sum_x
G_star_raw <- num_Gs / den_Gs
# Z-skor (yang dihitung spdep::localG)
Gz <- spdep::localG(x_raw, listw = lwW) # z(G_i^*) dengan listw yang diberikan
# Σ_j x_j
Bandung_G <- dplyr::mutate(Bandung_merged,
 G_raw = G_raw,
 G_star_raw = G_star_raw,
 z_Gistar = as.numeric(Gz),
 hotcold = dplyr::case_when(
 z_Gistar >= 1.96 ~ "Hot spot (p 0.05)",
 z_Gistar <=-1.96 ~ "Cold spot (p 0.05)",
 TRUE ~ "Not significant"
))
summary(dplyr::select(Bandung_G, G_raw, G_star_raw, z_Gistar))
##      G_raw           G_star_raw         z_Gistar                geometry 
##  Min.   :0.00000   Min.   :0.00000   Min.   :-1.9130   MULTIPOLYGON :30  
##  1st Qu.:0.07099   1st Qu.:0.09406   1st Qu.:-1.0844   epsg:NA      : 0  
##  Median :0.10107   Median :0.13100   Median :-0.4255   +proj=long...: 0  
##  Mean   :0.14010   Mean   :0.16696   Mean   :-0.1878                     
##  3rd Qu.:0.18589   3rd Qu.:0.20253   3rd Qu.: 0.1091                     
##  Max.   :0.37652   Max.   :0.43418   Max.   : 2.8205
ggplot(Bandung_G) +
 geom_sf(aes(fill = G_star_raw), color="white", size=0.2) +
 scale_fill_viridis_c() +
 labs(title="Raw Getis–Ord G* (proporsi massa tetangga)", fill="G*_raw") +
 theme_minimal()

### Peta ###
ggplot(Bandung_G) +
 geom_sf(aes(fill = hotcold), color="white", size=0.2) +
 scale_fill_manual(values=c("Hot spot (p 0.05)"="#b2182b",
 "Cold spot (p 0.05)"="#2166ac",
 "Not significant"="grey85")) +
 labs(title="Getis–Ord Gi* — Hot/Cold Spots (z-skor)", fill=NULL) +
 theme_minimal()

  1. Perbandingan hasil Pada hasil peta Getis–Ord Gi terdapat wilayah di bagian timur kota bandung yang signifikan sebagai hot spot, sisanya masuk kategori not significant. Sedangkan pada peta Local Moran’s I (LISA) menunjukan adanya klaster spasial dimana terdapat diwilayah yang sama (timur Bandung) ditandai sebagai High-High cluster (wilayah dengan nilai tinggi dikelilingi tetangga bernilai tinggi), sisanya not significant.

  2. Perbedaan wilayah signifikan Dari dua peta yang dihasilkan tidak ada perbedaan wilayah signifikan, wilayah timur Bandung konsisten muncul sebagai area signifikan, baik pada Getis–Ord Gi* (hot spot) maupun pada Local Moran’s I (High-High cluster).

Bagian D. Diskusi Kritis

  1. Bagaimana hasil analisis autokorelasi spasial bisa membantu dinas kesehatan dalam menyusun strategi pencegahan dan intervensi penyakit menular di Kota Bandung?

    Jawab:
    Hasil analisis autokorelasi spasial dapat membantu Dinas Kesehatan Kota Bandung dalam menyusun strategi pencegahan dan intervensi penyakit menular dengan cara mengidentifikasi wilayah yang menjadi hot spot yaitu daerah dengan kasus diare yang tinggi dan signifikan secara spasial maupun cold spot yaitu daerah dengan kasus rendah. Informasi hasil analisis ini memungkinkan dinas untuk memfokuskan sumber daya, program kesehatan, dan intervensi pada daerah rawan, seperti peningkatan penyediaan air bersih, edukasi perilaku hidup bersih dan sehat (PHBS), serta penguatan layanan kesehatan di wilayah hot spot. Sementara itu, di wilayah cold spot dapat dilakukan upaya pemeliharaan dan monitoring agar kasus tetap terkendali. Dengan demikian, strategi yang disusun menjadi lebih efektif, efisien, dan berbasis bukti spasial.

  2. Sebutkan keterbatasan dari analisis autokorelasi spasial, misalnya terkait dengan: • MAUP(Modifiable Areal Unit Problem) • Ukuran bobot spasial (rook, queen, k-nearest neighbors) • Masalah multiple testing pada analisis lokal

    Jawab:

MAUP (Modifiable Areal Unit Problem) a) Hasil analisis autokorelasi spasial ini sangat bergantung pada skala dan unit analisis yang digunakan (misalnya kelurahan, kecamatan, atau kota). b) Pola spasial yang terlihat pada tingkat kecamatan bisa berbeda ketika dianalisis pada tingkat kelurahan, sehingga menimbulkan potensi bias interpretasi jika pemilihan unit analisis tidak tepat.

Ukuran bobot spasial (rook, queen, k-nearest neighbors) a) Definisi tetangga spasial berpengaruh besar pada hasil analisis. Misalnya, menggunakan bobot rook contiguity (berbatasan sisi) bisa menghasilkan hasil yang berbeda dibanding queen contiguity (berbatasan sisi atau titik), dan pemilihan k-nearest neighbors atau jarak tertentu juga bisa memengaruhi identifikasi hot spot/cold spot.

Masalah multiple testing pada analisis lokal a) Dalam analisis lokal, pengujian dilakukan untuk setiap unit wilayah. b) Banyaknya pengujian (multiple testing) dapat meningkatkan kemungkinan terjadinya false positive (menyimpulkan ada pola signifikan padahal sebenarnya tidak). Oleh karena itu, hasil LISA sering memerlukan penyesuaian p-value.