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 autokorelasi spasial negatif.
  • Menurut saya Autokorelasi Spasial Positif adalah nilai suatu variabel pada suatu lokasi yang mirip dengan nilai lokasi di sekitarnya nilai yang mirip bisa tinggi dikelilingi tinggi atau rendah dikelilingi rendah. Sedangkan Autokorelasi Spasial Negatif adalah nilai suatu variabel pada suatu lokasi yang berbeda dengan nilai lokasi di sekitarnya.
  1. Berikan contoh sederhana dari fenomena yang bisa menunjukkan masing-masing pola tersebut.

  • Contoh Autokorelasi Spasial Positif (I > 0)

    Harga rumah di kawasan perumahan mewah biasanya merata sama sama mahal (Tinggi dikelilingi Tinggi) jika di kawasan perumahan sederhana biasanya merasa sama sama rendah (Rendah dikelilingi Rendah)

  • Contoh Autokorelasi Spasial Negatif (I < 0)

    Taman kota di tengah pusat kota biasanya di pusat kota memiliki suhu udara yang tinggi karena urban heat island. Sedangkan Taman kota memiliki suhu lebih rendah karena banyak vegetasi.

  1. Tuliskan rumus matematis dari :

  • Moran’s I

    Rumus Moran’s I dituliskan sebagai

    \[ 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

    Rumus Geary’s C adalah

    \[ 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 Ii

    \[ 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. \]

  • Get is-Ord Gi, dan Gi

    Tanpa memasukkan lokasi \(i\) (local \(G_i\)): \[ G_i(d) = \frac{\sum_{j \neq i} w_{ij}(d) x_j}{\sum_{j \neq i} x_j}. \]

    Memasukkan lokasi \(i\) (local \(G_i^*\)):

\[ 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. \]

  1. Jelaskan perbedaan utama antara ukuran global dan lokal

    Ukuran global mengukur tingkat autokorelasi spasial secara keseluruhan untuk seluruh wilayah studi untuk mengetahui kovarians antar lokasi dan perbedaan antar tetangga, sedangkan ukuran lokal mengukur autokorelasi spasial pada setiap lokasi atau unit wilayah untuk mendeteksi hot atau cold spots.

Bagian B. Analisis Data (Simulasi)

  1. Gunakan data spasial kecamatan di Kota Bandung (atau gunakan grid simulasi jika data asli tidak tersedia).

    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)
    ## 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)
    ## Warning: package 'ggplot2' was built under R version 4.3.3
    library(mapview)
    ## Warning: package 'mapview' was built under R version 4.3.3
    library(spData)

    Menggunakan data spasial kecamatan

    # Arahkan ke folder kerja Anda
setwd("/Users/ASUS/Downloads/Spatial")

# Peta administratif Indonesia level kecamatan (gadm36 level 3)
Indo_Kec <- readRDS('gadm36_IDN_3_sp.rds')
# Ambil hanya Kota Bandung
Bandung <- Indo_Kec[Indo_Kec$NAME_2 == "Kota Bandung",]
plot(Bandung)

    Menggunakan grid simulasi

    setwd("/Users/ASUS/Downloads/Spatial")

# ==== menggunakan file RDS ====
Indo_Kec <- readRDS("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)

# ==== Opsi B: DEMO grid jika data eksternal tidak ada ====
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
}
ggplot(Bandung_sf) + geom_sf(fill="grey90", color="white") +
labs(title="Kota Bandung - Unit (contoh)") + theme_minimal()

  2. 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

    library(openxlsx)
    ## Warning: package 'openxlsx' was built under R version 4.3.3
    setwd("/Users/ASUS/Downloads/Spatial")

# Pastikan nama file sesuai dengan yang ada di folder (Cikungunyah_sim.csv)
Cikungunyah <- read.xlsx("Cikungunyah_sim.xlsx", sep = ";")

head(Cikungunyah)

  3. Buat peta choropleth dari data simulasi tersebut

    # === Peta administratif Indonesia level kecamatan (gadm36 level 3) ===
Indo_Kec <- readRDS("gadm36_IDN_3_sp.rds")

# Ambil hanya Kota Bandung
Bandung <- Indo_Kec[Indo_Kec$NAME_2 == "Kota Bandung", ]

# Tambahkan id untuk join (30 kecamatan di Kota Bandung)
Bandung$id <- 1:30

# Ubah ke format sf
Bandung_sf <- st_as_sf(Bandung)

# Gabungkan dengan data kasus (pastikan objek Cikungunyah sudah ada di environment)
Bandung_merged <- Bandung_sf %>%
  left_join(Cikungunyah, by = "id")

# Ambil centroid tiap kecamatan untuk posisi label
centroids <- st_centroid(Bandung_merged)
    ## Warning: st_centroid assumes attributes are constant over geometries
    # === Plot peta jumlah kasus + label nama kecamatan ===
ggplot(Bandung_merged) +
  geom_sf(aes(fill = Cikungunyah), color = "white", size = 0.2) +
  geom_sf_text(data = centroids, aes(label = NAME_3), size = 2.5, color = "black") +
  scale_fill_gradient(low = "yellow", high = "red") +
  theme_bw() +
  theme(
    panel.grid.major = element_blank(),
    panel.grid.minor = element_blank(),
    axis.text.x = element_blank(),
    axis.text.y = element_blank(),
    legend.position = "right"
  ) +
  labs(
    title = "Peta Kasus Cikungunyah Kota Bandung",
    fill  = "Jumlah Kasus"
  )
    ## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data

  4. Apa pola spasial yang terlihat secara visual?

    Peta kasus Cikungunyah di Kota Bandung menunjukkan bahwa penyebaran penyakit tidak merata, melainkan terkonsentrasi di wilayah tengah kota seperti Batununggal dan Antapani yang mencatat jumlah kasus sangat tinggi. Sementara itu, wilayah pinggiran kota seperti Ujung Berung, Cibiru, maupun Bojongloa Kidul cenderung memiliki kasus lebih rendah. Pola spasial ini mengindikasikan adanya klaster kasus di pusat kota yang kemungkinan dipengaruhi oleh kepadatan penduduk dan aktivitas perkotaan yang lebih intens.

Bagian C. Pengukuran Autokorelasi

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

  • Berapa nilai Moran’s I?

    # === Tambah kolom rate10k ke dataset asli ===
Cikungunyah <- read.xlsx("Cikungunyah_sim.xlsx")

# Hitung rate per 10.000 penduduk
Cikungunyah$rate10k <- (Cikungunyah$Cikungunyah / Cikungunyah$Penduduk) * 10000

# Gabungkan dengan shapefile Bandung
Bandung_merged <- Bandung_sf %>%
  left_join(Cikungunyah, by = "id")

# Plot peta rate10k
ggplot(Bandung_merged) +
  geom_sf(aes(fill = rate10k), color="white", size=0.2) +
  scale_fill_viridis_c(option="magma") +
  labs(title="Rate Cikungunyah (per 10.000 Penduduk)", fill="Rate10k") +
  theme_minimal()

    nb <- spdep::poly2nb(sf::as_Spatial(Bandung_sf), queen = TRUE)
lwW <- spdep::nb2listw(nb, style="W")

moran_res <- spdep::moran.test(Bandung_merged$rate10k, lwW,
randomisation = TRUE, alternative = "two.sided")
moran_res
    ## 
##  Moran I test under randomisation
## 
## data:  Bandung_merged$rate10k  
## weights: lwW    
## 
## Moran I statistic standard deviate = 0.22129, p-value = 0.8249
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      -0.009263025      -0.034482759       0.012988897

  • Apakah signifikan secara statistik (uji permutasi)?

    Hasil uji Moran’s I menunjukkan nilai p-value = 0.8249 (> 0.05), sehingga secara statistik tidak signifikan. Artinya, kita gagal menolak hipotesis nol (tidak ada autokorelasi spasial). Dengan demikian, pola spasial rate Cikungunyah di Kota Bandung tidak menunjukkan bukti adanya keterkaitan antar kecamatan.

  • Apa artinya bagi pola spasial penyakit?

    Ketidaksignifikanan Moran’s I berarti distribusi kasus Cikungunyah bersifat acak atau tersebar merata antar kecamatan, tanpa kecenderungan membentuk kluster (hotspot) maupun pola pengelompokan tertentu. Dengan kata lain, penyakit ini tidak memperlihatkan konsentrasi spasial yang terdeteksi pada tingkat kecamatan di Kota Bandung.

  1. Hitung Geary’s C
geary_res <- geary.test(Bandung_merged$rate10k, lwW,
                        randomisation = TRUE, alternative = "two.sided")
geary_res
## 
##  Geary C test under randomisation
## 
## data:  Bandung_merged$rate10k 
## weights: lwW   
## 
## Geary C statistic standard deviate = 0.062913, p-value = 0.9498
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.99258161        1.00000000        0.01390381

  • Bagaimana perbandingannya dengan Moran’s I?

    Hasil uji Geary’s C menunjukkan nilai statistik 0.99 dengan p-value 0.9498 (> 0.05), yang juga tidak signifikan. Sama seperti Moran’s I, Geary’s C menolak adanya autokorelasi spasial pada data rate Cikungunyah di Kota Bandung. Dengan demikian, kedua ukuran konsisten menyimpulkan bahwa tidak ada pola kluster spasial yang berarti.

  • Jelaskan perbedaan sensitivitas kedua ukuran ini

    Moran’s I lebih sensitif dalam menangkap pola global atau kecenderungan umum keterkaitan antar wilayah, sedangkan Geary’s C lebih peka terhadap variasi lokal atau perbedaan nilai antar tetangga terdekat. Perbedaan sensitivitas ini menjelaskan mengapa dalam beberapa kasus keduanya bisa memberi hasil berbeda, meskipun pada data ini keduanya sama-sama tidak mendeteksi autokorelasi spasial.

  1. Hitung Local Moran’s I (LISA)
x <- scale(Bandung_merged$rate10k)[, 1]
lagx <- lag.listw(lwW, x)
lisa <- 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 <- 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 <- bind_cols(Bandung_merged, lisa_df) |>
  mutate(quad = ifelse(Pi.two.sided <= alpha, quad, "Not significant"))

ggplot(Bandung_LISA) +
  geom_sf(aes(fill = quad), color = "white", size = 0.2) +
  geom_sf_text(data = centroids, aes(label = NAME_3), size = 2.5, color = "black") +
  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_bw()
## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data

  • Apa interpretasi hasil ini untuk kasus penyakit menular?

    Hasil peta LISA menunjukkan bahwa sebagian besar kecamatan di Kota Bandung masuk kategori Not significant, sehingga tidak terdapat pola kluster spasial yang berarti. Namun, terdapat satu kecamatan dengan kategori Low-High (Outlier), yaitu wilayah dengan tingkat kasus rendah tetapi dikelilingi oleh tetangga dengan tingkat kasus lebih tinggi. Hal ini menandakan distribusi penyakit secara umum acak, namun ada perbedaan lokal yang menarik untuk ditelusuri lebih lanjut.

  1. Hitung Getis-Ord Gi*

  • Tentukan kecamatan yang termasuk hot spot dan cold spot

    # === Getis-Ord G dan G* ===

# Matriks bobot biner (untuk G dan G*)
lwB <- spdep::nb2listw(nb, style = "B")

# Variabel rate10k
x_raw <- Bandung_merged$rate10k
sum_x <- sum(x_raw)

# Matriks bobot spasial biner (B)
Wb <- spdep::listw2mat(lwB)

# --- Hitung G_i (tanpa i) ---
num_G <- as.numeric(Wb %*% x_raw)     # Σ_{j≠i} w_ij x_j
den_G <- (sum_x - x_raw)              # total tanpa i
G_raw <- num_G / den_G

# --- Hitung G_i* (dengan i) ---
Wb_star <- Wb
diag(Wb_star) <- 1                    # set w_ii = 1
num_Gs <- as.numeric(Wb_star %*% x_raw)
den_Gs <- sum_x
G_star_raw <- num_Gs / den_Gs

# --- Z-score Getis–Ord G* ---
Gz <- spdep::localG(x_raw, listw = lwW)  # z(G_i*)

# Gabungkan ke shapefile
Bandung_G <- Bandung_merged %>%
  mutate(
    G_raw = G_raw,
    G_star_raw = G_star_raw,
    z_Gistar = as.numeric(Gz),
    hotcold = case_when(
      z_Gistar >= 1.96 ~ "Hot spot (p<0.05)",
      z_Gistar <= -1.96 ~ "Cold spot (p<0.05)",
      TRUE ~ "Not significant"
    )
  )

# Ringkasan hasil
summary(dplyr::select(Bandung_G, G_raw, G_star_raw, z_Gistar))
    ##      G_raw           G_star_raw         z_Gistar                geometry 
##  Min.   :0.03344   Min.   :0.04891   Min.   :-1.6569   MULTIPOLYGON :30  
##  1st Qu.:0.11805   1st Qu.:0.13834   1st Qu.:-0.4793   epsg:NA      : 0  
##  Median :0.17843   Median :0.20489   Median : 0.1131   +proj=long...: 0  
##  Mean   :0.17208   Mean   :0.19928   Mean   : 0.1980                     
##  3rd Qu.:0.21809   3rd Qu.:0.26073   3rd Qu.: 0.8150                     
##  Max.   :0.32743   Max.   :0.34506   Max.   : 2.5231
    # --- Visualisasi Getis-Ord G* ---
ggplot(Bandung_G) +
  geom_sf(aes(fill = hotcold), color = "white", size = 0.2) +
  geom_sf_text(data = centroids, aes(label = NAME_3), size = 2.5, color = "black") +
  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-score)", fill = NULL) +
  theme_bw()
    ## Warning in st_point_on_surface.sfc(sf::st_zm(x)): st_point_on_surface may not
## give correct results for longitude/latitude data

  • Bandingkan Hasilnya dengan peta LISA

    Jika dibandingkan dengan peta LISA, pola yang muncul pada Getis–Ord Gi* lebih fokus pada identifikasi hot spot tunggal, sedangkan LISA mampu menunjukkan kombinasi klaster High-High dan Low-Low yang lebih menyebar. Artinya, LISA memberikan informasi spasial yang lebih detail mengenai hubungan lokal antarwilayah, sementara Gi* lebih menekankan intensitas konsentrasi kasus pada satu wilayah tertentu.

  • Apakah ada perbedaan wilayah yang ditandai sebagai klaster signifikan

    Ya, terdapat perbedaan. Pada hasil Gi*, hanya satu kecamatan yang terdeteksi sebagai hot spot signifikan dengan p<0.05, sementara pada LISA terdapat beberapa wilayah lain yang ditandai sebagai klaster signifikan (baik high maupun low). Hal ini menunjukkan bahwa kedua metode memiliki sensitivitas berbeda: LISA lebih sensitif dalam mendeteksi pola hubungan lokal, sedangkan Gi* lebih ketat dalam menandai hot spot tunggal.

Bagian D. Diskusi Kritis

Soal 1

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

Jawaban:

Analisis autokorelasi spasial, seperti Moran’s I atau LISA, membantu dinas kesehatan mengidentifikasi pola pengelompokan kasus penyakit menular di Kota Bandung. Dengan mengetahui lokasi hotspot, intervensi dapat difokuskan pada area berisiko tinggi, sehingga pencegahan, pengendalian vektor, dan distribusi sumber daya kesehatan menjadi lebih tepat sasaran. Hasil ini juga mendukung kebijakan berbasis bukti untuk mengurangi penyebaran penyakit. Menurut Anselin (1995), analisis lokal mampu mengungkap cluster signifikan, sementara Cliff dan Ord (1981) menekankan pentingnya mempertimbangkan faktor spasial agar tidak terjadi bias dalam interpretasi pola penyakit.

Soal 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

Jawaban:

  1. MAUP (Modifiable Areal Unit Problem):
    Hasil autokorelasi spasial sangat sensitif terhadap satuan wilayah analisis. Perubahan unit spasial (misalnya kelurahan vs kecamatan) dapat menghasilkan pola yang berbeda, sehingga interpretasi dapat menimbulkan bias. (Openshaw, 1984).

  2. Ukuran bobot spasial:
    Pemilihan matriks bobot spasial (rook, queen, atau k-nearest neighbors) akan memengaruhi hasil pengukuran autokorelasi. Misalnya, matriks queen contiguity dapat menghasilkan pola yang berbeda dibanding k-nearest neighbors. Tidak ada standar tunggal sehingga hasil bisa subjektif tergantung definisi kedekatan. (Anselin, 1988).

  3. Masalah multiple testing:
    Dalam analisis lokal (misalnya Local Moran’s I), dilakukan banyak uji statistik serentak pada setiap lokasi. Hal ini meningkatkan risiko false positives atau kesalahan tipe I. Oleh karena itu, diperlukan koreksi seperti Bonferroni atau False Discovery Rate (FDR). (Ord & Getis, 1995).