Modul Praktikum

Autokorelasi Spasial Univariat dan Bivariat — STA1352 Regresi Spasial

1 Mulai dari sebuah pertanyaan

Bayangkan kita diminta membaca peta kriminalitas suatu kota. Beberapa wilayah tampak memiliki angka tinggi dan letaknya berdekatan. Apakah pengelompokan itu cukup kuat untuk dibedakan dari susunan nilai yang acak? Setelah menjawabnya, kita memperluas pertanyaan: apakah kriminalitas suatu wilayah berkaitan dengan pendapatan wilayah di sekitarnya?

Kita akan mengikuti satu kasus dari awal sampai akhir: 49 lingkungan di Columbus, Ohio, menggunakan data tahun 1980 yang sudah tersedia dalam paket R. Ini contoh historis untuk belajar metode, bukan gambaran kondisi kota saat ini.

Alur berpikir praktikum
Pahami data → lihat peta → tentukan tetangga → buat bobot → uji pola global → cari pola lokal → perluas ke dua variabel → tulis kesimpulan.

install.packages(c("sf", "spdep", "spData", "ggplot2", "knitr", "rmarkdown"))

Internet diperlukan saat memasang paket. Datanya berasal dari paket lokal setelah instalasi. Instalasi sengaja tidak berjalan otomatis ketika Knit.

paket <- c("sf", "spdep", "spData", "ggplot2", "knitr", "rmarkdown")
belum_ada <- paket[!vapply(paket, requireNamespace, logical(1), quietly = TRUE)]
if (length(belum_ada) > 0) {
  stop("Pasang paket berikut melalui Console: ", paste(belum_ada, collapse = ", "))
}

library(sf)
library(spdep)
library(ggplot2)

B <- as.integer(params$nsim)
alpha <- params$alpha
stopifnot(B >= 99, alpha > 0, alpha < 1)
set.seed(2026)
spdep::set.coresOption(NULL)  # eksekusi berurutan; tidak memerlukan cluster

theme_set(theme_minimal(base_size = 12))
fmt <- function(x) formatC(as.numeric(x), digits = 3, format = "f")
fp <- function(p) format.pval(as.numeric(p), digits = 3, eps = 0.0001)
putusan <- function(p) {
  if (p < alpha) "menolak H0" else "belum dapat menolak H0"
}

library() membuka paket; B menyimpan banyaknya permutasi; alpha menjadi ambang pengujian. fmt() dan fp() hanya merapikan tampilan angka. Nilai alpha ditetapkan sebelum melihat hasil.

2 Mengenal data Columbus

Dataset Columbus merupakan data spasial mengenai 49 wilayah lingkungan (neighbourhoods) di Kota Columbus, Ohio, Amerika Serikat, pada tahun 1980. Data ini tersedia dalam paket R spData dan dapat digunakan untuk mempelajari pola suatu variabel serta hubungannya dengan kondisi wilayah di sekitarnya. Dokumentasi data Columbus.

Satu baris data mewakili satu wilayah, bukan satu penduduk atau satu rumah tangga. Setiap wilayah memiliki informasi karakteristik wilayah dan geometri berupa poligon yang menggambarkan batasnya. Jadi, kita dapat mengetahui bukan hanya nilai variabelnya, tetapi juga posisi wilayah dan siapa saja tetangganya.

2.1 Membaca data dan mengenali unit pengamatan

Dalam dataset Columbus, unit pengamatannya adalah wilayah lingkungan (neighbourhood). Artinya, setiap baris data mewakili satu wilayah dan memuat karakteristik wilayah tersebut, seperti angka kriminalitas (CRIME), pendapatan rumah tangga (INC), dan nilai rumah (HOVAL). Karena terdapat 49 wilayah yang diamati, dataset ini memiliki 49 baris.

Sebagai contoh, nilai CRIME, INC, dan HOVAL pada baris pertama semuanya merujuk pada wilayah yang sama. Baris kedua memuat informasi untuk wilayah lain, demikian seterusnya. Dengan demikian, analisis yang kita lakukan membandingkan karakteristik antarwilayah, bukan antarindividu.

Selain informasi dalam tabel, setiap baris juga terhubung dengan geometri berupa poligon, yaitu bentuk batas wilayah pada peta. Geometri ini memungkinkan kita menampilkan nilai variabel dalam peta dan menentukan wilayah mana yang saling bertetangga.

Langkah pertama praktikum adalah membaca data spasial Columbus dari paket spData dan menyimpannya dalam objek bernama col. Selanjutnya, kita memeriksa kelengkapan data, menambahkan nomor urut id untuk memudahkan penyebutan wilayah, dan menampilkan enam baris pertama untuk mengenali susunan datanya.

lokasi_data <- system.file("shapes/columbus.gpkg", package = "spData")
stopifnot(nzchar(lokasi_data))
col <- st_read(lokasi_data, quiet = TRUE)

# Dataset ini memakai koordinat lokal dengan satuan tidak diketahui.
# Simpan metadata aslinya; lepaskan label CRS untuk menggambar bidang lokal.
# Angka koordinat tidak diubah dan tidak direproyeksi ke bujur/lintang.
crs_asli <- st_crs(col)
col <- st_set_crs(col, NA)

# ID mengikuti urutan geometri, bukan urutan nilai CRIME.
col$id <- seq_len(nrow(col))
stopifnot(nrow(col) == 49, all(st_is_valid(col)))
stopifnot(all(is.finite(col$CRIME)), all(is.finite(col$INC)))

knitr::kable(
  head(st_drop_geometry(col)[, c("id", "NEIG", "CRIME", "INC", "HOVAL")]),
  digits = 2, caption = "Enam wilayah pertama; ID dipakai pada tabel dan peta modul."
)

Enam wilayah pertama; ID dipakai pada tabel dan peta modul.

id	NEIG	CRIME	INC	HOVAL
1	5	15.73	19.53	80.47
2	1	18.80	21.23	44.57
3	6	30.63	15.96	26.35
4	2	32.39	4.48	33.20
5	7	50.73	11.25	23.23
6	8	26.07	16.03	28.75

Kolom	Arti	Satuan
CRIME	Pencurian rumah dan kendaraan dalam wilayah	Per 1.000 rumah tangga
INC	Pendapatan rumah tangga	Ribuan dolar AS
HOVAL	Nilai rumah; dipakai pada latihan tambahan	Ribuan dolar AS
geometry	Bentuk dan posisi batas wilayah	Poligon
id	Nomor baris yang kita buat	Penanda, bukan variabel analisis

CRIME tidak mencakup semua jenis kejahatan. Istilah tinggi/rendah di modul bersifat relatif terhadap data ini. Metadata Columbus juga tidak menyediakan dasar yang cukup untuk menafsirkan jarak koordinat sebagai kilometer; karena itu kita memakai persinggungan batas wilayah. Metadata data Columbus.

Sebagai ilustrasi, CRIME = 30 berarti terdapat 30 kejadian pencurian rumah dan kendaraan per 1.000 rumah tangga. Angka ini merupakan ukuran relatif terhadap jumlah rumah tangga, bukan berarti 30% rumah tangga menjadi korban. Sementara itu, INC = 15 menunjukkan nilai pendapatan rumah tangga sebesar 15 ribu dolar AS pada ukuran yang dicatat dalam dataset.

Kita akan menggunakan data ini untuk menjawab dua pertanyaan secara bertahap:

1. Autokorelasi spasial univariat: apakah wilayah dengan CRIME tinggi cenderung bertetangga dengan wilayah yang CRIME-nya juga tinggi?
2. Asosiasi spasial bivariat: apakah wilayah dengan CRIME tinggi cenderung bertetangga dengan wilayah yang memiliki INC rendah?

Pada pertanyaan pertama, kita membandingkan variabel yang sama antarwilayah. Pada pertanyaan kedua, kita membandingkan dua variabel berbeda antarwilayah. Perbedaan inilah yang menjadi dasar perpindahan dari analisis univariat ke bivariat.

summary(st_drop_geometry(col)[, c("CRIME", "INC", "HOVAL")])
#>      CRIME              INC             HOVAL      
#>  Min.   : 0.1783   Min.   : 4.477   Min.   :17.90  
#>  1st Qu.:20.0485   1st Qu.: 9.963   1st Qu.:25.70  
#>  Median :34.0008   Median :13.380   Median :33.50  
#>  Mean   :35.1288   Mean   :14.375   Mean   :38.44  
#>  3rd Qu.:48.5855   3rd Qu.:18.324   3rd Qu.:43.30  
#>  Max.   :68.8920   Max.   :31.070   Max.   :96.40

2.2 Lihat peta sebelum menghitung indeks

Sebelum menghitung autokorelasi spasial, kita memetakan nilai CRIME untuk melihat persebaran kriminalitas antarwilayah. Warna pada peta menunjukkan besar nilainya: semakin mendekati kuning, semakin tinggi angka kriminalitas; semakin mendekati biru tua, semakin rendah.

# Dua peta terpisah menjaga legenda sesuai satuan masing-masing.
ggplot(col) +
  geom_sf(aes(fill = CRIME), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "C", name = "CRIME") +
  labs(title = "Di mana kriminalitas relatif tinggi?", subtitle = "Columbus, 1980") +
  coord_sf(datum = NA) +
  theme_void() + theme(legend.position = "right")

Peta awal: kedua variabel memiliki legenda dan satuan yang berbeda.

Cara membaca: warna menunjukkan besarnya nilai. Wilayah yang tampak berdekatan dan berwarna mirip memberi dugaan pengelompokan. Kita masih memerlukan aturan tetangga dan uji statistik untuk menilai dugaan tersebut.

Secara visual, wilayah dengan angka kriminalitas relatif tinggi tampak berdekatan di bagian tengah peta, sedangkan nilai yang lebih rendah banyak ditemukan di bagian tepi. Pola ini memberi dugaan awal adanya pengelompokan spasial, yang selanjutnya perlu diperiksa melalui uji autokorelasi spasial.

ggplot(col) +
  geom_sf(aes(fill = INC), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "D", name = "INC") +
  labs(title = "Bagaimana pola pendapatan rumah tangga?", subtitle = "Columbus, 1980") +
  coord_sf(datum = NA) +
  theme_void() + theme(legend.position = "right")

Peta menunjukkan bahwa pendapatan rumah tangga (INC) relatif rendah, ditandai warna ungu hingga biru, banyak ditemukan di bagian tengah. Pendapatan relatif tinggi, ditandai warna hijau hingga kuning, tampak pada beberapa wilayah di bagian tepi, terutama sisi kanan peta.

Dibandingkan dengan peta CRIME, terlihat kecenderungan bahwa wilayah berpendapatan rendah memiliki kriminalitas lebih tinggi. Namun, pengamatan visual ini belum membuktikan hubungan statistik, khususnya antara kriminalitas suatu wilayah dan pendapatan wilayah tetangganya.

Berhenti sejenak. Pilih satu bagian peta dengan CRIME relatif tinggi. Apakah tetangganya tampak serupa? Simpan dugaanmu untuk dibandingkan dengan hasil Moran.

3 Apa arti “bertetangga”?

3.1 Queen contiguity dan bobot spasial

Analisis spasial memerlukan definisi ketetanggaan untuk menentukan hubungan antarwilayah. Salah satu pendekatannya adalah Queen contiguity, yaitu dua wilayah dianggap bertetangga apabila berbagi sisi atau titik sudut. Sebagai pembanding, Rook contiguity hanya menganggap dua wilayah bertetangga apabila berbagi sisi.

Dalam praktikum ini, hubungan ketetanggaan ditentukan menggunakan Queen contiguity, kemudian dinyatakan dalam bobot spasial. Bobot tersebut menentukan kontribusi setiap wilayah tetangga dalam perhitungan statistik spasial.

nb <- poly2nb(col, queen = TRUE, row.names = as.character(col$id))
jumlah_tetangga <- card(nb)
summary(jumlah_tetangga)
#>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#>   2.000   3.000   4.000   4.816   6.000  10.000
stopifnot(all(jumlah_tetangga > 0))

lw_B <- nb2listw(nb, style = "B", zero.policy = FALSE)
lw_W <- nb2listw(nb, style = "W", zero.policy = FALSE)
W <- listw2mat(lw_W)

stopifnot(all(diag(W) == 0), all(abs(rowSums(W) - 1) < 1e-10))
data.frame(
  wilayah = nrow(col),
  pasangan_tetangga_unik = sum(jumlah_tetangga) / 2,
  minimum_tetangga = min(jumlah_tetangga),
  maksimum_tetangga = max(jumlah_tetangga)
)
#>   wilayah pasangan_tetangga_unik minimum_tetangga maksimum_tetangga
#> 1      49                    118                2                10

Berdasarkan Queen contiguity, terdapat 118 pasangan wilayah bertetangga pada 49 wilayah Columbus. Setiap wilayah memiliki 2–10 tetangga, dengan rata-rata sekitar 4,82 tetangga. Dengan demikian, seluruh wilayah memiliki tetangga dan tidak terdapat wilayah yang terisolasi. Setiap pasangan dihitung sekali; hubungan wilayah A–B sama dengan B–A.

Daftar ketetanggaan disimpan dalam objek nb, kemudian digunakan untuk membentuk dua jenis bobot:

• Bobot biner (lw_B): bernilai 1 untuk wilayah yang bertetangga dan 0 untuk wilayah yang tidak bertetangga.
• Bobot terstandardisasi per baris (lw_W): bobot setiap tetangga dibagi dengan jumlah tetangga wilayah asal, sehingga total bobot pada setiap baris menjadi 1. Sebagai contoh, jika suatu wilayah memiliki empat tetangga, setiap tetangganya memperoleh bobot ¼.

Matriks W merupakan bentuk matriks dari lw_W. Elemen diagonalnya bernilai nol karena wilayah tidak dianggap sebagai tetangga bagi dirinya sendiri. Pemeriksaan melalui stopifnot() memastikan bahwa kondisi tersebut terpenuhi sebelum analisis dilanjutkan.

titik <- st_coordinates(st_centroid(st_geometry(col)))
plot(st_geometry(col), col = "#f0f5f7", border = "#8195a1", main = "Siapa bertetangga dengan siapa?")
plot(nb, titik, add = TRUE, col = "#25858e", pch = 19, cex = 0.45)

Garis menunjukkan hubungan Queen antarwilayah, bukan jaringan jalan.

Gambar menunjukkan hubungan ketetanggaan antarwilayah berdasarkan Queen contiguity. Setiap titik mewakili pusat geometris (centroid) suatu wilayah, sedangkan garis menghubungkan dua wilayah yang berbagi sisi atau titik sudut.

Jumlah garis yang terhubung pada suatu titik menunjukkan banyaknya tetangga wilayah tersebut. Seluruh wilayah memiliki hubungan ketetanggaan sehingga tidak terdapat wilayah yang terisolasi. Garis pada gambar menunjukkan hubungan spasial, bukan jaringan jalan maupun kemiripan nilai kriminalitas.

3.2 Spatial lag: ringkasan nilai tetangga

Spatial lag digunakan untuk merangkum nilai suatu variabel di wilayah-wilayah tetangga. Pada contoh ini, kita menghitung rata-rata angka kriminalitas tetangga menggunakan bobot yang sudah dibuat.

Karena setiap tetangga mendapat bobot yang sama, perhitungannya cukup dengan merata-ratakan nilai CRIME seluruh tetangga. Nilai CRIME wilayah yang sedang diamati tidak ikut dihitung.

\[ (Wx)_i = \sum_j w_{ij}x_j. \]

x <- col$CRIME
y <- col$INC
n <- length(x)
lag_x <- lag.listw(lw_W, x)

i <- 1L
j <- nb[[i]]
knitr::kable(data.frame(
  id_tetangga = col$id[j], CRIME_tetangga = x[j], bobot = lw_W$weights[[i]]
), digits = 3, caption = paste("Tetangga wilayah ID", col$id[i]))

Tetangga wilayah ID 1

id_tetangga	CRIME_tetangga	bobot
2	18.802	0.5
3	30.627	0.5

c(CRIME_lokasi = x[i],
  rata_rata_manual = mean(x[j]),
  lag_dari_R = lag_x[i])
#>     CRIME_lokasi rata_rata_manual       lag_dari_R 
#>         15.72598         24.71427         24.71427
stopifnot(isTRUE(all.equal(mean(x[j]), lag_x[i])))

Wilayah ID 1 bertetangga dengan wilayah ID 2 dan ID 3. Angka kriminalitas kedua tetangganya adalah 18,802 dan 30,627, sehingga rata-ratanya sekitar 24,714 per 1.000 rumah tangga. Nilai rata-rata inilah yang disebut spatial lag CRIME untuk wilayah ID 1.

Sementara itu, angka kriminalitas wilayah ID 1 sendiri adalah 15,726 per 1.000 rumah tangga. Jadi, angka kriminalitas wilayah ID 1 lebih rendah daripada rata-rata tetangganya. Hasil perhitungan manual dan fungsi lag.listw() menunjukkan nilai yang sama.

4 Autokorelasi spasial global: satu variabel

Autokorelasi spasial global digunakan untuk melihat apakah nilai suatu variabel membentuk pola tertentu di seluruh wilayah yang diamati. Pada data Columbus, kita ingin mengetahui apakah wilayah dengan angka kriminalitas tinggi cenderung berdekatan dengan wilayah yang juga memiliki angka kriminalitas tinggi.

Pola hubungan antarwilayah dapat berupa:

• Autokorelasi positif: wilayah dengan nilai tinggi cenderung bertetangga dengan wilayah bernilai tinggi, begitu pula wilayah bernilai rendah dengan wilayah bernilai rendah.
• Autokorelasi negatif: wilayah dengan nilai tinggi cenderung bertetangga dengan wilayah bernilai rendah, dan sebaliknya.

Istilah global berarti hasil analisis menggambarkan pola secara keseluruhan, belum menunjukkan wilayah mana yang membentuk pola tersebut. Untuk menilai apakah pola yang terlihat didukung secara statistik, kita menggunakan uji autokorelasi spasial.

Beberapa metode yang dapat digunakan adalah Join Count, Moran’s I, Geary’s C, dan Getis–Ord G. Pemilihan metode disesuaikan dengan jenis data dan pertanyaan analisis, seperti dirangkum dalam tabel berikut.

Global berarti kita merangkum seluruh wilayah menjadi sebuah statistik. Berikut pertanyaan yang membedakan metodenya.

Metode	Jenis nilai	Pertanyaan utama	Bagian output yang dibaca
Join Count	Kategori	Apakah kategori sama sering bertetangga?	Jumlah pasangan dan p-value
Moran’s I	Numerik	Apakah nilai serupa mengelompok?	I, nilai harapan, p-value
Geary’s C	Numerik	Apakah selisih nilai tetangga kecil?	C, p-value
Getis–Ord G	Numerik nonnegatif	Apakah nilai tinggi/rendah mengelompok?	G relatif nilai harapan, Z, p-value

4.1 A. Join Count: mulai dari kategori

Join Count digunakan untuk mengetahui apakah wilayah dengan kategori yang sama cenderung bertetangga. Pada contoh ini, nilai CRIME dibagi menjadi dua kategori: Tinggi jika nilainya melebihi median dan Rendah jika nilainya kurang dari atau sama dengan median.

Selanjutnya, kita menghitung pasangan wilayah bertetangga yang sama-sama berkategori Tinggi atau sama-sama berkategori Rendah. Jumlah pasangan tersebut dibandingkan dengan jumlah yang diharapkan jika kategori tersebar secara acak.

Hipotesis yang digunakan untuk setiap kategori adalah:

• H₀: kategori tersebar secara acak, dengan jumlah wilayah pada setiap kategori tetap.
• H₁: pasangan wilayah dengan kategori yang sama lebih banyak daripada yang diharapkan jika tersebar secara acak.

Argumen alternative = "greater" digunakan karena kita ingin menguji apakah jumlah pasangan tersebut lebih banyak dari harapan.

kategori <- factor(
  ifelse(x > median(x), "Tinggi", "Rendah"),
  levels = c("Rendah", "Tinggi")
)
table(kategori)
#> kategori
#> Rendah Tinggi 
#>     25     24
jc <- joincount.test(kategori, lw_B, alternative = "greater", sampling = "nonfree")
jc
#> 
#>  Join count test under nonfree sampling
#> 
#> data:  kategori 
#> weights: lw_B 
#> 
#> Std. deviate for Rendah = 1.1164, p-value = 0.1321
#> alternative hypothesis: greater
#> sample estimates:
#> Same colour statistic           Expectation              Variance 
#>              35.00000              30.10204              19.24708 
#> 
#> 
#>  Join count test under nonfree sampling
#> 
#> data:  kategori 
#> weights: lw_B 
#> 
#> Std. deviate for Tinggi = 6.163, p-value = 3.57e-10
#> alternative hypothesis: greater
#> sample estimates:
#> Same colour statistic           Expectation              Variance 
#>              54.00000              27.69388              18.21936

ringkas_jc <- data.frame(
  kategori = levels(kategori),
  pasangan_teramati = vapply(jc, function(a) unname(a$estimate[1]), numeric(1)),
  harapan_H0 = vapply(jc, function(a) unname(a$estimate[2]), numeric(1)),
  p_value = vapply(jc, function(a) a$p.value, numeric(1))
)
knitr::kable(ringkas_jc, digits = 4, row.names = FALSE)

kategori	pasangan_teramati	harapan_H0	p_value
Rendah	35	30.1020	0.1321
Tinggi	54	27.6939	0.0000

Hasil pembagian berdasarkan median menghasilkan 25 wilayah berkategori Rendah dan 24 wilayah berkategori Tinggi. Pada taraf signifikansi 5%, hasil Join Count menunjukkan:

• Kategori Tinggi: terdapat 54 pasangan wilayah bertetangga, lebih banyak daripada harapan sekitar 27,694 pasangan. Nilai p sebesar 3,57 × 10⁻¹⁰, lebih kecil dari 0,05, sehingga H₀ ditolak. Artinya, wilayah dengan kriminalitas tinggi cenderung mengelompok.
• Kategori Rendah: terdapat 35 pasangan wilayah bertetangga, dibandingkan harapan sekitar 30,102 pasangan. Nilai p sebesar 0,1321, lebih besar dari 0,05, sehingga H₀ belum dapat ditolak. Bukti pengelompokan wilayah dengan kriminalitas rendah belum cukup.

Angka 0,0000 pada tabel merupakan hasil pembulatan p-value yang sangat kecil, bukan berarti p-value tepat nol. Kategori Tinggi dan Rendah pada analisis ini bersifat relatif terhadap median data.

4.2 B. Moran’s I: apakah nilai serupa mengelompok?

Moran’s I digunakan untuk mengetahui apakah wilayah dengan nilai yang mirip cenderung bertetangga. Pada analisis ini, kita menggunakan nilai CRIME asli tanpa membaginya menjadi kategori.

Pertanyaan yang ingin dijawab adalah: apakah wilayah dengan kriminalitas tinggi cenderung berdekatan dengan wilayah yang juga memiliki kriminalitas tinggi, begitu pula untuk nilai yang rendah?

Hipotesis yang digunakan adalah:

• H₀: nilai CRIME tersebar secara acak di antara wilayah.
• H₁: terdapat autokorelasi spasial positif, yaitu nilai CRIME yang mirip cenderung mengelompok.

Pengujian dilakukan melalui dua pendekatan. Fungsi moran.test() menggunakan pendekatan distribusi teoretis, sedangkan moran.mc() membandingkan hasil pengamatan dengan hasil pengacakan nilai CRIME antarlokasi. Argumen alternative = "greater" digunakan untuk menguji pengelompokan nilai yang serupa.

\[ I=\frac{n}{S_0}\frac{\sum_i\sum_j w_{ij}(x_i-\bar{x})(x_j-\bar{x})} {\sum_i(x_i-\bar{x})^2},\qquad S_0=\sum_i\sum_jw_{ij}. \]

moran_analitik <- moran.test(x, lw_W, alternative = "greater", randomisation = TRUE)
set.seed(2026)
moran_perm <- moran.mc(x, lw_W, nsim = B, alternative = "greater")

moran_analitik
#> 
#>  Moran I test under randomisation
#> 
#> data:  x  
#> weights: lw_W    
#> 
#> Moran I statistic standard deviate = 5.5894, p-value = 1.139e-08
#> alternative hypothesis: greater
#> sample estimates:
#> Moran I statistic       Expectation          Variance 
#>       0.500188557      -0.020833333       0.008689289
moran_perm
#> 
#>  Monte-Carlo simulation of Moran I
#> 
#> data:  x 
#> weights: lw_W  
#> number of simulations + 1: 10000 
#> 
#> statistic = 0.50019, observed rank = 10000, p-value = 1e-04
#> alternative hypothesis: greater
I_obs <- unname(moran_perm$statistic)
EI <- -1 / (n - 1)

Hasil analisis menunjukkan Moran’s I sebesar 0,5002, lebih tinggi daripada nilai harapannya jika pola acak, yaitu −0,0208. Hasil ini menunjukkan kecenderungan wilayah dengan angka kriminalitas yang mirip untuk bertetangga.

P-value uji analitik sebesar 1,139 × 10⁻⁸, sedangkan p-value dari 9.999 permutasi sebesar 0,0001. Keduanya lebih kecil dari 0,05, sehingga H₀ ditolak.

Dengan demikian, terdapat bukti autokorelasi spasial positif pada CRIME: secara keseluruhan, wilayah dengan kriminalitas tinggi cenderung berdekatan dengan wilayah bernilai tinggi, dan wilayah bernilai rendah dengan wilayah bernilai rendah. Hasil global ini belum menunjukkan lokasi klaster tertentu; untuk mengetahuinya, diperlukan analisis lokal.

4.2.1 Mengapa memakai permutasi?

Uji permutasi membantu menilai apakah pengelompokan yang ditemukan dapat dijelaskan oleh penempatan nilai secara acak. Caranya, nilai CRIME dipindahkan secara acak antarwilayah, kemudian Moran’s I dihitung kembali. Nilai datanya tetap sama, sedangkan letaknya berubah; hubungan ketetanggaan juga tetap dipertahankan.

Proses ini diulang sebanyak 9.999 kali. Hasilnya digunakan sebagai pembanding untuk mengetahui apakah Moran’s I dari data asli jauh lebih tinggi daripada hasil pengacakan.

sim_I <- moran_perm$res[seq_len(B)]
hist(sim_I, breaks = 30, col = "#acd4dc", border = "white",
     main = "Apakah pola teramati lazim pada peta acak?",
     xlab = "Moran's I dari permutasi",
     xlim = range(c(sim_I, I_obs)))
abline(v = I_obs, col = "#b83248", lwd = 3)
legend("topright", legend = paste("I teramati =", fmt(I_obs)),
       col = "#b83248", lwd = 3, bty = "n")

Garis merah adalah Moran teramati; histogram adalah hasil pengacakan nilai CRIME.

Histogram menunjukkan sebaran Moran’s I dari hasil pengacakan, yang sebagian besar berada di sekitar nol. Garis merah menunjukkan Moran’s I data asli, yaitu sekitar 0,500, yang berada jauh di sebelah kanan sebaran tersebut.

Artinya, pengelompokan pada data asli jauh lebih kuat daripada yang dihasilkan oleh pengacakan. Hasil ini sejalan dengan p-value 0,0001, sehingga pada taraf signifikansi 5% terdapat bukti bahwa wilayah dengan nilai CRIME yang mirip cenderung mengelompok.

4.3 C. Geary’s C: lihat perbedaan antartetangga

Geary’s C digunakan untuk melihat kemiripan nilai antarwilayah yang bertetangga berdasarkan selisih nilainya. Jika angka kriminalitas suatu wilayah dan tetangganya hampir sama, selisihnya kecil dan nilai Geary’s C cenderung rendah.

Nilai C dibaca dengan acuan berikut:

• C < 1: wilayah yang bertetangga cenderung memiliki nilai serupa.
• C > 1: wilayah yang bertetangga cenderung memiliki nilai berbeda.
• C mendekati 1: mendekati nilai yang diharapkan pada pola acak.

Pada contoh ini, H₀ menyatakan nilai CRIME tersebar secara acak, sedangkan H₁ menyatakan wilayah dengan nilai serupa cenderung mengelompok. P-value digunakan untuk menilai apakah bukti pengelompokan tersebut cukup kuat.

geary_hasil <- geary.test(x, lw_W, alternative = "greater", randomisation = TRUE)
geary_hasil
#> 
#>  Geary C test under randomisation
#> 
#> data:  x 
#> weights: lw_W   
#> 
#> Geary C statistic standard deviate = 4.7431, p-value = 1.053e-06
#> alternative hypothesis: Expectation greater than statistic
#> sample estimates:
#> Geary C statistic       Expectation          Variance 
#>       0.540528203       1.000000000       0.009384264
C_obs <- unname(geary_hasil$estimate[1])

Hasil analisis menunjukkan Geary’s C sebesar 0,541, lebih kecil dari nilai harapannya, yaitu 1. Artinya, angka kriminalitas antarwilayah yang bertetangga cenderung serupa.

P-value sebesar 1,053 × 10⁻⁶, lebih kecil dari 0,05, sehingga H₀ ditolak. Dengan demikian, terdapat bukti autokorelasi spasial positif pada CRIME. Kesimpulan ini sejalan dengan hasil Moran’s I sebelumnya, yaitu nilai kriminalitas yang serupa cenderung mengelompok.

Pada fungsi geary.test(), argumen alternative = "greater" digunakan untuk menguji autokorelasi positif, yang ditunjukkan oleh C lebih kecil dari 1. Jadi, argumen tersebut tidak berarti menguji apakah C lebih besar dari 1.

4.4 D. Getis–Ord G global: konsentrasi nilai tinggi

Getis–Ord G global digunakan untuk menilai pengelompokan nilai tinggi atau rendah secara keseluruhan. Pada contoh ini, kita ingin mengetahui apakah wilayah dengan angka kriminalitas tinggi cenderung saling berdekatan.

Analisis menggunakan nilai CRIME asli dan bobot biner, dengan hipotesis:

• H₀: nilai CRIME tersebar secara acak antarwilayah.
• H₁: nilai CRIME yang tinggi mengelompok lebih kuat daripada yang diharapkan pada pola acak.

Argumen alternative = "greater" digunakan karena pengujian difokuskan pada pengelompokan nilai tinggi.

stopifnot(all(x >= 0))
G_global <- globalG.test(x, lw_B, alternative = "greater")
G_global
#> 
#>  Getis-Ord global G statistic
#> 
#> data:  x 
#> weights: lw_B   
#> 
#> standard deviate = 4.6382, p-value = 1.758e-06
#> alternative hypothesis: greater
#> sample estimates:
#> Global G statistic        Expectation           Variance 
#>       1.278075e-01       1.003401e-01       3.507057e-05

Hasil analisis menunjukkan G sebesar 0,1278, lebih besar daripada nilai harapannya pada pola acak, yaitu 0,1003. Nilai Z sebesar 4,638 menunjukkan bahwa G teramati berada di atas nilai harapan tersebut.

P-value sebesar 1,758 × 10⁻⁶, lebih kecil dari 0,05, sehingga H₀ ditolak. Dengan demikian, terdapat bukti bahwa wilayah dengan angka kriminalitas tinggi cenderung saling berdekatan atau mengelompok.

Hasil ini menggambarkan pola secara keseluruhan. Untuk mengetahui wilayah mana yang membentuk pengelompokan tersebut, analisis perlu dilanjutkan dengan statistik lokal.

5 Autokorelasi spasial lokal: di mana polanya?

5.1 Local Moran’s I dan empat kuadran

Analisis global menunjukkan adanya pengelompokan kriminalitas secara keseluruhan. Selanjutnya, Local Moran’s I digunakan untuk menilai pola pada setiap wilayah dengan membandingkan angka kriminalitas wilayah tersebut dan tetangganya.

Sebagai langkah awal, nilai CRIME diubah menjadi skor baku (z). Nilai z positif berarti angka kriminalitas berada di atas rata-rata seluruh wilayah, sedangkan nilai z negatif berarti berada di bawah rata-rata. Kemudian, kita menghitung rata-rata skor baku tetangganya (spatial lag).

Berdasarkan kedua nilai tersebut, setiap wilayah dikelompokkan ke dalam empat kuadran: High–High (HH), Low–Low (LL), High–Low (HL), dan Low–High (LH). Pengelompokan ini membantu mengenali pola awal sebelum dilakukan uji signifikansi.

zx <- as.numeric(scale(x))
zy <- as.numeric(scale(y))
wzx <- lag.listw(lw_W, zx)
wzy <- lag.listw(lw_W, zy)

# Fungsi buatan modul; dipakai kembali untuk kasus bivariat.
kuadran <- function(nilai_lokasi, nilai_tetangga) {
  ifelse(nilai_lokasi == 0 | nilai_tetangga == 0, "Pada batas",
    ifelse(nilai_lokasi > 0 & nilai_tetangga > 0, "HH",
      ifelse(nilai_lokasi < 0 & nilai_tetangga < 0, "LL",
        ifelse(nilai_lokasi > 0 & nilai_tetangga < 0, "HL", "LH"))))
}
q_uni <- kuadran(zx, wzx)
table(q_uni)
#> q_uni
#> HH HL LH LL 
#> 21  3  5 20

Pada grafik, setiap titik mewakili satu wilayah. Sumbu horizontal menunjukkan skor baku CRIME wilayah tersebut, sedangkan sumbu vertikal menunjukkan rata-rata skor baku CRIME tetangganya. Garis nol memisahkan nilai di atas dan di bawah rata-rata.

Hasil pembagian kuadran menunjukkan:

Kuadran	Jumlah wilayah	Interpretasi
HH	21	Kriminalitas wilayah tinggi dan rata-rata kriminalitas tetangganya juga tinggi.
LL	20	Kriminalitas wilayah rendah dan rata-rata kriminalitas tetangganya juga rendah.
HL	3	Kriminalitas wilayah tinggi, tetapi rata-rata kriminalitas tetangganya rendah.
LH	5	Kriminalitas wilayah rendah, tetapi rata-rata kriminalitas tetangganya tinggi.

Sebagian besar wilayah berada pada kuadran HH dan LL, sehingga nilai kriminalitas wilayah cenderung sejalan dengan nilai tetangganya. Pola ini sesuai dengan garis yang menanjak dan hasil Moran’s I global yang positif.

Namun, posisi pada kuadran belum menunjukkan bahwa pola suatu wilayah signifikan. Uji Local Moran diperlukan untuk menilai apakah pola tersebut cukup kuat dibandingkan dengan hasil pengacakan.

ggplot(data.frame(zx, wzx, kuadran = q_uni), aes(zx, wzx, color = kuadran)) +
  geom_hline(yintercept = 0, color = "grey65") +
  geom_vline(xintercept = 0, color = "grey65") +
  geom_point(size = 2.5) +
  geom_abline(intercept = 0, slope = I_obs, color = "#183f56", linewidth = 0.8) +
  scale_color_manual(values = c(HH = "#b83248", LL = "#2878ad", HL = "#e8a23d", LH = "#58b6b1")) +
  labs(title = "Moran scatterplot: CRIME dengan CRIME tetangga",
       x = "z CRIME lokasi", y = "Wz CRIME tetangga", color = "Kuadran")

5.1.1 Jalankan uji lokal

Setelah mengenali kuadran setiap wilayah, kita melakukan uji Local Moran’s I untuk menilai apakah pola tersebut signifikan secara statistik.

Pengujian dilakukan dengan mempertahankan nilai CRIME wilayah yang sedang diperiksa, lalu mengacak nilai dari wilayah lain untuk mengisi posisi tetangganya. Hasil pengamatan dibandingkan dengan hasil pengacakan tersebut. Uji dua sisi digunakan agar dapat mendeteksi pola nilai yang serupa maupun berbeda dengan tetangga.

Karena pengujian dilakukan pada banyak wilayah sekaligus, p-value juga disesuaikan menggunakan metode Benjamini–Hochberg (BH). Penyesuaian ini membantu mengurangi temuan signifikan yang muncul secara kebetulan akibat banyaknya pengujian.

lisa <- localmoran_perm(
  x, lw_W, nsim = B, alternative = "two.sided",
  mlvar = FALSE, iseed = 2026, no_repeat_in_row = TRUE
)
colnames(lisa)
#> [1] "Ii"                 "E.Ii"               "Var.Ii"             "Z.Ii"              
#> [5] "Pr(z != E(Ii))"     "Pr(z != E(Ii)) Sim" "Pr(folded) Sim"     "Skewness"          
#> [9] "Kurtosis"
p_uni <- lisa[, "Pr(z != E(Ii)) Sim"]
p_uni_BH <- p.adjust(p_uni, method = "BH")
stopifnot(max(abs(lisa[, "Ii"] - zx * wzx)) < 1e-8)

hasil_uni <- data.frame(
  id = col$id, CRIME = x, lag_CRIME = lag_x,
  Ii = lisa[, "Ii"], kuadran = q_uni,
  p_permutasi = p_uni, p_BH = p_uni_BH
)
knitr::kable(head(hasil_uni[order(hasil_uni$p_permutasi), ], 8), digits = 4,
             caption = "Delapan p-value lokal terkecil; ini bukan peringkat risiko wilayah.")

Delapan p-value lokal terkecil; ini bukan peringkat risiko wilayah.

	id	CRIME	lag_CRIME	Ii	kuadran	p_permutasi	p_BH
24	24	38.2979	54.9056	0.2239	HH	0.0006	0.0229
16	16	54.8387	52.3655	1.2135	HH	0.0010	0.0229
29	29	60.7504	51.7907	1.5249	HH	0.0014	0.0229
25	25	61.2992	49.9713	1.3874	HH	0.0032	0.0392
32	32	19.1456	13.8547	1.2146	LL	0.0050	0.0490
28	28	56.9198	47.4242	0.9570	HH	0.0088	0.0719
37	37	42.4451	50.9942	0.4146	HH	0.0114	0.0798
36	36	14.3056	18.7549	1.2179	LL	0.0136	0.0833

Tabel menampilkan delapan wilayah dengan p-value permutasi terkecil. Kolom Ii menunjukkan nilai Local Moran, kuadran menunjukkan jenis pola, sedangkan p_permutasi dan p_BH menunjukkan p-value sebelum dan sesudah penyesuaian BH.

Pada taraf signifikansi 5%, lima wilayah dalam tabel tetap signifikan setelah penyesuaian BH:

• ID 24, 16, 29, dan 25 memiliki pola HH. Angka kriminalitas wilayah tersebut berada di atas rata-rata, begitu pula rata-rata angka kriminalitas tetangganya.
• ID 32 memiliki pola LL. Angka kriminalitas wilayah tersebut berada di bawah rata-rata, begitu pula rata-rata angka kriminalitas tetangganya.

Sementara itu, ID 28, 37, dan 36 signifikan berdasarkan p-value awal, tetapi tidak lagi signifikan setelah penyesuaian BH karena p_BH > 0,05. Artinya, bukti pola lokal pada ketiga wilayah tersebut belum cukup setelah memperhitungkan banyaknya pengujian.

5.1.2 Banyak uji memerlukan perhatian

Uji Local Moran dilakukan pada 49 wilayah. Ketika banyak pengujian dilakukan sekaligus, sebagian hasil dapat terlihat signifikan hanya karena kebetulan. Oleh karena itu, p-value disesuaikan menggunakan metode Benjamini–Hochberg (BH) untuk membantu mengendalikan temuan tersebut.

Kita kemudian membandingkan dua peta: peta berdasarkan p-value awal dan peta berdasarkan p-value setelah penyesuaian BH. Wilayah diberi warna sesuai kuadrannya jika p-value yang digunakan lebih kecil dari 0,05. Wilayah lainnya ditampilkan dengan warna abu-abu. kontrol yang lebih konservatif terhadap ketergantungan umum.

warna_pola <- c(
  HH = "#b83248", LL = "#2878ad", HL = "#e8a23d", LH = "#58b6b1",
  "Belum signifikan" = "#e1e6e9", "Pada batas" = "#adb8bf"
)
label_signifikan <- function(q, p) {
  factor(ifelse(p < alpha, q, "Belum signifikan"), levels = names(warna_pola))
}
peta_pola <- function(q, p, judul, subjudul) {
  tmp <- col
  tmp$pola <- label_signifikan(q, p)
  ggplot(tmp) +
    geom_sf(aes(fill = pola), color = "white", linewidth = 0.3) +
    scale_fill_manual(values = warna_pola, drop = FALSE, name = "Pola") +
    labs(title = judul, subtitle = subjudul) +
    coord_sf(datum = NA) +
  theme_void() + theme(legend.position = "bottom", plot.title = element_text(size = 14))
}

peta_pola(q_uni, p_uni, "Local Moran: CRIME", "Eksplorasi | p permutasi mentah")

peta_pola(q_uni, p_uni_BH, "Local Moran: CRIME", "Pembanding | p disesuaikan BH")

Pada peta pertama, terdapat 15 wilayah dengan pola lokal yang signifikan berdasarkan p-value awal. Pola HH berwarna merah terutama terlihat di bagian tengah, sedangkan pola LL berwarna biru terlihat pada beberapa wilayah di bagian tepi.

Setelah penyesuaian BH, terdapat 5 wilayah yang tetap signifikan, yaitu 4 wilayah berpola HH dan 1 wilayah berpola LL. Artinya, bukti pengelompokan pada kelima wilayah tersebut tetap cukup kuat setelah memperhitungkan banyaknya pengujian.

Wilayah yang berubah menjadi abu-abu tidak berarti kehilangan pola atau mengalami perubahan kriminalitas. Yang berubah adalah keputusan statistiknya: bukti yang tersedia belum cukup untuk menyatakan pola lokalnya signifikan setelah penyesuaian BH.

5.2 Getis–Ord Gi: konsentrasi tinggi dan rendah secara lokal

Getis–Ord Gi digunakan untuk mencari wilayah yang lingkungan tetangganya memiliki konsentrasi nilai tinggi (hotspot) atau nilai rendah (coldspot). Pada contoh ini, kita memeriksa konsentrasi angka kriminalitas di sekitar setiap wilayah.

Fungsi localG() menghasilkan nilai Z. Z positif menunjukkan kecenderungan konsentrasi nilai tinggi, sedangkan Z negatif menunjukkan kecenderungan konsentrasi nilai rendah. Namun, penetapan hotspot atau coldspot juga harus mempertimbangkan p-value.

Pengujian menggunakan pendekatan normal dua sisi, kemudian p-value disesuaikan dengan metode BH karena banyak wilayah diuji sekaligus.

Gi_z <- as.numeric(localG(x, lw_B))
Gi_p <- 2 * pnorm(-abs(Gi_z))   # pendekatan normal, dua sisi
Gi_p_BH <- p.adjust(Gi_p, method = "BH")

col$Gi_pola <- factor(
  ifelse(Gi_p_BH >= alpha, "Belum signifikan",
         ifelse(Gi_z > 0, "Hotspot", "Coldspot")),
  levels = c("Hotspot", "Coldspot", "Belum signifikan")
)
ggplot(col) +
  geom_sf(aes(fill = Gi_pola), color = "white", linewidth = 0.3) +
  scale_fill_manual(values = c(Hotspot = "#b83248", Coldspot = "#2878ad",
                              "Belum signifikan" = "#e1e6e9"), drop = FALSE) +
  labs(title = "Getis–Ord Gi: konsentrasi CRIME di tetangga",
       subtitle = "Uji normal dua sisi; p disesuaikan BH", fill = "Pola") +
  coord_sf(datum = NA) +
  theme_void() + theme(legend.position = "bottom")

knitr::kable(head(data.frame(id = col$id, Z_Gi = Gi_z, p = Gi_p, p_BH = Gi_p_BH)), digits = 4)

id	Z_Gi	p	p_BH
1	-0.9378	0.3483	0.5147
2	-0.9857	0.3243	0.4966
3	-0.7184	0.4725	0.6093
4	-0.0667	0.9469	0.9522
5	1.0471	0.2950	0.4819
6	0.4547	0.6493	0.8158

Pada taraf signifikansi 5%, terdapat 4 wilayah yang teridentifikasi sebagai hotspot setelah penyesuaian BH. Wilayah tersebut ditampilkan dengan warna merah di bagian tengah peta. Artinya, lingkungan tetangga keempat wilayah itu memiliki konsentrasi angka kriminalitas tinggi yang signifikan.

Tidak terdapat coldspot yang signifikan setelah penyesuaian BH. Wilayah berwarna abu-abu menunjukkan bahwa bukti konsentrasi nilai tinggi maupun rendah belum cukup, bukan berarti wilayah tersebut pasti memiliki kriminalitas rendah.

Tabel hanya menampilkan enam wilayah pertama. Seluruhnya memiliki p_BH > 0,05, sehingga belum dapat dinyatakan sebagai hotspot atau coldspot meskipun nilai Z-nya positif atau negatif.

Perbedaan Gi dan Gi* terletak pada nilai yang disertakan dalam perhitungan:

• Gi hanya menggunakan nilai wilayah-wilayah tetangga.
• Gi* menggunakan nilai wilayah-wilayah tetangga beserta nilai wilayah yang sedang diperiksa.

Fungsi include.self(nb) menambahkan wilayah itu sendiri ke dalam daftar tetangganya untuk menghitung Gi. Karena nilai yang disertakan berbeda, hasil Z pada kedua metode juga dapat berbeda. Untuk menentukan hotspot atau coldspot berdasarkan Gi, p-value dan penyesuaian BH perlu dihitung kembali menggunakan hasil Gi*.

lw_star <- nb2listw(include.self(nb), style = "B")
Gi_star_z <- as.numeric(localG(x, lw_star))
head(data.frame(id = col$id, Z_Gi = Gi_z, Z_Gi_star = Gi_star_z))
#>   id        Z_Gi   Z_Gi_star
#> 1  1 -0.93780765 -1.43277965
#> 2  2 -0.98565912 -1.34000810
#> 3  3 -0.71841891 -0.77199868
#> 4  4 -0.06666232 -0.13173336
#> 5  5  1.04713603  1.28536280
#> 6  6  0.45474575  0.06864804

Hasil kedua versi dapat berbeda karena definisi lingkungannya berbeda. Uji normal di contoh Gi juga berbeda dari inferensi permutasi Local Moran; alternatif permutasi Gi tersedia pada localG_perm(). Rujukan Gi dan Gi*.

6 Jembatan ke dua variabel

Pada analisis sebelumnya, kita membandingkan angka kriminalitas suatu wilayah dengan angka kriminalitas tetangganya. Analisis tersebut disebut univariat karena menggunakan satu variabel, yaitu CRIME.

Selanjutnya, kita mempelajari hubungan spasial antara dua variabel, yaitu kriminalitas (CRIME) dan pendapatan rumah tangga (INC). Analisis ini disebut bivariat, dengan pertanyaan:

Apakah wilayah dengan angka kriminalitas tinggi cenderung bertetangga dengan wilayah yang memiliki pendapatan rumah tangga rendah?

Dalam Moran bivariat, kita membandingkan CRIME suatu wilayah dengan rata-rata INC tetangganya. Sementara itu, korelasi Pearson membandingkan CRIME dan INC dalam wilayah yang sama.

Untuk menjaga kejelasan analisis, kita menetapkan X = CRIME dan Y = INC. Perbedaan pasangan yang dibandingkan pada setiap metode dirangkum dalam tabel berikut.

Analisis	Pasangan yang dibandingkan	Kalimat pertanyaan
Pearson	\(X_i\) dengan \(Y_i\)	Apakah dua variabel berkaitan dalam wilayah yang sama?
Moran univariat	\(X_i\) dengan \((WX)_i\)	Apakah nilai X mirip dengan X tetangga?
Moran bivariat	\(X_i\) dengan \((WY)_i\)	Apakah X berkaitan dengan Y tetangga?
Lee’s L	Pola lag X dengan pola lag Y	Apakah kedua variabel memiliki pola lingkungan yang sejalan?

6.1 Pearson sebagai pembanding deskriptif

Sebelum menganalisis hubungan spasial dua variabel, kita terlebih dahulu melihat hubungan antara pendapatan rumah tangga (INC) dan angka kriminalitas (CRIME) dalam wilayah yang sama menggunakan korelasi Pearson.

Koefisien Pearson menunjukkan arah dan kekuatan hubungan linear antara kedua variabel. Nilai positif menunjukkan hubungan searah, sedangkan nilai negatif menunjukkan hubungan berlawanan arah. Analisis ini belum memperhitungkan hubungan ketetanggaan antarwilayah.

r_xy <- cor(x, y)
r_xy
#> [1] -0.6955898

ggplot(col, aes(INC, CRIME)) +
  geom_point(color = "#25858e", size = 2.5) +
  geom_smooth(method = "lm", formula = y ~ x, se = FALSE, color = "#183f56") +
  labs(title = "Pendapatan dan kriminalitas dalam wilayah yang sama",
       subtitle = paste("Korelasi Pearson =", fmt(r_xy)),
       x = "INC | ribuan USD", y = "CRIME | per 1.000 rumah tangga")

Nilai korelasi Pearson sebesar −0,696 menunjukkan hubungan linear negatif antara pendapatan dan kriminalitas. Artinya, wilayah dengan pendapatan rumah tangga lebih tinggi cenderung memiliki angka kriminalitas lebih rendah. Kecenderungan ini terlihat dari sebaran titik dan garis yang menurun pada grafik.

Setiap titik mewakili pasangan nilai INC dan CRIME dalam satu wilayah yang sama. Hasil ini digunakan sebagai gambaran awal, belum menunjukkan hubungan dengan wilayah tetangga dan tidak membuktikan sebab-akibat. Selanjutnya, Moran bivariat digunakan untuk membandingkan kriminalitas suatu wilayah dengan pendapatan wilayah tetangganya.

7 Moran bivariat global

7.1 Apa yang dihitung?

Moran bivariat global digunakan untuk melihat hubungan antara nilai suatu variabel di sebuah wilayah dan nilai variabel lain pada tetangganya. Dalam contoh ini, kita membandingkan kriminalitas suatu wilayah (CRIME) dengan rata-rata pendapatan tetangganya (INC).

Kedua variabel terlebih dahulu diubah menjadi skor baku karena memiliki satuan berbeda. Selanjutnya, perhitungan dilakukan dalam dua arah: CRIME → INC tetangga dan INC → CRIME tetangga. Tanda panah menunjukkan urutan variabel yang dibandingkan, bukan hubungan sebab-akibat.

Dengan W yang jumlah setiap barisnya satu:

\[ I_{X\to Y} = \frac{\sum_i z_{Xi}(Wz_Y)_i}{\sum_i z_{Xi}^2}. \]

# Fungsi ringkas buatan modul; khusus bobot W tanpa wilayah terisolasi.
I_bivariat <- function(z_asal, z_tetangga, listw) {
  sum(z_asal * lag.listw(listw, z_tetangga)) / sum(z_asal^2)
}
I_xy <- I_bivariat(zx, zy, lw_W)
I_yx <- I_bivariat(zy, zx, lw_W)

knitr::kable(data.frame(
  ukuran = c("Pearson: CRIME dan INC satu wilayah",
             "Moran bivariat: CRIME -> INC tetangga",
             "Moran bivariat: INC -> CRIME tetangga"),
  nilai = c(r_xy, I_xy, I_yx)
), digits = 4)

ukuran	nilai
Pearson: CRIME dan INC satu wilayah	-0.6956
Moran bivariat: CRIME -> INC tetangga	-0.4286
Moran bivariat: INC -> CRIME tetangga	-0.4516

Hasil perhitungan menunjukkan:

• Pearson = −0,6956: wilayah dengan pendapatan lebih tinggi cenderung memiliki kriminalitas lebih rendah dalam wilayah yang sama.
• Moran CRIME → INC tetangga = −0,4286: wilayah dengan kriminalitas tinggi cenderung memiliki tetangga dengan rata-rata pendapatan rendah, dan sebaliknya.
• Moran INC → CRIME tetangga = −0,4516: wilayah dengan pendapatan tinggi cenderung memiliki tetangga dengan rata-rata kriminalitas rendah, dan sebaliknya.

Kedua nilai Moran bivariat berbeda karena pertukaran variabel mengubah pasangan yang dibandingkan, sementara bobot dari satu wilayah ke tetangganya belum tentu sama dengan bobot arah sebaliknya.

Nilai negatif menunjukkan kecenderungan hubungan yang berlawanan arah. Namun, signifikansinya masih perlu diperiksa melalui uji permutasi. Nilai indeks tersebut juga tidak dapat diartikan sebagai persentase pengaruh. Rujukan definisi Moran bivariat.

7.2 Uji permutasi: nyatakan apa yang diacak

Untuk menilai apakah hubungan spasial yang ditemukan signifikan, kita melakukan uji permutasi. Nilai pendapatan (INC) diacak antarwilayah, sedangkan nilai kriminalitas (CRIME) dan bobot ketetanggaan tetap.

Setiap pengacakan menghasilkan nilai Moran bivariat baru. Proses ini diulang 9.999 kali, kemudian hasilnya dibandingkan dengan Moran bivariat dari data asli.

Hipotesis yang digunakan adalah:

• H₀: penempatan nilai INC antarwilayah bersifat acak, dengan CRIME dan bobot spasial tetap.
• H₁: hubungan yang diamati lebih kuat, baik ke arah positif maupun negatif, daripada yang diharapkan dari pengacakan tersebut.

Karena kedua arah diperiksa, pengujian menggunakan nilai absolut indeks untuk menghitung p-value dua sisi.

set.seed(2027)
I_xy_sim <- replicate(B, {
  zy_acak <- sample(zy, replace = FALSE)
  I_bivariat(zx, zy_acak, lw_W)
})
p_xy <- (1 + sum(abs(I_xy_sim) >= abs(I_xy))) / (B + 1)
data.frame(I_CRIME_INC = I_xy, p_dua_sisi = p_xy, permutasi = B)
#>   I_CRIME_INC p_dua_sisi permutasi
#> 1  -0.4286322      1e-04      9999

Moran bivariat CRIME terhadap INC tetangga sebesar −0,4286, dengan p-value 0,0001. Karena p-value lebih kecil dari 0,05, H₀ ditolak.

Hasil ini mendukung hubungan spasial yang berlawanan arah: wilayah dengan kriminalitas tinggi cenderung memiliki tetangga dengan rata-rata pendapatan rendah, sedangkan wilayah dengan kriminalitas rendah cenderung memiliki tetangga dengan rata-rata pendapatan tinggi.

Kesimpulan tersebut berlaku terhadap pengacakan INC yang digunakan. Pengujian ini belum mengendalikan pola spasial pendapatan maupun hubungan pendapatan–kriminalitas dalam wilayah yang sama, dan tidak membuktikan sebab-akibat.

7.3 Bivariate Moran scatterplot

Bivariate Moran scatterplot memperlihatkan hubungan antara kriminalitas suatu wilayah dan rata-rata pendapatan tetangganya. Setiap titik mewakili satu wilayah.

Sumbu horizontal menunjukkan skor baku CRIME wilayah tersebut, sedangkan sumbu vertikal menunjukkan rata-rata skor baku INC tetangganya. Garis nol membagi grafik menjadi empat kuadran. Huruf pertama menunjukkan kondisi kriminalitas wilayah, sedangkan huruf kedua menunjukkan kondisi pendapatan tetangganya.

q_bv <- kuadran(zx, wzy)
df_bv <- data.frame(id = col$id, zx = zx, wzy = wzy, kuadran = q_bv)
ggplot(df_bv, aes(zx, wzy, color = kuadran)) +
  geom_hline(yintercept = 0, color = "grey65") +
  geom_vline(xintercept = 0, color = "grey65") +
  geom_point(size = 2.6) +
  geom_abline(intercept = 0, slope = I_xy, color = "#183f56", linewidth = 0.9) +
  scale_color_manual(values = warna_pola, drop = FALSE) +
  labs(title = "Moran bivariat: kriminalitas wilayah dan pendapatan tetangga",
       subtitle = paste("Kemiringan garis melalui nol =", fmt(I_xy)),
       x = "z CRIME lokasi", y = "Wz INC tetangga", color = "Kuadran")

Kuadran	Interpretasi spesifik kasus
HH	CRIME wilayah tinggi, rata-rata INC tetangga tinggi
LL	CRIME wilayah rendah, rata-rata INC tetangga rendah
HL	CRIME wilayah tinggi, rata-rata INC tetangga rendah
LH	CRIME wilayah rendah, rata-rata INC tetangga tinggi

Sebagian besar titik berada pada kuadran HL dan LH:

• HL (kanan bawah): kriminalitas wilayah di atas rata-rata, sedangkan rata-rata pendapatan tetangganya di bawah rata-rata.
• LH (kiri atas): kriminalitas wilayah di bawah rata-rata, sedangkan rata-rata pendapatan tetangganya di atas rata-rata.

Beberapa wilayah berada pada kuadran HH, yaitu kriminalitas wilayah dan rata-rata pendapatan tetangganya sama-sama tinggi, atau LL, yaitu keduanya sama-sama rendah. Tinggi dan rendah di sini dibandingkan dengan rata-rata masing-masing variabel di seluruh wilayah.

Garis yang menurun memiliki kemiringan sekitar −0,429, sesuai dengan Moran bivariat yang telah dihitung. Pola ini menunjukkan kecenderungan hubungan berlawanan arah antara kriminalitas wilayah dan pendapatan tetangganya.

Warna kuadran menunjukkan jenis pola, belum menunjukkan signifikansi setiap wilayah. Untuk menilainya, analisis dilanjutkan dengan Local Moran bivariat.

8 Local Moran bivariat

8.1 Dari pola keseluruhan ke wilayah tertentu

Local Moran bivariat digunakan untuk mengetahui wilayah mana yang menunjukkan hubungan antara kriminalitas wilayah tersebut dan rata-rata pendapatan tetangganya.

Pengujian dilakukan pada setiap wilayah dengan mempertahankan nilai CRIME wilayah tersebut, lalu membandingkan hasilnya dengan pengacakan nilai INC dari wilayah lain. Karena banyak wilayah diuji sekaligus, p-value disesuaikan menggunakan metode Benjamini–Hochberg (BH).

local_bv <- localmoran_bv(
  x, y, lw_W, nsim = B, scale = TRUE,
  alternative = "two.sided", iseed = 2028,
  no_repeat_in_row = TRUE
)
colnames(local_bv)
#> [1] "Ibvi"                 "E.Ibvi"               "Var.Ibvi"            
#> [4] "Z.Ibvi"               "Pr(z != E(Ibvi))"     "Pr(z != E(Ibvi)) Sim"
#> [7] "Pr(folded) Sim"

p_bv <- local_bv[, "Pr(z != E(Ibvi)) Sim"]
p_bv_BH <- p.adjust(p_bv, method = "BH")
stopifnot(max(abs(local_bv[, "Ibvi"] - zx * wzy)) < 1e-8)

hasil_bv <- data.frame(
  id = col$id, CRIME = x, INC = y,
  lag_INC = lag.listw(lw_W, y),
  Ibvi = local_bv[, "Ibvi"], kuadran = q_bv,
  p_permutasi = p_bv, p_BH = p_bv_BH
)
knitr::kable(head(hasil_bv[order(hasil_bv$p_permutasi), ], 10), digits = 4,
             caption = "Sepuluh p-value terkecil untuk CRIME lokasi dengan INC tetangga.")

Sepuluh p-value terkecil untuk CRIME lokasi dengan INC tetangga.

	id	CRIME	INC	lag_INC	Ibvi	kuadran	p_permutasi	p_BH
32	32	19.1456	18.942	26.0663	-1.9581	LH	0.0004	0.0196
8	8	38.4259	11.337	8.5385	-0.2016	HL	0.0008	0.0196
12	12	56.7057	10.048	9.0632	-1.2010	HL	0.0026	0.0425
40	40	16.2413	29.833	20.7966	-1.2710	LH	0.0060	0.0735
47	47	27.8229	18.950	26.0200	-0.8915	LH	0.0126	0.1225
16	16	54.8387	7.625	10.5606	-0.7878	HL	0.0150	0.1225
41	41	18.9051	22.207	22.5750	-1.3941	LH	0.0202	0.1263
28	28	56.9198	7.856	10.9571	-0.7804	HL	0.0224	0.1263
15	15	48.5855	9.873	10.0793	-0.6057	HL	0.0232	0.1263
14	14	57.0661	9.963	10.0772	-0.9880	HL	0.0280	0.1336

Pada taraf signifikansi 5%, wilayah ID 32, 8, dan 12 tetap signifikan setelah penyesuaian BH, karena nilai p_BH masing-masing kurang dari 0,05.

• Wilayah ID 32 berpola LH: kriminalitasnya di bawah rata-rata, sedangkan rata-rata pendapatan tetangganya di atas rata-rata.
• Wilayah ID 8 dan 12 berpola HL: kriminalitasnya di atas rata-rata, sedangkan rata-rata pendapatan tetangganya di bawah rata-rata.

Tinggi dan rendah ditentukan berdasarkan rata-rata masing-masing variabel di seluruh wilayah. Untuk membaca kondisi pendapatan tetangga, gunakan kolom lag_INC; kolom INC menunjukkan pendapatan wilayah itu sendiri.

Wilayah lain dalam tabel memiliki p_BH lebih dari 0,05 sehingga bukti pola lokalnya belum cukup setelah penyesuaian. Hasil ini menunjukkan keterkaitan antarwilayah, tetapi tidak membuktikan bahwa pendapatan tetangga menyebabkan tingkat kriminalitas suatu wilayah. Dokumentasi dan contoh fungsi.

8.2 Peta dan satu contoh interpretasi

Peta berikut menunjukkan lokasi pola hubungan antara kriminalitas suatu wilayah dan rata-rata pendapatan tetangganya. Peta pertama menggunakan p-value sebelum penyesuaian, sedangkan peta kedua menggunakan p-value setelah penyesuaian BH untuk memperhitungkan banyaknya pengujian. Perbandingan keduanya membantu kita melihat pola yang tetap signifikan setelah penyesuaian.

peta_pola(q_bv, p_bv, "Moran lokal bivariat: CRIME -> INC tetangga",
          "Eksplorasi | p permutasi mentah; huruf pertama CRIME, huruf kedua lag INC")

peta_pola(q_bv, p_bv_BH, "Moran lokal bivariat: CRIME -> INC tetangga",
          "Pembanding | p disesuaikan BH")

knitr::kable(as.data.frame(table(
  Kuadran = q_bv, Lolos_BH = ifelse(p_bv_BH < alpha, "Ya", "Belum")
)), caption = "Jumlah wilayah menurut kuadran dan hasil penyesuaian BH.")

Jumlah wilayah menurut kuadran dan hasil penyesuaian BH.

Kuadran	Lolos_BH	Freq
HH	Belum	4
HL	Belum	18
LH	Belum	18
LL	Belum	6
HH	Ya	0
HL	Ya	2
LH	Ya	1
LL	Ya	0

Pada taraf signifikansi 5%, terdapat 14 wilayah yang signifikan sebelum penyesuaian BH. Setelah penyesuaian, 3 wilayah tetap signifikan, yaitu:

• Dua wilayah berpola HL (oranye): kriminalitas wilayah di atas rata-rata, sedangkan rata-rata pendapatan tetangganya di bawah rata-rata.
• Satu wilayah berpola LH (hijau kebiruan): kriminalitas wilayah di bawah rata-rata, sedangkan rata-rata pendapatan tetangganya di atas rata-rata.

Tidak ada wilayah berpola HH atau LL yang tetap signifikan setelah penyesuaian BH. Wilayah berwarna abu-abu menunjukkan bahwa bukti statistiknya belum cukup untuk menyatakan adanya pola lokal yang signifikan.

Berkurangnya wilayah berwarna terjadi karena penyesuaian dalam menilai signifikansi; nilai kriminalitas dan pendapatannya tetap sama.

# Ambil satu wilayah untuk latihan membaca, bukan untuk mengubah aturan pengujian.
contoh <- hasil_bv[which.min(hasil_bv$p_permutasi), ]
knitr::kable(contoh, digits = 4, row.names = FALSE)

id	CRIME	INC	lag_INC	Ibvi	kuadran	p_permutasi	p_BH
32	19.1456	18.942	26.0663	-1.9581	LH	4e-04	0.0196

Wilayah ID 32 memiliki CRIME sekitar 19,146 per 1.000 rumah tangga, sedangkan rata-rata INC tetangganya sekitar 26,066 ribu USD. Wilayah ini termasuk pola LH karena kriminalitasnya di bawah rata-rata seluruh wilayah, sementara rata-rata pendapatan tetangganya di atas rata-rata.

Nilai Local Moran bivariat sebesar −1,9581, dengan p-value permutasi 0,0004 dan p-value setelah penyesuaian BH 0,0196. Karena p-value BH kurang dari 0,05, pola LH tersebut tetap signifikan.

Pada peta lokasi, warna oranye menandai wilayah ID 32, sedangkan hijau muda menandai tetangganya, yaitu ID 20, 23, 40, dan 41. Nilai INC keempat tetangga inilah yang digunakan untuk menghitung lag_INC. Warna pada peta lokasi menunjukkan peran wilayah, bukan kategori HL atau LH.

lokator <- col
lokator$peran <- "Wilayah lain"
lokator$peran[nb[[contoh$id]]] <- "Tetangga contoh"
lokator$peran[contoh$id] <- "Wilayah contoh"
label_id <- data.frame(id = col$id, X = titik[, 1], Y = titik[, 2])
ggplot(lokator) +
  geom_sf(aes(fill = peran), color = "white", linewidth = 0.4) +
  geom_text(data = label_id, aes(X, Y, label = id), size = 2.5, color = "#183446") +
  scale_fill_manual(values = c("Wilayah contoh" = "#e8a23d",
                              "Tetangga contoh" = "#9edbd6",
                              "Wilayah lain" = "#e8edf0")) +
  labs(title = paste("Lokasi wilayah contoh: ID", contoh$id), fill = "Peran") +
  coord_sf(datum = NA) +
  theme_void() + theme(legend.position = "bottom")

Temukan wilayah contoh dan tetangganya; label menunjukkan ID modul.

9 Pengayaan: dua sudut pandang bivariat lain

Bagian ini menghubungkan praktikum dengan materi Geary bivariat dan Lee’s L dalam slide. Fokuskan pada pertanyaan dan interpretasinya; penurunan rumus dapat dibaca setelah kelas.

9.1 Geary bivariat: selisih lintas variabel

Geary bivariat dalam modul ini digunakan untuk melihat seberapa berbeda nilai kriminalitas suatu wilayah dan pendapatan wilayah tetangganya. Karena satuannya berbeda, CRIME dan INC terlebih dahulu diubah menjadi skor baku agar dapat dibandingkan.

Perhitungan menggunakan kuadrat selisih kedua skor tersebut pada pasangan wilayah yang bertetangga. Nilai yang lebih kecil menunjukkan skor yang lebih mirip. Hasilnya kemudian dibandingkan dengan hasil pengacakan lokasi nilai INC, sementara CRIME dan hubungan ketetanggaan tetap.

Dengan skor baku sampel, rumus pada bahan kuliah menyederhana menjadi:

\[ C_{XY}=\frac{1}{2S_0}\sum_i\sum_j w_{ij}(z_{Xi}-z_{Yj})^2. \]

C_bivariat <- function(za, zb, W) {
  sum(W * outer(za, zb, "-")^2) / (2 * sum(W))
}
C_xy <- C_bivariat(zx, zy, W)

set.seed(2029)
C_xy_sim <- replicate(B, C_bivariat(zx, sample(zy), W))
p_C_bawah <- (1 + sum(C_xy_sim <= C_xy)) / (B + 1)
p_C_atas <- (1 + sum(C_xy_sim >= C_xy)) / (B + 1)
p_C_xy <- min(1, 2 * min(p_C_bawah, p_C_atas))

knitr::kable(data.frame(
  Cxy = C_xy, rerata_permutasi = mean(C_xy_sim),
  harapan_teoretis = (n - 1) / n, p_dua_sisi = p_C_xy
), digits = 4)

Cxy	rerata_permutasi	harapan_teoretis	p_dua_sisi
1.5029	0.9801	0.9796	2e-04

Nilai Geary bivariat yang diperoleh adalah 1,5029, lebih besar daripada rata-rata hasil pengacakan sebesar 0,9801. Nilai p-value dua sisi sebesar 0,0002, lebih kecil dari 0,05, sehingga H0 ditolak.

Artinya, perbedaan skor baku kriminalitas wilayah dan pendapatan tetangganya lebih besar daripada yang biasanya muncul pada susunan acak. Hasil ini sejalan dengan temuan Moran bivariat sebelumnya yang menunjukkan kecenderungan hubungan berlawanan arah.

Pada perhitungan ini, nilai pembanding saat pengacakan sekitar 0,98. Oleh karena itu, interpretasi didasarkan pada hasil pengacakan dan p-value, bukan langsung memakai patokan angka 1 dari Geary univariat.

9.2 Lee’s L: pola lingkungan kedua variabel

Lee’s L digunakan untuk melihat hubungan antara rata-rata kriminalitas tetangga dan rata-rata pendapatan tetangga di sekitar setiap wilayah. Dengan demikian, kedua variabel sama-sama diringkas berdasarkan kondisi wilayah tetangganya.

Pengujian dilakukan dengan mengacak lokasi pasangan CRIME–INC secara bersama-sama, sementara hubungan ketetanggaan tetap. Cara ini mempertahankan pasangan nilai kriminalitas dan pendapatan dari wilayah yang sama, tetapi mengubah susunannya pada peta. Untuk W pada contoh ini:

\[ L=\frac{\sum_i (Wz_X)_i(Wz_Y)_i} {\sqrt{\sum_i z_{Xi}^2\sum_i z_{Yi}^2}}. \]

L_manual <- sum(wzx * wzy) / sqrt(sum(zx^2) * sum(zy^2))
L_paket <- lee(x, y, lw_W, n = n)$L
stopifnot(isTRUE(all.equal(as.numeric(L_paket), L_manual)))

set.seed(2030)
L_uji <- lee.mc(x, y, lw_W, nsim = B, alternative = "two.sided")
L_uji
#> 
#>  Monte-Carlo simulation of Lee's L
#> 
#> data:  x ,  y 
#> weights: lw_W  
#> number of simulations + 1: 10000 
#> 
#> statistic = -0.47031, observed rank = 1, p-value = 2e-04
#> alternative hypothesis: two.sided

Nilai Lee’s L sebesar −0,4703 menunjukkan bahwa pola lingkungan kedua variabel cenderung berlawanan arah. Lingkungan dengan rata-rata kriminalitas tetangga yang tinggi cenderung memiliki rata-rata pendapatan tetangga yang rendah, dan sebaliknya.

Hasil uji menghasilkan p-value 0,0002, lebih kecil dari 0,05, sehingga H0 ditolak. Artinya, pola yang diamati berbeda secara signifikan dari pola yang biasanya muncul ketika lokasi pasangan CRIME–INC diacak.

Nilai −0,4703 menunjukkan arah dan besarnya statistik hubungan spasial, bukan persentase pengaruh. Hasil ini juga tidak membuktikan hubungan sebab-akibat antara pendapatan dan kriminalitas.

9.3 Latihan penerapan

Latihan ini menggunakan data Columbus dan bobot ketetanggaan Queen yang telah dibuat. Gunakan taraf signifikansi 5% dan jelaskan setiap hasil dengan bahasa sendiri.

9.3.1 1. Memahami nilai wilayah dan nilai tetangganya

Pilih wilayah ID 16, kemudian:

1. Tampilkan ID seluruh tetangganya.
2. Tampilkan nilai CRIME wilayah tersebut dan nilai CRIME setiap tetangganya.
3. Hitung rata-rata CRIME tetangganya secara manual, lalu bandingkan dengan hasil lag.listw().
4. Jelaskan apakah kriminalitas wilayah ID 16 lebih tinggi atau lebih rendah daripada rata-rata kriminalitas tetangganya.

Hasil yang dikumpulkan: tabel tetangga, hasil perhitungan, dan interpretasi 2–3 kalimat.

9.3.2 2. Membaca pola kriminalitas secara keseluruhan

Gunakan hasil Moran’s I dan Geary’s C untuk variabel CRIME.

1. Buat tabel berisi nama metode, nilai statistik, dan p-value. Untuk Moran’s I, gunakan p-value permutasi.
2. Tentukan keputusan pengujian masing-masing metode pada taraf 5%.
3. Jelaskan apakah kedua metode memberikan kesimpulan yang sejalan.
4. Jelaskan mengapa Moran’s I yang positif dan Geary’s C yang kurang dari 1 dapat menunjukkan pola yang sama.

Hasil yang dikumpulkan: satu tabel perbandingan dan satu paragraf penjelasan.

9.3.3 3. Menilai pengaruh penyesuaian BH

Gunakan hasil Local Moran univariat untuk CRIME.

1. Hitung jumlah wilayah yang signifikan sebelum dan setelah penyesuaian BH.
2. Sebutkan ID dan kuadran wilayah yang tetap signifikan setelah penyesuaian.
3. Pilih satu wilayah yang signifikan sebelum penyesuaian, tetapi tidak lagi signifikan setelah BH. Tuliskan kedua p-value-nya.
4. Jelaskan mengapa perubahan status tersebut tidak berarti nilai kriminalitas wilayahnya berubah.

Hasil yang dikumpulkan: tabel ringkas dan interpretasi 3–4 kalimat.

9.3.4 4. Membandingkan dua pola lokal bivariat

Gunakan hasil Local Moran bivariat CRIME terhadap INC tetangga untuk wilayah ID 8 dan ID 32.

1. Buat tabel berisi ID, CRIME, INC, lag_INC, kuadran, dan p-value BH.
2. Jelaskan arti kuadran masing-masing wilayah.
3. Tentukan apakah pola keduanya signifikan setelah penyesuaian BH.
4. Jelaskan mengapa pembacaan huruf kedua pada kuadran menggunakan lag_INC.
5. Tuliskan satu kesimpulan tentang hubungan kriminalitas wilayah dan pendapatan tetangganya tanpa menyatakan hubungan sebab-akibat.

Hasil yang dikumpulkan: satu tabel dan satu paragraf perbandingan.

9.3.5 5. Memeriksa apakah pilihan tetangga mengubah kesimpulan

Bentuk ketetanggaan Rook dengan mengganti queen = TRUE menjadi queen = FALSE pada poly2nb(). Buat kembali bobot dengan style = "W".

1. Periksa apakah ada wilayah yang tidak memiliki tetangga.
2. Bandingkan jumlah pasangan tetangga unik antara Queen dan Rook.
3. Hitung kembali Moran’s I global untuk CRIME menggunakan bobot Rook.
4. Jalankan uji permutasi dengan jumlah pengacakan dan arah pengujian yang sama seperti pada bobot Queen.
5. Bandingkan nilai indeks, p-value, dan kesimpulan kedua bobot.

Hasil yang dikumpulkan: tabel perbandingan Queen–Rook dan jawaban singkat: Apakah kesimpulan pola kriminalitas tetap sama ketika definisi tetangga berubah?

Petunjuk: gunakan card(nb) untuk menghitung jumlah tetangga. Jumlah pasangan tetangga unik dapat dihitung dengan sum(card(nb)) / 2.

9.3.6 Ketentuan jawaban

Sertakan sintaks yang digunakan dan output yang relevan. Setiap interpretasi perlu menyebutkan variabel, wilayah atau cakupan analisis, serta bukti yang mendukung kesimpulan. Bedakan hasil yang menggambarkan seluruh wilayah dengan hasil untuk wilayah tertentu.