Pemetaan dan Analisis Spasial Kasus Diare di Kota Bandung Menggunakan Moran’s I dan LISA

Bagian A. Konsep Teori

1. Pengertian Autokorelasi Spasial Positif & Negatif

Autokorelasi spasial positif berarti nilai variabel di lokasi-lokasi yang berdekatan cenderung mirip satu sama lain, misalnya wilayah dengan angka kasus tinggi dikelilingi oleh wilayah lain yang juga tinggi (clustering). Autokorelasi spasial negatif berarti nilai yang berdekatan cenderung berbeda pola seperti checkerboard, di mana tinggi-dan-rendah saling berdekatan. (Referensi rumus dan konsep: materi kursus Spatial Statistics).

2. Contoh Fenomena

Positif: Penyakit menular (misal: diare) menyebar melalui kontak lokal sehingga kecamatan berdekatan menunjukkan tingkat kasus tinggi.

Negatif: Zonasi penggunaan lahan yang saling bergantian (mis industri vs perumahan) memberi nilai variabel yang berseberangan pada tetangga sehingga terlihat pola negatif.

3. Rumus Matematis

• Moran’s I (global)

\[ I = \frac{N}{S_0} \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 (global)

\[ C = \frac{(N-1)\sum_{i=1}^{N}\sum_{j=1}^{N} w_{ij}(x_i - x_j)^2} {2 S_0 \sum_{i=1}^{N}(x_i - \bar{x})^2} \]

• Local Moran’s \(I_i\) (LISA)

\[ I_i = \frac{(x_i - \bar{x})}{m^2} \sum_{j=1}^{N} w_{ij}(x_j - \bar{x}), \qquad m^2 = \frac{1}{N} \sum_{k=1}^{N}(x_k - \bar{x})^2 \]

• 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} \]

\[ G_i^*(d) = \frac{\sum_{j=1}^{N} w_{ij}(d) \, x_j} {\sum_{j=1}^{N} x_j} \]

Tabel berikut merangkum perbedaan utama antara ukuran global dan lokal dalam autokorelasi spasial:

4. Perbedaan Ukuran Global dan Lokal

## Warning: package 'knitr' was built under R version 4.3.3

## Warning: package 'kableExtra' was built under R version 4.3.3

Aspek	Ukuran.Global	Ukuran.Lokal
Cakupan	Menghasilkan satu nilai ringkasan untuk seluruh wilayah studi. Menunjukkan kecenderungan umum apakah wilayah dengan nilai tinggi berdekatan dengan wilayah tinggi lainnya, atau sebaliknya.	Menghasilkan nilai statistik untuk setiap unit wilayah (misalnya kecamatan). Menunjukkan apakah suatu unit merupakan bagian dari cluster atau outlier dibandingkan tetangganya.
Tujuan Utama	Menentukan ada atau tidaknya autokorelasi spasial secara keseluruhan di wilayah penelitian.	Mengidentifikasi lokasi spesifik yang membentuk pola spasial, seperti hotspot (nilai tinggi dikelilingi tinggi), coldspot (nilai rendah dikelilingi rendah), atau outlier (tinggi di antara rendah, atau rendah di antara tinggi).
Contoh Ukuran	Moran’s I, Geary’s C	Local Moran’s Ii, Getis–Ord Gi, Getis–Ord Gi*
Interpretasi	Jika signifikan, artinya terdapat pola spasial global (misalnya clustering atau pola papan catur). Namun, tidak memberikan informasi lokasi spesifik.	Jika signifikan, dapat menunjukkan unit wilayah tertentu yang membentuk cluster atau outlier. Misalnya, kecamatan A signifikan High-High, kecamatan B signifikan Low-Low.
Sensitivitas	Lebih menekankan pada kecenderungan rata-rata global. Kurang sensitif terhadap variasi kecil di tingkat lokal.	Lebih sensitif terhadap variasi spasial lokal dan heterogenitas antar wilayah. Mampu mendeteksi pola yang tidak tampak pada ukuran global.
Keterbatasan	Tidak dapat menunjukkan lokasi spesifik cluster; hanya menyatakan ada atau tidaknya autokorelasi secara umum.	Karena dilakukan banyak uji statistik (satu per unit wilayah), rentan terhadap masalah multiple testing sehingga interpretasi p-value perlu hati-hati.

Bagian B. Analisis Data (Simulasi)

1. Data Spasial

Peta Kota Bandung diambil dari shapefile resmi, dan untuk keperluan analisis spasial, dibuat grid simulasi berbentuk persegi yang terdiri dari 30 unit wilayah sebagai representasi kecamatan. Grid ini akan diplot di atas peta asli.

library(sf)
library(dplyr)
library(ggplot2)

setwd("D:/Semester 5/Spasial/Spatial")
Indo_Kec <- readRDS("gadm36_IDN_3_sp.rds")
Bandung <- subset(Indo_Kec, NAME_2 %in% c("Kota Bandung", "Bandung"))
if(nrow(Bandung) == 0){
  stop("Data Kota Bandung tidak ditemukan di RDS!")
}
Bandung$id <- seq_len(nrow(Bandung))

Bandung_sf <- sf::st_as_sf(Bandung)

if(!exists("Bandung_sf")) {
  bb <- sf::st_as_sfc(
    sf::st_bbox(
      c(xmin = 107.55, ymin = -6.98,
        xmax = 107.72, ymax = -6.85),
      crs = 4326
    )
  )
  
  grid <- sf::st_make_grid(bb, n = c(6, 5)) |>
    sf::st_as_sf() |>
    dplyr::mutate(id = dplyr::row_number())
  
  Bandung_sf <- grid
}

Bandung_sf <- sf::st_make_valid(Bandung_sf)
Bandung_sf$id <- if(!"id" %in% names(Bandung_sf)) seq_len(nrow(Bandung_sf)) else Bandung_sf$id
ggplot(Bandung_sf) + geom_sf(fill="grey90", color="white") +
labs(title="Kota Bandung") + theme_minimal()

2. Data Simulasi Kasus Penyakit

Kita buat data simulasi jumlah kasus diare per 10.000 penduduk. Kita menggunakan distribusi Poisson dengan rata-rata berbeda antar wilayah. Jumlah penduduk tiap wilayah dibuat acak (3.000–8.000).

library(sf)
library(dplyr)
library(ggplot2)
library(viridis)

set.seed(123)  # biar reproducible

# Baca Shapefile Level 3 (Kecamatan)
adm3 <- st_read("C:/Users/Fauzan L/Downloads/gadm41_IDN_shp/gadm41_IDN_3.shp", quiet = TRUE)

# Filter Khusus Kota Bandung
bandung <- adm3 %>% filter(NAME_2 == "Kota Bandung")
n_units <- nrow(bandung)   # jumlah kecamatan
n_units

## [1] 30

set.seed(2025)

# Populasi Acak Per Kecamatan (3.000 - 8.000 Penduduk)
bandung <- bandung %>%
  mutate(pop = round(runif(n(), 3000, 8000)))

# Buat Lambda Berbeda Tiap Kecamatan
lambda_base <- runif(n_units, 1, 3)
lambda_base[1:5] <- lambda_base[1:5] + 4    
lambda_base[(n_units-4):n_units] <- lambda_base[(n_units-4):n_units] - 0.8  
lambda_base <- pmax(lambda_base, 0.1)

# Simulasi Kasus
cases <- rpois(n_units, lambda = lambda_base * (bandung$pop/10000))
rate10k <- (cases / bandung$pop) * 10000

bandung <- bandung %>%
  mutate(lambda = lambda_base,
         cases = cases,
         rate10k = rate10k)

# Ringkasan
bandung %>% st_drop_geometry() %>% 
  select(pop, lambda, cases, rate10k) %>% summary()

##       pop           lambda           cases          rate10k      
##  Min.   :3119   Min.   :0.9038   Min.   :0.000   Min.   : 0.000  
##  1st Qu.:4754   1st Qu.:1.3547   1st Qu.:0.000   1st Qu.: 0.000  
##  Median :5546   Median :2.0610   Median :1.000   Median : 1.709  
##  Mean   :5666   Mean   :2.5646   Mean   :1.367   Mean   : 2.333  
##  3rd Qu.:6865   3rd Qu.:2.8484   3rd Qu.:2.000   3rd Qu.: 3.096  
##  Max.   :7762   Max.   :6.3238   Max.   :7.000   Max.   :12.565

Peta Choropleth dari Data Simulasi

ggplot(bandung) +
  geom_sf(aes(fill = rate10k), color = "white", size = 0.3) +
  scale_fill_viridis_c(option = "plasma", name = "Kasus per 10k") +
  labs(
    title = "Peta Simulasi Kasus Diare per 10.000 Penduduk",
    subtitle = "Data simulasi untuk 30 kecamatan di Kota Bandung"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 14, face = "bold"),
    plot.subtitle = element_text(size = 11)
  )

Pola Spasial (Visual)

Dari peta simulasi di atas, kita bisa mengamati:

• Wilayah dengan nilai tinggi (warna lebih gelap) terkonsentrasi di beberapa bagian grid.
• Ada indikasi pola cluster (wilayah bertetangga dengan rate tinggi).
• Ada juga wilayah dengan nilai rendah yang dikelilingi wilayah tinggi (indikasi outlier).

Bagian C : Pengukuran Autokorelasi

Hitung Moran’s I untuk data simulasi kasus diare di Kota Bandung.

library(spdep)     # untuk moran, geary, localmoran
library(sf)
library(dplyr)
library(ggplot2)
library(viridis)

# Neighbors & Weights
nb <- poly2nb(bandung)              # daftar tetangga antar kecamatan
lw <- nb2listw(nb, style = "W")     # bobot spasial (row-standardized)

# 1. Moran’s I
moran <- moran.test(bandung$rate10k, lw, randomisation = TRUE)
moran

## 
##  Moran I test under randomisation
## 
## data:  bandung$rate10k  
## weights: lw    
## 
## Moran I statistic standard deviate = 2.1334, p-value = 0.01645
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.19011328       -0.03448276        0.01108287

Interpretasi:

• Nilai positif (0.19) menunjukkan ada kecenderungan autokorelasi spasial positif → kecamatan dengan angka kasus tinggi cenderung berdekatan dengan kecamatan lain yang juga tinggi (dan sebaliknya untuk rendah).
• Signifikan (p=0.016) pola ini bukan acak, ada indikasi nyata klaster penyakit.
• Ada wilayah-wilayah yang membentuk klaster risiko penyakit menular, sehingga intervensi kesehatan harus memperhatikan sebaran spasial ini, bukan hanya rata-rata kota secara keseluruhan.

# 2. Geary’s C 
geary <- geary.test(bandung$rate10k, lw, randomisation = TRUE)
geary

## 
##  Geary C test under randomisation
## 
## data:  bandung$rate10k 
## weights: lw   
## 
## Geary C statistic standard deviate = 0.69337, p-value = 0.244
## alternative hypothesis: Expectation greater than statistic
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.91233812        1.00000000        0.01598416

Interpretasi:

• Nilai mendekati 1 → menunjukkan pola relatif acak (tidak terlalu klaster).
• Tidak signifikan (p=0.244) → tidak ada bukti kuat adanya autokorelasi spasial menurut ukuran ini.

**Perbedaan dengan Moran’s I:**

Moran’s I lebih sensitif terhadap **pola global** (klaster luas)Geary’s C lebih sensitif terhadap **perbedaan lokal antar tetangga langsung**.

Jadi, hasil ini menunjukkan pola klaster mungkin ada di tingkat global, tapi perbedaan lokal antar kecamatan tidak terlalu kuat.

# 3. Local Moran’s I (LISA)
local_moran <- localmoran(bandung$rate10k, lw)

bandung$Ii     <- local_moran[,1]  # nilai statistik
bandung$Ii_p   <- local_moran[,5]  # p-value
bandung$lag_rate <- lag.listw(lw, bandung$rate10k)  # rata-rata tetangga

# Kategorisasi cluster
bandung$quad_sig <- NA
bandung$quad_sig[bandung$rate10k >= mean(bandung$rate10k) &
                 bandung$lag_rate >= mean(bandung$lag_rate) &
                 bandung$Ii_p <= 0.05] <- "High-High"

bandung$quad_sig[bandung$rate10k <= mean(bandung$rate10k) &
                 bandung$lag_rate <= mean(bandung$lag_rate) &
                 bandung$Ii_p <= 0.05] <- "Low-Low"

bandung$quad_sig[bandung$rate10k >= mean(bandung$rate10k) &
                 bandung$lag_rate <= mean(bandung$lag_rate) &
                 bandung$Ii_p <= 0.05] <- "High-Low"

bandung$quad_sig[bandung$rate10k <= mean(bandung$rate10k) &
                 bandung$lag_rate >= mean(bandung$lag_rate) &
                 bandung$Ii_p <= 0.05] <- "Low-High"

# Peta LISA
ggplot(bandung) +
  geom_sf(aes(fill = quad_sig), color="white") +
  scale_fill_manual(values=c("High-High"="red",
                           "Low-Low"="blue",
                           "High-Low"="orange",
                           "Low-High"="lightblue"),
                  na.value="grey90",
                  name="Cluster LISA")

Interpretasi:

• Kecamatan dengan kategori High-High = cluster kasus tinggi, hotspot penyakit.
• Kecamatan dengan Low-Low = cluster rendah, wilayah sehat.
• High-Low dan Low-High = outlier (anomali), artinya wilayah dengan nilai berbeda dibandingkan tetangganya.
• Hasil ini membantu mengidentifikasi wilayah fokus intervensi kesehatan

# 4. Getis–Ord Gi* 
gi <- localG(bandung$rate10k, lw)
bandung$GiZ <- as.numeric(gi)

# Hotspot/Coldspot
bandung$Gi_cat <- ifelse(bandung$GiZ > 1.65, "Hotspot",
                  ifelse(bandung$GiZ < -1.65, "Coldspot", "Not Sig"))

# Peta Getis–Ord Gi*
ggplot(bandung) +
  geom_sf(aes(fill = Gi_cat), color="white") +
  scale_fill_manual(values=c("Hotspot"="red",
                             "Coldspot"="blue",
                             "Not Sig"="grey90"),
                    name="Getis-Ord Gi*") +
  labs(title="Peta Hotspot/Coldspot (Getis–Ord Gi*)",
       subtitle="Simulasi kasus diare per 10.000 penduduk") +
  theme_minimal()

Interpretasi:
- Kecamatan dengan Hotspot, signifikan kasus tinggi yang dikelilingi tinggi.
- Kecamatan dengan Coldspot, signifikan rendah yang dikelilingi rendah.
- Perbandingan dengan LISA:
- LISA bisa mendeteksi outlier (High-Low, Low-High).
- Gi* hanya fokus pada hotspot/coldspot, jadi hasilnya bisa lebih sederhana tetapi lebih fokus untuk pemetaan wilayah rawan penyakit.

Bagian D. Diskusi Kritis

1. Pemanfaatan hasil autokorelasi spasial untuk strategi kesehatan

Hasil analisis autokorelasi spasial, baik ukuran global (Moran’s I, Geary’s C) maupun lokal (Local Moran’s I, Getis–Ord Gi*), sangat bermanfaat bagi dinas kesehatan dalam menyusun strategi pencegahan dan intervensi penyakit menular di Kota Bandung.

• Ukuran global memberikan gambaran umum apakah pola penyebaran penyakit bersifat acak atau membentuk klaster. Jika terdeteksi adanya klaster signifikan, maka dinas kesehatan perlu meningkatkan kewaspadaan karena penyebaran penyakit berpotensi saling memengaruhi antar wilayah.
• **Ukuran lokal (LISA dan Gi*) membantu mengidentifikasi hotspot penyakit, yaitu kecamatan dengan tingkat kasus tinggi yang dikelilingi wilayah dengan kasus tinggi pula. Informasi ini penting untuk menentukan prioritas intervensi**, misalnya penambahan fasilitas kesehatan, distribusi obat, program sanitasi, atau penyuluhan kesehatan masyarakat pada wilayah rawan.
• Sebaliknya, deteksi coldspot atau wilayah sehat dapat menjadi contoh praktik baik (best practice) yang bisa direplikasi ke wilayah lain.
  Dengan demikian, analisis spasial ini mendukung strategi kesehatan yang lebih efektif, efisien, dan berbasis bukti (evidence-based policy).

---

2. Keterbatasan analisis autokorelasi spasial

Meskipun berguna, terdapat beberapa keterbatasan yang perlu diperhatikan:

• MAUP (Modifiable Areal Unit Problem):
  Hasil analisis sangat dipengaruhi oleh skala dan cara pengelompokan wilayah. Misalnya, analisis di tingkat kecamatan bisa menghasilkan pola berbeda dibandingkan analisis di tingkat kelurahan. Hal ini dapat menyebabkan bias interpretasi jika tidak dipertimbangkan dengan hati-hati.

• Ukuran bobot spasial:
  Hasil autokorelasi bergantung pada definisi tetangga (rook, queen, atau k-nearest neighbors). Pemilihan bobot yang berbeda bisa menghasilkan hasil yang berbeda pula. Oleh karena itu, analisis sebaiknya diuji dengan beberapa definisi bobot spasial untuk memastikan konsistensi hasil.

• Masalah multiple testing pada analisis lokal:
  Dalam analisis LISA, setiap unit wilayah diuji secara terpisah. Hal ini menimbulkan risiko false positive (menyatakan ada klaster padahal tidak ada) karena banyaknya uji statistik yang dilakukan secara simultan. Oleh karena itu, interpretasi hasil harus hati-hati, misalnya dengan menerapkan koreksi p-value atau mengutamakan pola yang konsisten secara spasial.

---

Kesimpulan:
Analisis autokorelasi spasial memberikan wawasan penting dalam memahami pola penyebaran penyakit menular. Namun, penggunaannya harus mempertimbangkan keterbatasan metodologis agar strategi kesehatan yang dirumuskan benar-benar akurat dan efektif.