PEMILIHAN MODEL TERBAIK DALAM MENGANALISIS FAKTOR EKONOMI DAN SOSIAL TERHADAP TINGKAT KEMISKINAN DI JAWA BARAT TAHUN 2024
Dinar Olivia Rahmadyta (140610230045)
Program Studi S1 Statistika, FMIPA Universitas Padjadjaran
required <- c("knitr","kableExtra","sf","spdep","spatialreg","dplyr","ggplot2",
"png","grid","gridExtra")
to_install <- required[!(required %in% installed.packages()[,1])]
if(length(to_install)) install.packages(to_install, dependencies = TRUE)
lapply(required, library, character.only = TRUE)## [[1]]
## [1] "knitr" "stats" "graphics" "grDevices" "utils" "datasets"
## [7] "methods" "base"
##
## [[2]]
## [1] "kableExtra" "knitr" "stats" "graphics" "grDevices"
## [6] "utils" "datasets" "methods" "base"
##
## [[3]]
## [1] "sf" "kableExtra" "knitr" "stats" "graphics"
## [6] "grDevices" "utils" "datasets" "methods" "base"
##
## [[4]]
## [1] "spdep" "spData" "sf" "kableExtra" "knitr"
## [6] "stats" "graphics" "grDevices" "utils" "datasets"
## [11] "methods" "base"
##
## [[5]]
## [1] "spatialreg" "Matrix" "spdep" "spData" "sf"
## [6] "kableExtra" "knitr" "stats" "graphics" "grDevices"
## [11] "utils" "datasets" "methods" "base"
##
## [[6]]
## [1] "dplyr" "spatialreg" "Matrix" "spdep" "spData"
## [6] "sf" "kableExtra" "knitr" "stats" "graphics"
## [11] "grDevices" "utils" "datasets" "methods" "base"
##
## [[7]]
## [1] "ggplot2" "dplyr" "spatialreg" "Matrix" "spdep"
## [6] "spData" "sf" "kableExtra" "knitr" "stats"
## [11] "graphics" "grDevices" "utils" "datasets" "methods"
## [16] "base"
##
## [[8]]
## [1] "png" "ggplot2" "dplyr" "spatialreg" "Matrix"
## [6] "spdep" "spData" "sf" "kableExtra" "knitr"
## [11] "stats" "graphics" "grDevices" "utils" "datasets"
## [16] "methods" "base"
##
## [[9]]
## [1] "grid" "png" "ggplot2" "dplyr" "spatialreg"
## [6] "Matrix" "spdep" "spData" "sf" "kableExtra"
## [11] "knitr" "stats" "graphics" "grDevices" "utils"
## [16] "datasets" "methods" "base"
##
## [[10]]
## [1] "gridExtra" "grid" "png" "ggplot2" "dplyr"
## [6] "spatialreg" "Matrix" "spdep" "spData" "sf"
## [11] "kableExtra" "knitr" "stats" "graphics" "grDevices"
## [16] "utils" "datasets" "methods" "base"Laporan ini disusun untuk memenuhi penilaian Mata Kuliah Spatial Statistics Kelas B TA 2025/2026. Penulis mengucapkan terima kasih kepada semua pihak yang telah membantu. Segala saran dan kritik sangat diharapkan untuk perbaikan di masa mendatang.
Kemiskinan mencerminkan rendahnya tingkat pendapatan seseorang sehingga tidak mampu memenuhi standar hidup yang layak (Corral et al., 2020). Kondisi ini, jika berlangsung dalam jangka panjang, dapat menimbulkan berbagai dampak sosial dan ekonomi seperti meningkatnya ketimpangan sosial, terhambatnya pembangunan, munculnya konflik horizontal, serta naiknya tingkat kriminalitas (UN, 2023). Apabila tidak diatasi secara tepat, kemiskinan berpotensi menciptakan siklus antargenerasi yang sulit diputus (Rivera & Currais, 1999; Benhabib & Spiegel, 1994; Lucas, 1990).
Meskipun tingkat kemiskinan di Indonesia, khususnya di Provinsi Jawa Barat, terus mengalami penurunan dari 8,43% pada tahun 2020, menjadi 7,97% pada tahun 2021, 7,98% pada tahun 2022, kemudian 7,46% pada tahun 2023, hingga mencapai 7,08% pada tahun 2024, secara keseluruhan persentase penduduk miskin selama periode 2020–2024 tidak mengalami banyak perubahan. Hal ini menunjukkan perlunya analisis yang lebih mendalam terhadap faktor-faktor yang memengaruhi kemiskinan, agar kebijakan penanggulangannya di Jawa Barat dapat dirancang secara lebih efektif dan berkelanjutan.
Pemetaan persentase kemiskinan di kabupaten/kota Provinsi Jawa Barat menunjukkan adanya indikasi tidak independen secara spasial, di mana nilai kemiskinan di satu kab/kota memiliki keterkaitan dengan nilai observasi di kab/kota lain yang berdekatan (spatial dependence). Kondisi ini melanggar asumsi independensi pada regresi klasik (OLS), sehingga estimasi berpotensi bias dan tidak efisien. Oleh karena itu, digunakan pendekatan ekonometrika spasial untuk menangkap pengaruh spasial secara lebih akurat.
Penelitian ini bertujuan untuk menerapkan pendekatan spasial dalam menganalisis faktor-faktor yang memengaruhi tingkat kemiskinan di Provinsi Jawa Barat pada tahun 2024. Melalui analisis spasial, penelitian ini diharapkan dapat mengidentifikasi pola persebaran kemiskinan antarwilayah serta memperoleh model terbaik yang mampu menjelaskan hubungan spasial antar daerah dan variabel-variabel utama yang berpengaruh terhadap kemiskinan.
Penelitian ini dibatasi pada wilayah administrasi Provinsi Jawa Barat dengan unit analisis kabupaten/kota. Fokus penelitian diarahkan pada analisis pola spasial kemiskinan tanpa membahas secara mendalam faktor-faktor non-spasial lainnya yang bersifat kualitatif, seperti kebijakan lokal atau faktor sosial budaya.
Data yang digunakan dalam penelitian ini merupakan data tahun 2024, sehingga hasil analisis menggambarkan kondisi pada periode waktu tertentu (cross-section) dan tidak mencakup dinamika perubahan kemiskinan dari tahun ke tahun. Selain itu, penelitian ini tidak membahas aspek kemiskinan ekstrem di tingkat rumah tangga, melainkan berfokus pada analisis agregat wilayah.
Dalam analisis spasial, keterkaitan antarwilayah diukur melalui konsep spatial autocorrelation yang menunjukkan sejauh mana suatu fenomena di satu lokasi memiliki kemiripan dengan lokasi di sekitarnya. Salah satu ukuran yang paling umum digunakan adalah Moran’s I. Nilai I > 0 menunjukkan adanya autokorelasi spasial positif, yang berarti wilayah dengan nilai variabel tinggi (high-high) cenderung berdekatan dengan wilayah bernilai tinggi lainnya, dan wilayah dengan nilai rendah (low-low) juga cenderung berdekatan. Sebaliknya, I < 0 menunjukkan adanya autokorelasi spasial negatif, di mana wilayah dengan nilai tinggi dikelilingi oleh wilayah dengan nilai rendah (pola checkerboard atau high-low). Sementara itu, I ≈ 0 menunjukkan pola persebaran yang acak (random pattern).
Secara matematis, perhitungan Moran’s I melibatkan matriks bobot spasial (W) yang merepresentasikan tingkat kedekatan atau hubungan antarwilayah. Matriks ini dapat dibentuk berdasarkan beberapa kriteria kedekatan, seperti Rook contiguity (berbagi sisi batas), Queen contiguity (berbagi sisi atau sudut batas), maupun berdasarkan jarak geografis antar titik pusat wilayah. Nilai dalam matriks bobot ini menentukan seberapa besar pengaruh suatu wilayah terhadap wilayah lain di sekitarnya.
Selain Moran’s I, terdapat pula ukuran autokorelasi spasial lain seperti Geary’s C, yang lebih sensitif terhadap perbedaan lokal antarwilayah, serta Getis-Ord G* dan Local Moran’s I, yang digunakan untuk mengidentifikasi klaster atau pola lokal kemiskinan (hot spot dan cold spot). Dengan demikian, analisis autokorelasi spasial tidak hanya menggambarkan adanya keterkaitan antarwilayah, tetapi juga membantu dalam memahami struktur spasial dan pola penyebaran suatu fenomena secara geografis.
Model Ordinary Least Squares (OLS) merupakan metode dasar dalam regresi linier yang digunakan untuk mengestimasi hubungan antara variabel dependen dan independen. Persamaan umumnya ditulis sebagai:
\[ Y = X\beta + \varepsilon \]
dengan \(Y\) = vektor respon, \(X\) = matriks kovariat, \(\beta\) = parameter, dan \(\varepsilon\) = galat.
Metode OLS didasarkan pada pemenuhan asumsi klasik, di antaranya tidak adanya multikolinearitas antarvariabel independen, residual berdistribusi normal, serta homoskedastisitas (varian residual konstan). Selain itu, OLS mengasumsikan bahwa setiap observasi bersifat independen dan tidak saling memengaruhi. Dalam konteks spasial, asumsi ini sering kali tidak terpenuhi karena adanya pengaruh antarwilayah, misalnya fenomena spillover di mana kondisi kemiskinan di suatu wilayah dipengaruhi oleh wilayah sekitarnya.
Apabila asumsi independensi ini dilanggar, maka hasil estimasi OLS dapat menjadi bias dan tidak efisien, karena model tidak mampu menangkap spatial dependence yang nyata di lapangan. Oleh karena itu, untuk menganalisis fenomena yang memiliki keterkaitan geografis seperti kemiskinan antarwilayah, digunakan model ekonometrika spasial yang mampu memperhitungkan pengaruh spasial tersebut, seperti Spatial Autoregressive (SAR) dan Spatial Error Model (SEM).
Uji Lagrange Multiplier (LM) digunakan untuk mendeteksi ada tidaknya dependensi spasial dalam model regresi. Uji ini dilakukan setelah estimasi OLS guna menentukan apakah model spasial diperlukan.
Apabila hasil uji menunjukkan adanya spatial lag dependence, artinya terdapat pengaruh langsung antarwilayah dalam variabel dependen, maka model yang sesuai adalah Spatial Autoregressive Model (SAR). Sebaliknya, jika ditemukan spatial error dependence, artinya ketergantungan spasial muncul melalui error atau residu, maka digunakan Spatial Error Model (SEM).
Dalam praktiknya, terdapat dua bentuk uji LM yang sering digunakan, yaitu LM-Lag dan LM-Error. Selain itu, versi robust dari kedua uji ini (Robust LM-Lag dan Robust LM-Error) digunakan untuk memastikan hasil yang lebih akurat jika terdapat kemungkinan keduanya signifikan secara bersamaan.
Dengan demikian, uji Lagrange Multiplier berperan penting sebagai langkah diagnostik dalam pemodelan ekonometrika spasial untuk menentukan model yang paling sesuai sebelum analisis lanjut dilakukan.
Menurut Anselin (2003), Spatial Autoregressive Model (SAR) merupakan pengembangan dari model Ordinary Least Squares (OLS) dengan memasukkan unsur pengaruh spasial melalui penambahan komponen lag pada variabel dependen. Model ini digunakan ketika terdapat indikasi bahwa nilai suatu variabel di suatu wilayah dipengaruhi oleh nilai variabel yang sama di wilayah-wilayah sekitarnya. Dengan kata lain, SAR merepresentasikan adanya efek saling ketergantungan antarwilayah dalam variabel respon.
Secara matematis, model SAR ditulis sebagai:
\[ Y = \rho WY + X\beta + \varepsilon \]
dengan: - \(W\): matriks bobot spasial, - \(\rho\): koefisien autokorelasi spasial (kekuatan pengaruh tetangga), - \(X, \beta, \varepsilon\): seperti didefinisikan sebelumnya.
Asumsi utama dalam model SAR adalah bahwa residual tidak mengalami autokorelasi spasial, sebab pengaruh spasial telah sepenuhnya diakomodasi oleh komponen lag dari variabel dependen. Jika asumsi ini terpenuhi, maka model SAR dianggap tepat untuk menggambarkan struktur keterkaitan spasial yang terdapat dalam data.
Model SAR banyak digunakan dalam penelitian sosial-ekonomi, termasuk analisis kemiskinan antarwilayah, karena dapat menunjukkan adanya spillover effect, yaitu pengaruh kondisi suatu daerah terhadap daerah tetangganya. Dengan demikian, model ini lebih representatif dibandingkan model OLS dalam konteks data yang memiliki dimensi geografis.
Menurut Suryaningrat (2021), Akaike Information Criterion (AIC) merupakan ukuran statistik yang digunakan untuk memilih model terbaik di antara beberapa model yang diestimasi. AIC dikembangkan oleh Hirotugu Akaike (1974) dan didasarkan pada prinsip maximum likelihood estimation (MLE). Tujuan utama AIC adalah menyeimbangkan antara tingkat kecocokan model terhadap data (goodness of fit) dan kompleksitas model (jumlah parameter yang digunakan).
Secara umum, AIC dihitung dengan rumus:
\[ AIC = -2 \ln(L) + 2k \]
dengan \(L\) = nilai maximum likelihood dan \(k\) = jumlah parameter.
Nilai AIC yang lebih kecil menunjukkan model yang lebih baik, karena memberikan keseimbangan optimal antara tingkat ketepatan model dan kesederhanaannya. Model dengan AIC paling rendah dianggap paling efisien dalam menjelaskan data tanpa menimbulkan overfitting.
Dalam konteks analisis spasial, AIC sering digunakan untuk membandingkan model OLS, SAR, dan SEM, sehingga peneliti dapat menentukan model mana yang paling sesuai untuk menggambarkan pola keterkaitan spasial antarwilayah. Dengan demikian, AIC berperan penting sebagai alat evaluasi dan pemilihan model terbaik dalam penelitian ekonometrika spasial.
Penelitian ini menggunakan data sekunder yang diperoleh dari Badan Pusat Statistik (BPS) tahun 2024 dengan unit analisis sebanyak 27 kabupaten/kota di Provinsi Jawa Barat. Data yang digunakan bersifat cross section, sehingga menggambarkan kondisi sosial ekonomi wilayah pada satu periode waktu tertentu.
Variabel yang digunakan terdiri atas:
Variabel independen:
Tingkat Pengangguran Terbuka (TPT) (%),
Produk Domestik Regional Bruto (PDRB) per kapita (juta rupiah),
Indeks Pembangunan Manusia (IPM), dan
Gini Ratio.
Keempat variabel independen tersebut dipilih karena secara teoritis memiliki keterkaitan erat dengan tingkat kesejahteraan masyarakat. Tingkat pengangguran berpotensi meningkatkan kemiskinan karena berkurangnya sumber pendapatan rumah tangga. PDRB per kapita menggambarkan kapasitas ekonomi suatu daerah; semakin tinggi PDRB, umumnya semakin rendah tingkat kemiskinan. IPM mencerminkan kualitas hidup masyarakat dari aspek pendidikan, kesehatan, dan daya beli, sedangkan Gini Ratio mengukur ketimpangan distribusi pendapatan yang juga dapat memengaruhi tingkat kemiskinan.
Selain data atribut, penelitian ini juga menggunakan data spasial berupa peta batas administrasi kabupaten/kota di Jawa Barat. Data tersebut diperoleh dari geoBoundaries Indonesia Administrative Level 2 (ADM2) melalui situs https://www.geoboundaries.org dalam format shapefile (.shp). Data spasial ini digunakan untuk membentuk matriks bobot spasial (W) yang menggambarkan hubungan kedekatan antarwilayah berdasarkan batas administratif, yang kemudian menjadi dasar dalam analisis spasial.
Analisis dilakukan menggunakan pendekatan statistik spasial dengan bantuan perangkat lunak R. Pendekatan ini bertujuan untuk:
Mengidentifikasi pola persebaran spasial tingkat kemiskinan di Jawa Barat.
Menganalisis faktor-faktor sosial ekonomi yang memengaruhi kemiskinan.
Menentukan model spasial terbaik yang dapat menggambarkan keterkaitan antarwilayah secara signifikan.
Langkah pertama dalam analisis adalah melakukan analisis deskriptif statistik terhadap seluruh variabel untuk memberikan gambaran umum kondisi sosial ekonomi di masing-masing kabupaten/kota. Selanjutnya dilakukan uji autokorelasi spasial global menggunakan Moran’s I dan Geary’s C untuk mengetahui ada tidaknya keterkaitan spasial antarwilayah. Jika hasil menunjukkan autokorelasi positif yang signifikan, maka dilanjutkan dengan analisis autokorelasi spasial lokal (LISA) seperti Local Moran’s I dan Getis-Ord 𝐺∗ untuk mendeteksi klaster daerah dengan tingkat kemiskinan tinggi (hot spots) atau rendah (cold spots).
Tahap berikutnya dilakukan analisis regresi OLS (Ordinary Least Squares) sebagai model dasar untuk menguji hubungan antarvariabel. Hasil dari model OLS selanjutnya diuji dengan Lagrange Multiplier (LM Test) untuk menentukan apakah terdapat dependensi spasial pada variabel dependen (spatial lag dependence) atau pada residual (spatial error dependence).
Apabila hasil LM menunjukkan adanya efek spasial, maka dilakukan estimasi model ekonometrika spasial, yaitu:
Spatial Autoregressive Model (SAR) jika terdapat spatial lag dependence;
Spatial Error Model (SEM) jika terdapat spatial error dependence; atau
Model lain seperti Spatial Durbin Model (SDM) dan Spatial Durbin Error Model (SDEM) jika pengaruh spasial juga terdapat pada variabel independen.
Pemilihan model terbaik dilakukan berdasarkan nilai Akaike Information Criterion (AIC), tingkat signifikansi parameter spasial (ρ atau λ), serta kesesuaian teori ekonomi spasial yang mendukung hubungan antarwilayah.
Secara umum, tahapan penelitian ini dapat dijelaskan melalui alur berikut:
Pengumpulan Data Mengumpulkan data atribut sosial ekonomi dari BPS dan data spasial dari geoBoundaries.
Analisis Statistik Deskriptif Menganalisis karakteristik umum variabel penelitian, termasuk nilai rata-rata, minimum, maksimum, dan penyebaran data.
Estimasi Model OLS dan Uji Asumsi Klasik Melakukan regresi linier dasar serta menguji asumsi normalitas, multikolinearitas, dan homoskedastisitas.
Uji Autokorelasi Spasial (Moran’s I dan Geary’s C) Menguji apakah terdapat ketergantungan spasial antarwilayah.
Uji Lagrange Multiplier (LM) Menentukan model spasial yang paling sesuai (SAR, SEM, atau SDM).
Estimasi Model Spasial Mengestimasi parameter model spasial dan menilai signifikansi koefisien spasial.
Pemilihan Model Terbaik Membandingkan hasil estimasi berdasarkan nilai AIC dan tingkat signifikansi.
Interpretasi dan Penarikan Kesimpulan Menginterpretasikan hasil model untuk menjelaskan pengaruh antarvariabel dan hubungan spasial antarwilayah.
Analisis deskriptif dilakukan untuk memberikan gambaran umum mengenai karakteristik data yang digunakan dalam penelitian ini. Analisis ini bertujuan untuk memahami kondisi sosial ekonomi di setiap kabupaten/kota di Provinsi Jawa Barat pada tahun 2024, sebelum dilakukan analisis spasial yang lebih lanjut.
Tabel berikut menyajikan statistik deskriptif dari masing-masing variabel yang digunakan, meliputi nilai rata-rata, minimum, maksimum, dan simpangan baku:
# Membuat tabel deskriptif
library(knitr)
deskriptif <- data.frame(
Variabel = c("Kemiskinan (%)", "TPT (%)", "PDRB (Juta Rp)", "IPM", "Gini Ratio"),
Rata_rata = c(8.01, 6.52, 53.11, 74.68, 0.37),
Minimum = c(2.34, 1.58, 24.91, 68.89, 0.29),
Maksimum = c(11.93, 8.97, 147.08, 83.75, 0.48),
Simpangan_Baku = c(2.65, 1.71, 33.26, 4.33, 0.05)
)
kable(deskriptif, caption = "Statistik Deskriptif Variabel Penelitian Tahun 2024")| Variabel | Rata_rata | Minimum | Maksimum | Simpangan_Baku |
|---|---|---|---|---|
| Kemiskinan (%) | 8.01 | 2.34 | 11.93 | 2.65 |
| TPT (%) | 6.52 | 1.58 | 8.97 | 1.71 |
| PDRB (Juta Rp) | 53.11 | 24.91 | 147.08 | 33.26 |
| IPM | 74.68 | 68.89 | 83.75 | 4.33 |
| Gini Ratio | 0.37 | 0.29 | 0.48 | 0.05 |
Berdasarkan tabel tersebut, rata-rata tingkat kemiskinan di Jawa Barat tahun 2024 sebesar 8,01%, dengan rentang antara 2,34% hingga 11,93%. Hal ini menunjukkan adanya ketimpangan tingkat kemiskinan antarwilayah, di mana beberapa kabupaten memiliki tingkat kemiskinan yang jauh lebih tinggi dibandingkan kota-kota besar seperti Bandung atau Bekasi.
Tingkat Pengangguran Terbuka (TPT) menunjukkan rata-rata sebesar 6,52%, dengan variasi yang cukup tinggi antarwilayah (simpangan baku 1,71). Wilayah dengan TPT tinggi umumnya adalah daerah dengan tingkat industrialisasi menengah dan ketergantungan pada sektor informal.
Sementara itu, Produk Domestik Regional Bruto (PDRB) per kapita memiliki rata-rata Rp 53,11 juta dengan simpangan baku Rp 33,26 juta, yang menunjukkan kesenjangan ekonomi antarwilayah. Beberapa daerah perkotaan seperti Kota Bekasi, Kota Bandung, dan Kabupaten Karawang memiliki PDRB per kapita yang jauh lebih tinggi dibandingkan kabupaten di wilayah selatan Jawa Barat yang relatif tertinggal secara ekonomi.
Rata-rata Indeks Pembangunan Manusia (IPM) di Jawa Barat sebesar 74,68, menunjukkan kategori pembangunan manusia yang cukup tinggi. Namun, perbedaan nilai IPM antarwilayah masih cukup signifikan — dengan nilai minimum 68,89 dan maksimum 83,75 yang mengindikasikan masih adanya ketimpangan dalam akses pendidikan, kesehatan, dan standar hidup layak.
Sedangkan Gini Ratio rata-rata sebesar 0,37, mencerminkan tingkat ketimpangan pendapatan yang tergolong sedang. Nilai Gini tertinggi sebesar 0,48 menunjukkan bahwa di beberapa wilayah, distribusi pendapatan belum merata dan kesenjangan sosial ekonomi masih cukup menonjol.
| Moran.s.I | Z | p | Ket |
|---|---|---|---|
| 0.283 | 2.3277 | 0.009964 | Signifikan |
Nilai Moran’s I sebesar 0.283 dengan nilai p = 0.009964 menunjukkan adanya autokorelasi spasial positif. Artinya, wilayah dengan tingkat kemiskinan tinggi cenderung berdekatan dengan wilayah yang juga memiliki kemiskinan tinggi (high-high), dan demikian pula sebaliknya.
geary <- data.frame(
Variabel = c("Kemiskinan", "TPT", "PDRB", "IPM", "GINI"),
GearyC = c(0.410, 0.465, 0.660, 0.607, 0.504)
)
kable(geary, caption = "Nilai Geary’s C Tiap Variabel")| Variabel | GearyC |
|---|---|
| Kemiskinan | 0.410 |
| TPT | 0.465 |
| PDRB | 0.660 |
| IPM | 0.607 |
| GINI | 0.504 |
Nilai Geary’s C untuk kemiskinan sebesar 0.410 (< 1) juga menunjukkan adanya autokorelasi spasial positif. Hal ini memperkuat hasil Moran’s I bahwa wilayah-wilayah dengan kondisi sosial ekonomi mirip (tingkat kemiskinan tinggi/rendah) cenderung saling berdekatan. Variabel PDRB, IPM, dan Gini Ratio memiliki nilai C mendekati 0.5–0.6, menandakan hubungan spasial yang lebih lemah, namun masih ada kecenderungan keterkaitan antarwilayah.
Hasil analisis LISA (Local Moran’s I) mengungkapkan pola klaster yang spesifik:
Wilayah selatan Jawa Barat (Garut, Tasikmalaya, Ciamis) membentuk klaster high-high kemiskinan, menandakan daerah dengan kemiskinan tinggi dikelilingi oleh daerah dengan kemiskinan tinggi pula.
Wilayah utara seperti Bekasi, Karawang, dan Purwakarta membentuk klaster low-low, mencerminkan daerah dengan kondisi ekonomi relatif lebih baik. Pola ini memperlihatkan adanya segregasi spasial ekonomi antara daerah selatan dan utara provinsi, yang dapat disebabkan oleh perbedaan struktur ekonomi dan akses pembangunan.
Analisis Getis–Ord G* menunjukkan daerah dengan hot spot kemiskinan di bagian selatan (Garut–Tasikmalaya–Ciamis), dan cold spot di wilayah industri seperti Bekasi dan Karawang. Temuan ini selaras dengan hasil LISA, menegaskan adanya polaritas pembangunan di Jawa Barat, di mana kemajuan ekonomi lebih terpusat di kawasan metropolitan dan industri, sementara daerah selatan tertinggal.
\[ Y = 44.770 - 0.0879(TPT) - 0.0046(PDRB) - 0.5521(IPM) + 14.2366(GINI) \]
Secara umum, variabel IPM dan Gini Ratio berpengaruh signifikan terhadap tingkat kemiskinan. IPM berpengaruh negatif, artinya peningkatan kualitas pembangunan manusia mampu menekan kemiskinan. Sebaliknya, Gini Ratio berpengaruh positif, yang menunjukkan semakin tinggi ketimpangan pendapatan, semakin tinggi pula tingkat kemiskinan.
linear <- data.frame(
Statistik_Uji = "Ramsey RESET",
F = 0.4760,
p = 0.6281,
Ket = "Hubungan Linier"
)
kable(linear)| Statistik_Uji | F | p | Ket |
|---|---|---|---|
| Ramsey RESET | 0.476 | 0.6281 | Hubungan Linier |
Nilai p sebesar 0.6281 > 0.05 menandakan hubungan antarvariabel bersifat linier, sehingga model OLS dapat digunakan sebagai dasar analisis.
norm <- data.frame(
Statistik_Uji = "Shapiro-Wilk",
W = 0.95524,
p = 0.2865,
Ket = "Normal"
)
kable(norm)| Statistik_Uji | W | p | Ket |
|---|---|---|---|
| Shapiro-Wilk | 0.95524 | 0.2865 | Normal |
Residual berdistribusi normal (p > 0.05), sehingga asumsi normalitas terpenuhi.
homo <- data.frame(
Statistik_Uji = "Breusch-Pagan Test",
BP = 2.56,
p = 0.6339,
Ket = "Tidak heteroskedastis"
)
kable(homo)| Statistik_Uji | BP | p | Ket |
|---|---|---|---|
| Breusch-Pagan Test | 2.56 | 0.6339 | Tidak heteroskedastis |
Nilai p > 0.05 menunjukkan tidak ada heteroskedastisitas; variasi residual antarwilayah homogen.
resid <- data.frame(
Statistik_Uji = "Moran’s I Residual",
I = 2.3277,
p = 0.015,
Ket = "Signifikan"
)
kable(resid)| Statistik_Uji | I | p | Ket |
|---|---|---|---|
| Moran’s I Residual | 2.3277 | 0.015 | Signifikan |
Adanya autokorelasi spasial pada residual membuktikan bahwa model OLS belum menangkap efek spasial, sehingga perlu dilanjutkan ke model spasial (SAR/SEM).
lmtest <- data.frame(
Statistik_Uji = c("LM Error", "LM Lag", "Robust LM Error", "Robust LM Lag", "LM SARMA"),
Value = c(3.4944, 8.8489, 1.0849, 6.4394, 9.9338),
p = c(0.0616, 0.0029, 0.2976, 0.0112, 0.0070),
Ket = c("ns", "sig.", "ns", "sig.", "sig.")
)
kable(lmtest)| Statistik_Uji | Value | p | Ket |
|---|---|---|---|
| LM Error | 3.4944 | 0.0616 | ns |
| LM Lag | 8.8489 | 0.0029 | sig. |
| Robust LM Error | 1.0849 | 0.2976 | ns |
| Robust LM Lag | 6.4394 | 0.0112 | sig. |
| LM SARMA | 9.9338 | 0.0070 | sig. |
Berdasarkan hasil uji Lagrange Multiplier (LM), diketahui bahwa nilai LM Lag (p = 0.0029) dan Adjusted LM Lag (p = 0.0112) signifikan pada α = 0.05, sedangkan LM Error dan Adjusted LM Error tidak signifikan. Hal ini menunjukkan adanya pengaruh spasial pada variabel dependen (lag spasial), bukan pada error-nya. Dengan demikian, model spasial yang paling tepat digunakan adalah model Spatial Autoregressive (SAR). Hasil ini juga diperkuat oleh nilai LM SARMA (p = 0.007) yang signifikan, menegaskan adanya dependensi spasial secara keseluruhan dalam model.
model <- data.frame(
Model = c("OLS", "SAR", "SDM", "GNS", "SAC/SARAR", "SEM", "SDEM"),
LogLik = c(-50.19, -43.976, -40.163, -39.951, -43.449, -46.03, -40.166),
AIC = c(112.38, 101.953, 102.327, 103.903, 102.898, 106.06, 102.332),
Ket = c(
"Model dasar tanpa efek spasial.",
"AIC terendah, terbaik menangkap efek spasial.",
"Mempertimbangkan lag X dan Y.",
"Model kompleks dengan banyak parameter.",
"Gabungan lag dan error.",
"Efek spasial hanya pada error.",
"Mirip SDM tapi efek hanya pada error."
)
)
kable(model)| Model | LogLik | AIC | Ket |
|---|---|---|---|
| OLS | -50.190 | 112.380 | Model dasar tanpa efek spasial. |
| SAR | -43.976 | 101.953 | AIC terendah, terbaik menangkap efek spasial. |
| SDM | -40.163 | 102.327 | Mempertimbangkan lag X dan Y. |
| GNS | -39.951 | 103.903 | Model kompleks dengan banyak parameter. |
| SAC/SARAR | -43.449 | 102.898 | Gabungan lag dan error. |
| SEM | -46.030 | 106.060 | Efek spasial hanya pada error. |
| SDEM | -40.166 | 102.332 | Mirip SDM tapi efek hanya pada error. |
Berdasarkan hasil perbandingan berbagai model spasial, model Spatial Autoregressive (SAR) memiliki nilai AIC paling rendah (101.95) dibandingkan model lainnya. Hal ini menunjukkan bahwa SAR memberikan keseimbangan terbaik antara ketepatan model dan kompleksitas parameter. Model SAR juga didukung oleh hasil uji Lagrange Multiplier (LM) yang menunjukkan bahwa efek spasial signifikan terjadi pada lag variabel dependen, bukan pada error. Oleh karena itu, pemilihan model SAR dianggap paling tepat untuk menggambarkan hubungan spasial antarwilayah.
Berdasarkan hasil pengujian autokorelasi spasial, diketahui bahwa terdapat dependensi spasial baik pada variabel X maupun Y. Oleh karena itu, model Spatial Durbin Model (SDM) cocok digunakan karena mampu mengakomodasi pengaruh langsung dan tidak langsung antarwilayah.
Dengan estimasi parameter model SDM menggunakan metode Maximum Likelihood Estimation (MLE), persamaan model dapat ditulis sebagai berikut:
\[ \text{Kemiskinan}_i = 0.55\rho\sum_{j=1}^{n}w_{ij}y_j + 39.14 + 0.26(TPT)_i - 0.0033(PDRB)_i - 0.42(IPM)_i + 18.98(GINI)_i + 0.46\sum_{j=1}^{n}w_{ij}(TPT)_j + 0.01\sum_{j=1}^{n}w_{ij}(PDRB)_j - 0.06\sum_{j=1}^{n}w_{ij}(IPM)_j - 13.86\sum_{j=1}^{n}w_{ij}(GINI)_j \]
| Statistik_Uji | rho | LR | p | Ket |
|---|---|---|---|---|
| Likelihood Ratio | 0.55505 | 8.029 | 0.0046 | sig. |
Nilai ρ sebesar 0.55505 dengan p = 0.0046 menunjukkan bahwa terdapat spatial dependence yang signifikan. Artinya, kemiskinan di suatu wilayah dipengaruhi oleh kemiskinan di wilayah sekitarnya. Jika kemiskinan meningkat di satu kabupaten, maka wilayah tetangga cenderung mengalami peningkatan serupa. Selain itu, pola yang sama juga ditemukan pada variabel sosial ekonomi lain seperti TPT, PDRB, IPM, dan Gini Ratio, yang memiliki efek tidak hanya lokal tetapi juga lintas wilayah.
Dampak Langsung dan Tidak Langsung
impact_sdm <- data.frame(
Variabel = c("IPM", "Gini Ratio"),
Jenis_Impact = c("Direct", "Direct"),
Z = c(-4.626, 1.819),
P_value = c(0.000, 0.0689),
Ket = c("sig", "sig (α = 0.01)")
)
kable(impact_sdm, caption = "Dampak Variabel dalam Model SDM")| Variabel | Jenis_Impact | Z | P_value | Ket |
|---|---|---|---|---|
| IPM | Direct | -4.626 | 0.0000 | sig |
| Gini Ratio | Direct | 1.819 | 0.0689 | sig (α = 0.01) |
IPM memiliki dampak langsung negatif dan signifikan terhadap kemiskinan. Artinya, peningkatan kualitas pembangunan manusia berperan besar dalam menekan kemiskinan di daerah tersebut.
Gini Ratio berpengaruh positif dan signifikan pada α = 0.01, yang berarti ketimpangan pendapatan masih menjadi faktor penting yang memperburuk kemiskinan.
Tidak terdapat efek spillover yang kuat ke wilayah tetangga, menandakan pengaruh variabel lebih dominan secara lokal.
asumsi_sdm <- data.frame(
Uji = c("Residual Autokorelasi", "Normalitas", "Homoskedastisitas"),
Statistik = c(0.56865, 0.98449, 9.3304),
P_value = c(0.2848, 0.9466, 0.3152),
Ket = c("tidak terdapat autokorelasi spasial", "berdistribusi normal", "tidak terdapat heteroskedastisitas")
)
kable(asumsi_sdm, caption = "Hasil Uji Asumsi Model SDM")| Uji | Statistik | P_value | Ket |
|---|---|---|---|
| Residual Autokorelasi | 0.56865 | 0.2848 | tidak terdapat autokorelasi spasial |
| Normalitas | 0.98449 | 0.9466 | berdistribusi normal |
| Homoskedastisitas | 9.33040 | 0.3152 | tidak terdapat heteroskedastisitas |
Semua hasil uji menunjukkan p-value > 0.05, sehingga model SDM memenuhi seluruh asumsi dasar: residual tidak terautokorelasi, berdistribusi normal, dan bersifat homoskedastis. Dengan demikian, model SDM dapat dianggap valid dan reliabel.
Berdasarkan hasil perbandingan model, Spatial Autoregressive (SAR) terpilih sebagai model terbaik karena mampu menjelaskan adanya spatial dependence antar kabupaten/kota dan secara teoritis sesuai dengan hasil uji Lagrange Multiplier.
Persamaan model SAR adalah sebagai berikut:
\[ y_i = \rho \sum_{j=1}^{n} w_{ij}y_j + 28.449 + 0.1132(TPT)_i - 0.0087(PDRB)_i - 0.4132(IPM)_i + 12.705(GINI)_i \]
uji_sar <- data.frame(
Statistik_Uji = "Likelihood Ratio",
rho = 0.65269,
LR = 12.428,
p = 0.00042,
Ket = "sig."
)
kable(uji_sar, caption ="Hasil Estimasi Model SAR")| Statistik_Uji | rho | LR | p | Ket |
|---|---|---|---|---|
| Likelihood Ratio | 0.65269 | 12.428 | 0.00042 | sig. |
Nilai ρ sebesar 0.65269 signifikan pada α = 0.01, yang menunjukkan adanya pengaruh kuat antarwilayah (spatial dependence). Artinya, kemiskinan di satu daerah berpotensi memengaruhi daerah di sekitarnya secara langsung. Wilayah yang berdekatan secara geografis memiliki dinamika kemiskinan yang saling terhubung.
Dampak Langsung dan Tidak Langsung
impact_sar <- data.frame(
Variabel = c("IPM", "IPM", "IPM"),
Jenis_Impact = c("Direct", "Indirect", "Total"),
Koefisien = c(-0.4816, -0.7083, -1.1899),
P_value = c(0.0000, 0.504, 0.339),
Ket = c("sig", "not sig", "not sig")
)
kable(impact_sar, caption = "Dampak Langsung dan Tidak Langsung Variabel IPM dalam Model SAR")| Variabel | Jenis_Impact | Koefisien | P_value | Ket |
|---|---|---|---|---|
| IPM | Direct | -0.4816 | 0.000 | sig |
| IPM | Indirect | -0.7083 | 0.504 | not sig |
| IPM | Total | -1.1899 | 0.339 | not sig |
Direct effect: IPM berpengaruh negatif signifikan terhadap kemiskinan — semakin tinggi IPM, semakin rendah tingkat kemiskinan.
Indirect effect: pengaruh IPM terhadap wilayah tetangga tidak signifikan, menandakan belum adanya spillover effect pembangunan manusia lintas daerah.
Total effect: efek keseluruhan masih didominasi oleh pengaruh lokal.
asumsi_sar <- data.frame(
Uji = c("Residual Autokorelasi", "Normalitas", "Homoskedastisitas"),
Statistik = c(-0.111, 0.95524, 3.0169),
P_value = c(0.7077, 0.2865, 0.555),
Ket = c("not sig.", "sig.", "not sig.")
)
kable(asumsi_sar, caption = "Uji Asumsi Model SAR")| Uji | Statistik | P_value | Ket |
|---|---|---|---|
| Residual Autokorelasi | -0.11100 | 0.7077 | not sig. |
| Normalitas | 0.95524 | 0.2865 | sig. |
| Homoskedastisitas | 3.01690 | 0.5550 | not sig. |
Model SAR memenuhi asumsi klasik spasial: residual acak (tidak autokorelasi), data berdistribusi normal, dan varian residual homogen. Plot residual terhadap nilai prediksi juga menunjukkan pola acak di sekitar nol, yang menegaskan bahwa hubungan antarvariabel bersifat linier dan stabil.
bobot_sar <- data.frame(
Bobot = c("Queen", "Rook", "KNN-4"),
rho = c(0.6527, 0.6527, 0.4871),
AIC = c(101.9525, 101.9525, 105.5565),
LogLik = c(-43.97625, -43.97625, -45.77824)
)
kable(bobot_sar, caption = "Uji Sensitivitas Matriks Bobot Spasial pada Model SAR")| Bobot | rho | AIC | LogLik |
|---|---|---|---|
| Queen | 0.6527 | 101.9525 | -43.97625 |
| Rook | 0.6527 | 101.9525 | -43.97625 |
| KNN-4 | 0.4871 | 105.5565 | -45.77824 |
Uji sensitivitas menunjukkan bahwa penggunaan Queen dan Rook contiguity memberikan hasil identik, baik dari sisi nilai ρ maupun AIC. Ini menandakan model SAR stabil terhadap pilihan matriks bobot spasial. Sementara itu, hasil dari KNN-4 sedikit berbeda, dengan AIC lebih tinggi, menandakan bahwa kedekatan spasial di Jawa Barat lebih baik dijelaskan oleh hubungan administratif antarwilayah (berbagi batas) dibandingkan oleh jarak fisik murni.
Berdasarkan hasil prediksi model SDM dan SAR, keduanya menunjukkan pola konsentrasi kemiskinan yang relatif tinggi di wilayah Bogor, Sukabumi, dan Cianjur bagian utara, serta Subang dan Indramayu di bagian utara Jawa Barat. Sebaliknya, daerah dengan tingkat kemiskinan rendah umumnya berada di Bandung Raya (Kota Bandung, Kota Cimahi, Kabupaten Bandung, dan Kabupaten Bandung Barat) serta wilayah timur seperti Sumedang, Majalengka, dan Kuningan.
Jika dibandingkan, model SAR (Spatial Autoregressive) memberikan hasil prediksi yang lebih halus dan konsisten di seluruh kabupaten/kota dibandingkan model SDM (Spatial Durbin Model) yang masih menunjukkan variasi ekstrem antarwilayah. Misalnya, pada model SDM terdapat prediksi kemiskinan yang terlalu tinggi di sebagian kecil wilayah utara dan terlalu rendah di selatan. Hal ini mengindikasikan bahwa model SAR lebih mampu menangkap keterkaitan spasial antarwilayah secara proporsional tanpa menghasilkan prediksi yang berlebihan.
Hasil ini juga diperkuat oleh analisis residual. Pada model SDM, masih terdapat pola spasial pada residual terutama di Bekasi dan Bogor bagian utara (yang cenderung underfit atau terprediksi terlalu rendah) serta Tasikmalaya dan Ciamis (yang overfit atau terprediksi terlalu tinggi). Hal ini menunjukkan bahwa SDM belum sepenuhnya mengakomodasi pengaruh spasial antarwilayah.
Sebaliknya, pada model SAR, distribusi residual tampak lebih acak dan tidak menunjukkan pola spasial yang jelas, dengan nilai residual yang relatif kecil di hampir seluruh kabupaten/kota. Hal ini menunjukkan bahwa model SAR lebih baik dalam menjelaskan variasi spasial kemiskinan, menghasilkan prediksi yang lebih stabil dan akurat, serta menggambarkan hubungan spasial antarwilayah dengan lebih realistis.
Dengan demikian, dapat disimpulkan bahwa model Spatial Autoregressive (SAR) merupakan model paling sesuai untuk menggambarkan pola kemiskinan di Provinsi Jawa Barat karena mampu menjelaskan keterkaitan spasial secara menyeluruh antara wilayah berpendapatan tinggi dan wilayah dengan tingkat kemiskinan yang lebih besar.
Berdasarkan hasil perbandingan model, model Spatial Autoregressive (SAR) terpilih sebagai model terbaik dengan nilai AIC terendah (101.95) dan hasil prediksi yang stabil. Model SAR untuk kemiskinan di Jawa Barat dapat dituliskan sebagai:
\[ y_i = \rho \sum_{j=1}^{n} w_{ij} y_j + 28.449 + 0.1132(TPT)_i - 0.0087(PDRB)_i - 0.4132(IPM)_i + 12.705(GINI)_i \]
Model ini paling sesuai untuk menganalisis faktor sosial ekonomi yang
memengaruhi kemiskinan di Jawa Barat tahun 2024. Hasil menunjukkan
bahwa: 1. IPM berpengaruh negatif dan signifikan
terhadap kemiskinan, artinya peningkatan pembangunan manusia dapat
menurunkan tingkat kemiskinan.
2. Gini Ratio berpengaruh positif signifikan pada taraf
signifikansi 0.01, menandakan bahwa ketimpangan pendapatan memperburuk
kondisi kemiskinan.
3. TPT dan PDRB per kapita tidak
berpengaruh signifikan setelah memperhitungkan efek spasial.
4. Secara spasial, pola kemiskinan menunjukkan bahwa wilayah
selatan Jawa Barat (seperti Garut, Tasikmalaya, dan
Ciamis) cenderung memiliki kemiskinan lebih tinggi dibandingkan wilayah
utara dan barat yang lebih maju secara ekonomi.
Dengan demikian, model SAR berhasil menangkap spatial dependence yang nyata antar kabupaten/kota dan memberikan gambaran yang lebih akurat tentang penyebaran kemiskinan di Jawa Barat.
Kebijakan Pembangunan Manusia: Pemerintah Provinsi Jawa Barat perlu menempatkan Indeks Pembangunan Manusia (IPM) sebagai fokus utama kebijakan pengentasan kemiskinan. Peningkatan IPM melalui tiga dimensi utama yaitu pendidikan, kesehatan, dan standar hidup layak — terbukti efektif dalam menekan angka kemiskinan.
Pengurangan Ketimpangan Ekonomi: Berdasarkan variabel Gini Ratio, kebijakan ekonomi sebaiknya diarahkan pada peningkatan pemerataan pendapatan dan kesempatan kerja agar jarak antara kelompok kaya dan miskin tidak semakin lebar. Program kewirausahaan dan pelatihan kerja bagi kelompok berpendapatan rendah dapat menjadi strategi efektif.
Fokus Wilayah Selatan Jawa Barat:
Daerah seperti Garut, Tasikmalaya, dan Ciamis memerlukan perhatian lebih
dalam program pembangunan infrastruktur, pendidikan, dan pemberdayaan
ekonomi, mengingat daerah-daerah tersebut konsisten menunjukkan tingkat
kemiskinan tinggi dan ketertinggalan ekonomi.
Penelitian Lanjutan: Untuk penelitian berikutnya, disarankan menambahkan variabel-variabel sosial ekonomi lainnya seperti akses air bersih, listrik, perumahan layak, dan kepadatan penduduk agar analisis lebih komprehensif dan representatif terhadap kondisi riil masyarakat.
Kesimpulan Akhir:
Model Spatial Autoregressive (SAR) dapat dinyatakan sebagai model paling
tepat dalam menggambarkan keterkaitan spasial antarwilayah serta
faktor-faktor utama yang memengaruhi kemiskinan di Provinsi Jawa Barat.
Hasil ini diharapkan dapat menjadi dasar empiris dalam perumusan
kebijakan pembangunan yang lebih inklusif dan berkeadilan
antarwilayah.
link dashboard (https://dinarolivia.shinyapps.io/dashboard/) script analysis :
library(sf); library(dplyr); library(stringr) library(spdep);
library(spatialreg) library(tmap); library(broom); library(ggplot2)
library(psych); library(readxl); library(viridis)
library(scales); library(tibble); library(purrr) library(forcats);
library(htmltools); library(lmtest)
indo_adm2 <- st_read(“data/geoBoundaries-IDN-ADM2.shp”, quiet = TRUE)
kab_jabar <- c(“Bandung”,“Bandung Barat”,“Bekasi”,“Bogor”,“Ciamis”,“Cianjur”,“Cirebon”, “Garut”,“Indramayu”,“Karawang”,“Kuningan”,“Majalengka”,“Pangandaran”, “Purwakarta”,“Subang”,“Sukabumi”,“Sumedang”,“Tasikmalaya”, “Kota Bandung”,“Kota Bekasi”,“Kota Bogor”,“Kota Cimahi”, “Kota Cirebon”,“Kota Depok”,“Kota Sukabumi”,“Kota Tasikmalaya”, “Kota Banjar”)
indo_adm2 <- indo_adm2 %>% mutate(shapeName_clean = str_to_lower(str_squish(shapeName)))
kab_jabar_clean <- str_to_lower(kab_jabar)
jabar <- indo_adm2 %>% filter(shapeName_clean %in% kab_jabar_clean)
if (nrow(jabar) == 0) { message(“⚠️ Tidak ada match exact; cek kembali nama kabupaten/kota di shapefile.”) print(unique(indo_adm2$shapeName)[1:50]) stop(“Periksa & adaptasi nama kab/kota, lalu jalankan ulang.”) }
attr <- read_xlsx(“data/DATA SPASIAL.xlsx”)
print(names(attr))
attr <- attr %>% rename(kabupaten = starts_with(“KAB”)) %>% mutate(kab_clean = str_to_lower(str_squish(kabupaten)))
jabar <- jabar %>% left_join(attr, by = c(“shapeName_clean” = “kab_clean”))
glimpse(jabar)
jabar <- jabar %>% left_join(attr, by = c(“shapeName_clean” = “kab_clean”)) jabar
data_spasial <- jabar %>% sf::st_drop_geometry() %>%
dplyr::select(KEMISKINAN, TPT, PDRB, IPM, GINI)
psych::describe(data_spasial)
ggplot(data_spasial, aes(x = KEMISKINAN)) + geom_histogram(bins = 10, fill = “skyblue”, color = “black”) + labs(title = “Distribusi Tingkat Kemiskinan”, x = “Kemiskinan (%)”, y = “Frekuensi”)
ggplot(data_spasial, aes(x = TPT)) + geom_histogram(bins = 10, fill = “lightgreen”, color = “black”) + labs(title = “Distribusi Tingkat Pengangguran Terbuka”, x = “TPT (%)”, y = “Frekuensi”)
ggplot(data_spasial, aes(x = PDRB)) + geom_histogram(bins = 10, fill = “lightcoral”, color = “black”) + labs(title = “Distribusi PDRB per Kapita”, x = “PDRB (juta rupiah)”, y = “Frekuensi”)
ggplot(data_spasial, aes(x = IPM)) + geom_histogram(bins = 10, fill = “gold”, color = “black”) + labs(title = “Distribusi Indeks Pembangunan Manusia”, x = “IPM”, y = “Frekuensi”)
ggplot(data_spasial, aes(x = GINI)) + geom_histogram(bins = 10, fill = “plum”, color = “black”) + labs(title = “Distribusi Rasio Gini”, x = “GINI”, y = “Frekuensi”)
ggplot(data_spasial, aes(y = KEMISKINAN)) + geom_boxplot(fill = “skyblue”, color = “black”) + labs(title = “Boxplot Tingkat Kemiskinan”, y = “Kemiskinan (%)”, x = ““) + theme_minimal()
ggplot(data_spasial, aes(y = TPT)) + geom_boxplot(fill = “lightgreen”, color = “black”) + labs(title = “Boxplot Tingkat Pengangguran Terbuka”, y = “TPT (%)”, x = ““) + theme_minimal()
ggplot(data_spasial, aes(y = PDRB)) + geom_boxplot(fill = “lightcoral”, color = “black”) + labs(title = “Boxplot PDRB per Kapita”, y = “PDRB (juta rupiah)”, x = ““) + theme_minimal()
ggplot(data_spasial, aes(y = IPM)) + geom_boxplot(fill = “gold”, color = “black”) + labs(title = “Boxplot Indeks Pembangunan Manusia”, y = “IPM”, x = ““) + theme_minimal()
ggplot(data_spasial, aes(y = GINI)) + geom_boxplot(fill = “plum”, color = “black”) + labs(title = “Boxplot Rasio Gini”, y = “GINI”, x = ““) + theme_minimal()
ggplot(data_spasial, aes(x = TPT, y = KEMISKINAN)) + geom_point(color = “steelblue”, size = 3) + geom_smooth(method = “lm”, se = FALSE, color = “darkred”, linetype = “dashed”) + labs(title = “Hubungan antara Kemiskinan dan TPT”, x = “Tingkat Pengangguran Terbuka (%)”, y = “Kemiskinan (%)”) + theme_minimal()
ggplot(data_spasial, aes(x = PDRB, y = KEMISKINAN)) + geom_point(color = “seagreen”, size = 3) + geom_smooth(method = “lm”, se = FALSE, color = “darkred”, linetype = “dashed”) + labs(title = “Hubungan antara Kemiskinan dan PDRB per Kapita”, x = “PDRB per Kapita (juta rupiah)”, y = “Kemiskinan (%)”) + theme_minimal()
ggplot(data_spasial, aes(x = IPM, y = KEMISKINAN)) + geom_point(color = “goldenrod”, size = 3) + geom_smooth(method = “lm”, se = FALSE, color = “darkred”, linetype = “dashed”) + labs(title = “Hubungan antara Kemiskinan dan IPM”, x = “Indeks Pembangunan Manusia”, y = “Kemiskinan (%)”) + theme_minimal()
ggplot(data_spasial, aes(x = GINI, y = KEMISKINAN)) + geom_point(color = “purple”, size = 3) + geom_smooth(method = “lm”, se = FALSE, color = “darkred”, linetype = “dashed”) + labs(title = “Hubungan antara Kemiskinan dan Rasio Gini”, x = “Rasio Gini”, y = “Kemiskinan (%)”) + theme_minimal()
theme_map <- function() { theme_minimal(base_size = 11) + theme( axis.title = element_blank(), axis.text = element_blank(), panel.grid = element_blank(), legend.position = “right”, plot.title = element_text(face = “bold”) ) }
p_kemiskinan <- ggplot(jabar) + geom_sf(aes(fill = KEMISKINAN), color = “white”, size = 0.2) + scale_fill_viridis_c(name = “Kemiskinan (%)”, option = “C”, direction = -1) + labs(title = “Penyebaran Kemiskinan (%) – Jawa Barat”) + coord_sf() + theme_map()
p_tpt <- ggplot(jabar) + geom_sf(aes(fill = TPT), color = “white”, size = 0.2) + scale_fill_viridis_c(name = “TPT (%)”, option = “C”, direction = -1) + labs(title = “Penyebaran TPT (%) – Jawa Barat”) + coord_sf() + theme_map()
p_pdrb <- ggplot(jabar) + geom_sf(aes(fill = PDRB), color = “white”, size = 0.2) + scale_fill_viridis_c( name = “PDRB per Kapita(juta rupiah, log10)”, option = “C”, direction = -1, trans = “log10”, labels = label_number(accuracy = 0.1) ) + labs(title = “Penyebaran PDRB per Kapita – Jawa Barat”) + coord_sf() + theme_map()
p_ipm <- ggplot(jabar) + geom_sf(aes(fill = IPM), color = “white”, size = 0.2) + scale_fill_viridis_c(name = “IPM”, option = “C”, direction = -1) + labs(title = “Penyebaran IPM – Jawa Barat”) + coord_sf() + theme_map()
p_gini <- ggplot(jabar) + geom_sf(aes(fill = GINI), color = “white”, size = 0.2) + scale_fill_viridis_c(name = “Rasio Gini”, option = “C”, direction = -1) + labs(title = “Penyebaran Rasio Gini – Jawa Barat”) + coord_sf() + theme_map()
p_kemiskinan p_tpt p_pdrb p_ipm p_gini
nb_queen <- poly2nb(jabar, queen = TRUE) lw_queen <- nb2listw(nb_queen, style = “W”, zero.policy = TRUE)
vars <- c(“KEMISKINAN”, “TPT”, “PDRB”, “IPM”, “GINI”)
safe_moran_row_both <- function(y, varname, nb_full, nsim = 999, seed = 77) { # Tangani NA: subset data & struktur tetangga ok <- is.finite(y) y_use <- y[ok] nb_sub <- if (all(ok)) nb_full else spdep::subset.nb(nb_full, ok) lw_sub <- nb2listw(nb_sub, style = “W”, zero.policy = TRUE)
mt <- moran.test(y_use, lw_sub, zero.policy = TRUE) est <- mt$estimate
I_asym <- suppressWarnings(as.numeric(est[[“Moran I statistic”]])) EI <- suppressWarnings(as.numeric(est[[“Expectation”]])) VarI <- est[[“Variance”]] if (is.null(VarI)) VarI <- est[[“variance of Moran’s I”]] VarI <- suppressWarnings(as.numeric(VarI)) Z_asym <- suppressWarnings(as.numeric(mt\(statistic)) SD_asym <- if (!is.null(VarI) && is.finite(VarI)) sqrt(VarI) else NA_real_ p_asym <- mt\)p.value
set.seed(seed) mmc <- moran.mc(y_use, lw_sub, nsim = nsim, zero.policy = TRUE, alternative = “two.sided”) I_perm <- suppressWarnings(as.numeric(mmc\(statistic)) p_perm <- mmc\)p.value
tibble( Variabel = varname, Moran_I = round(I_asym, 6), Expected_I = round(EI, 6), SD = round(SD_asym, 6), Z_value = round(Z_asym, 6), P_asym = signif(p_asym, 6), I_perm = round(I_perm, 6), P_perm = signif(p_perm, 6), nsim = nsim ) }
hasil_moran_both <- map_dfr(vars, ~ safe_moran_row_both(jabar[[.x]], .x, nb_queen))
hasil_moran_both
safe_geary_row <- function(y, varname, nb_full) { # Tangani NA: subset data & struktur tetangga ok <- is.finite(y) if (!all(ok)) { # subset nb (fungsi subset.nb ada di spdep) nb_sub <- spdep::subset.nb(nb_full, ok) } else { nb_sub <- nb_full } lw_sub <- nb2listw(nb_sub, style = “W”, zero.policy = TRUE)
gt <- geary.test(y[ok], lw_sub, randomisation = TRUE, zero.policy = TRUE) est <- gt$estimate
C_asym <- suppressWarnings( as.numeric(est[[“Geary C statistic”]]) ) if (is.na(C_asym)) { C_asym <- suppressWarnings(as.numeric(est[[“Geary’s C”]])) } E_C <- suppressWarnings(as.numeric(est[[“Expectation”]])) if (is.na(E_C)) E_C <- 1 # default teoritis Geary’s C
Z_asym <- suppressWarnings(as.numeric(gt\(statistic)) p_asym <- gt\)p.value
set.seed(77) gmc <- geary.mc(y[ok], lw_sub, nsim = 999, zero.policy = TRUE) C_perm <- suppressWarnings(as.numeric(gmc\(statistic)) p_perm <- gmc\)p.value
tibble( Variabel = varname, Geary_C = round(C_asym, 6), Expected_C = round(E_C, 6), Z_value = round(Z_asym, 6), P_asym = signif(p_asym, 6), C_perm = round(C_perm, 6), P_perm = signif(p_perm, 6) ) }
hasil_geary <- map_dfr(vars, ~ safe_geary_row(jabar[[.x]], .x, nb_queen))
hasil_geary
safe_lisa <- function(y, varname, sf_obj, lw, sig_thr = 0.05) { ok <- is.finite(y) z <- as.numeric(scale(y[ok])) lag_z <- lag.listw(lw, z, zero.policy = TRUE)
lisa_df <- as.data.frame(localmoran(z, lw, zero.policy = TRUE))
# Hitung p-value dua sisi (karena versi spdep bisa beda nama) if (“Pr(z > 0)” %in% names(lisa_df)) { p_right <- lisa_df[,“Pr(z > 0)”] p_left <- 1 - p_right p_two <- 2pmin(p_right, p_left) } else if (“Pr(z <= 0)” %in% names(lisa_df)) { p_left <- lisa_df[,“Pr(z <= 0)”] p_right <- 1 - p_left p_two <- 2pmin(p_right, p_left) } else if (“Pr(z != E(Ii))” %in% names(lisa_df)) { p_two <- lisa_df[,“Pr(z != E(Ii))”] } else { zval <- lisa_df[,“Z.Ii”] p_two <- 2*pnorm(-abs(zval)) }
# Klasifikasi kuadran quad_raw <- ifelse(z >= 0 & lag_z >= 0, “High-High”, ifelse(z < 0 & lag_z < 0, “Low-Low”, ifelse(z >= 0 & lag_z < 0, “High-Low”, “Low-High”))) cluster <- ifelse(p_two <= sig_thr, quad_raw, “Not Significant”)
sf_sub <- sf_obj[ok, ] %>% mutate( z_score = z, lag_z = lag_z, Ii = lisa_df\(Ii, Z_Ii = lisa_df\)Z.Ii, p_lisa = p_two, LISA_cat = factor(cluster, levels = c(“High-High”,“Low-Low”,“High-Low”,“Low-High”,“Not Significant”)), Variabel = varname )
# Tabel ringkasan kategori tab <- sf_sub %>% st_drop_geometry() %>% group_by(LISA_cat) %>% summarise( Kabupaten = paste(kabupaten, collapse = “,”), Jumlah = n(), .groups = “drop” )
list(sf = sf_sub, tabel = tab) }
lisa_results <- map(vars, ~ safe_lisa(jabar[[.x]], .x, jabar, lw_queen))
lisa_tables <- map2_df(lisa_results, vars, ~ mutate(.x$tabel, Variabel = .y)) lisa_tables
plots <- map2(lisa_results, vars, ~ { ggplot(.x$sf) + geom_sf(aes(fill = LISA_cat), color = “white”, size = 0.2) + scale_fill_manual(values = c( “High-High” = “#d7191c”, “Low-Low” = “#2c7bb6”, “High-Low” = “#fdae61”, “Low-High” = “#abd9e9”, “Not Significant” = “grey80” ), drop = FALSE, name = “Cluster”) + labs(title = paste0(“Peta Cluster LISA -”, .y), subtitle = “Metode Local Moran’s I (Queen Contiguity)”) + theme_minimal() + theme(axis.text = element_blank(), axis.title = element_blank(), panel.grid = element_blank()) })
print(lisa_tables) plots[[1]] # tampilkan peta KEMISKINAN plots[[2]] # tampilkan peta TPT plots[[3]] # tampilkan peta PDRB plots[[4]] # tampilkan peta IPM plots[[5]] # tampilkan peta GINI
cols_gstar <- c(“Hot spot” = “#D7191C”, “Cold spot” = “#1B9E77”, “Not significant” = “grey99”)
compute_gistar <- function(sf_obj, varname, lw) { x_raw <- sf_obj[[varname]] ok <- is.finite(x_raw)
# Subset bobot & data bila ada NA nb_sub <- if (all(ok)) lw\(neighbours else spdep::subset.nb(lw\)neighbours, ok) lw_sub <- nb2listw(nb_sub, style = “W”, zero.policy = TRUE)
# Matriks bobot (untuk G, G*) Wb <- spdep::listw2mat(lw_sub) sum_x <- sum(x_raw[ok])
# — G_i (tanpa i) — num_G <- as.numeric(Wb %*% x_raw[ok]) den_G <- (sum_x - x_raw[ok]) G_raw <- num_G / den_G
# — G_i* (dengan i): set w_ii = 1 — Wb_star <- Wb; diag(Wb_star) <- 1 num_Gs <- as.numeric(Wb_star %*% x_raw[ok]) den_Gs <- sum_x G_star_raw <- num_Gs / den_Gs
# — Z-score Gi* (uji lokal) — Gz <- spdep::localG(x_raw[ok], listw = lw_sub, zero.policy = TRUE)
# Susun vektor full-length (NA di posisi non-ok) z_Gistar_full <- rep(NA_real_, nrow(sf_obj)); z_Gistar_full[ok] <- as.numeric(Gz) G_raw_full <- rep(NA_real_, nrow(sf_obj)); G_raw_full[ok] <- G_raw G_star_raw_full <- rep(NA_real_, nrow(sf_obj)); G_star_raw_full[ok] <- G_star_raw
# Kategori GSTAR_cat <- case_when( z_Gistar_full >= 1.96 ~ “Hot spot”, z_Gistar_full <= -1.96 ~ “Cold spot”, TRUE ~ “Not significant” )
# Gabungkan ke sf sf_out <- sf_obj |> mutate( !!paste0(“G_raw_”, varname) := G_raw_full, !!paste0(“Gstar_”, varname) := G_star_raw_full, !!paste0(“z_Gi_”, varname) := z_Gistar_full, !!paste0(“GSTAR_”, varname) := factor(GSTAR_cat, levels = c(“Hot spot”,“Cold spot”,“Not significant”)), Variabel = varname )
# Tabel ringkasan hotspot/coldspot tab <- sf_out |> st_drop_geometry() |> transmute( Variabel = varname, Kabupaten = kabupaten, Kategori = .data[[paste0(“GSTAR_”, varname)]], Gi_star = round(.data[[paste0(“z_Gi_”, varname)]], 3) ) |> filter(Kategori %in% c(“Hot spot”,“Cold spot”)) |> arrange(Kategori, dplyr::desc(if_else(Kategori == “Hot spot”, Gi_star, -Gi_star)))
list(sf = sf_out, table = tab) }
gistar_results <- map(vars, ~ compute_gistar(jabar, .x, lw_queen))
gistar_tables <- map_dfr(gistar_results, “table”) gistar_tables
cols_gstar <- c( “Hot spot” = “#D7191C”, “Cold spot” = “#1B9E77”, “Not significant” = “#EAEAEA” # bukan putih; abu muda biar batas terlihat )
sf_plot <- gistar_results[[which(vars == “KEMISKINAN”)]]$sf batas_prov <- sf::st_as_sf(sf::st_cast(sf::st_union(sf_plot), “MULTILINESTRING”))
fill_col <- paste0(“GSTAR_”, “KEMISKINAN”) sf_plot <- sf_plot |> mutate( !!fill_col := as.character(.data[[fill_col]]), !!fill_col := ifelse(is.na(.data[[fill_col]]), “Not significant”, .data[[fill_col]]), !!fill_col := forcats::fct(.data[[fill_col]], levels = c(“Hot spot”,“Cold spot”,“Not significant”)) )
gistar_plots <- purrr::map2(gistar_results, vars, ~{ sf_plot <- .x$sf fill_col <- paste0(“GSTAR_”, .y)
batas_prov <- sf::st_as_sf(sf::st_cast(sf::st_union(sf_plot), “MULTILINESTRING”))
sf_plot <- sf_plot |> mutate( !!fill_col := as.character(.data[[fill_col]]), !!fill_col := ifelse(is.na(.data[[fill_col]]), “Not significant”, .data[[fill_col]]), !!fill_col := forcats::fct(.data[[fill_col]], levels = c(“Hot spot”,“Cold spot”,“Not significant”)) )
ggplot(sf_plot) + geom_sf(aes(fill = .data[[fill_col]]), color = “#FFFFFF”, size = 0.3) + geom_sf(data = sf::st_boundary(sf_plot), color = “#6B7280”, size = 0.35, fill = NA) + geom_sf(data = batas_prov, color = “#374151”, size = 0.6, fill = NA) + scale_fill_manual(values = cols_gstar, limits = c(“Hot spot”,“Cold spot”,“Not significant”), drop = FALSE, name = “Kategori”) + labs(title = paste0(“Getis–Ord G* –”, .y)) + coord_sf(expand = FALSE) + theme_minimal(base_size = 12) + theme(axis.title = element_blank(), axis.text = element_blank(), panel.grid = element_blank()) }) gistar_plots[[1]] # KEMISKINAN gistar_plots[[2]] # TPT gistar_plots[[3]] # PDRB gistar_plots[[4]] # IPM gistar_plots[[5]] # GINI
ols_jabar <- lm(KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar) summary(ols_jabar)
car::vif(ols_jabar) # Terpenuhi < 5
shapiro.test(residuals(ols_jabar)) # Terpenuhi terima H0
lmtest::bptest(ols_jabar) # Terpenuhi terima H0
moran.test(residuals(ols_jabar), lw_queen) moran.mc(residuals(ols_jabar), lw_queen, nsim = 999)
lm.LMtests(ols_jabar, lw_queen, test = “all”)
sdm_jabar <- lagsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, Durbin = TRUE, # aktifkan Spatial Lag of X method = “eigen”, # metode umum untuk estimasi SDM zero.policy = TRUE )
summary(sdm_jabar)
impacts_sdm <- impacts(sdm_jabar, listw = lw_queen, R = 999) summary(impacts_sdm, zstats = TRUE, short = FALSE)
gns_jabar <- sacsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, Durbin = TRUE, # aktifkan Spatial Lag of X method = “eigen”, # metode numerik untuk dekomposisi spasial zero.policy = TRUE )
summary(gns_jabar) impacts_gns <- impacts(gns_jabar, listw = lw_queen, R = 999) summary(impacts_gns, zstats = TRUE, short = TRUE)
sac_jabar <- sacsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, method = “eigen”, # metode numerik dekomposisi eigen zero.policy = TRUE )
summary(sac_jabar)
sar_jabar <- lagsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, method = “eigen”, zero.policy = TRUE )
summary(sar_jabar)
imp_sar <- impacts( sar_jabar, listw = lw_queen, R = 500, # jumlah replikasi simulasi (biasanya 500–1000) zero.policy = TRUE )
summary(imp_sar, zstats = TRUE, short = TRUE)
y <- jabar$KEMISKINAN
R2_sar <- 1 - (sar_jabar$s2 / var(y)) R2_sar
fitMeasures(sar_jabar)
slx_jabar <- lmSLX( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, Durbin = TRUE, # aktifkan lag pada variabel X (WX) zero.policy = TRUE )
summary(slx_jabar)
impacts_slx <- impacts(slx_jabar, listw = lw_queen, R = 999) summary(impacts_slx, zstats = TRUE, short = FALSE)
sem_jabar <- errorsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, method = “eigen”, zero.policy = TRUE )
summary(sem_jabar)
sdem_jabar <- errorsarlm( KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen, Durbin = TRUE, # aktifkan lag pada variabel X method = “eigen”, zero.policy = TRUE )
summary(sdem_jabar)
compare_stats <- data.frame( Model = c(“OLS”,“SLX”,“SAR”,“SDM”,“GNS”,“SAC/SARAR”,“SEM”,“SDEM”), LogLik = c(logLik(ols_jabar), logLik(slx_jabar), logLik(sar_jabar), logLik(sdm_jabar), logLik(gns_jabar), logLik(sac_jabar), logLik(sem_jabar), logLik(sdem_jabar)), AIC = c(AIC(ols_jabar), AIC(slx_jabar), AIC(sar_jabar), AIC(sdm_jabar), AIC(gns_jabar), AIC(sac_jabar), AIC(sem_jabar), AIC(sdem_jabar)) ) compare_stats
jabar_pred <- jabar %>% mutate( y_hat = as.numeric(fitted(sar_jabar)), # prediksi SAR resid = as.numeric(residuals(sar_jabar)), # residual z_resid = as.numeric(scale(resid)) # residual distandarkan (opsional) )
p_pred <- ggplot(jabar_pred) + geom_sf(aes(fill = y_hat), color = “white”, size = 0.2) + scale_fill_viridis_c(name = “Prediksi(%)”, option = “C”, direction = -1) + labs(title = “Peta Prediksi (SAR) – Kemiskinan Jawa Barat”) + theme_minimal(base_size = 12) + theme(axis.title = element_blank(), axis.text = element_blank(), panel.grid = element_blank()) p_pred
p_kemiskinan
lim_resid <- max(abs(jabar_pred$z_resid), na.rm = TRUE)
p_resid <- ggplot(jabar_pred) + geom_sf(aes(fill = z_resid), color = “white”, size = 0.2) + scale_fill_gradient2( name = “Residual (z)”, low = “#2c7bb6”, mid = “white”, high = “#d7191c”, midpoint = 0, limits = c(-lim_resid, lim_resid), oob = squish ) + labs(title = “Peta Residual Standar (SAR) – Kemiskinan Jawa Barat”, subtitle = “Merah: terprediksi terlalu rendah (underfit) • Biru: terprediksi terlalu tinggi (overfit)”) + theme_minimal(base_size = 12) + theme(axis.title = element_blank(), axis.text = element_blank(), panel.grid = element_blank()) p_resid
moran_res_sar <- moran.test(jabar_pred$resid, lw_queen, zero.policy = TRUE) moran_res_sar
jabar_pred <- jabar_pred %>% mutate(res_cat = cut(z_resid, breaks = c(-Inf, -1.96, -1, 1, 1.96, Inf), labels = c(“<= -1.96”, “(-1.96,-1]”, “(-1,1]”, “(1,1.96]”, “>= 1.96”), include.lowest = TRUE))
p_resid_disc <- ggplot(jabar_pred) + geom_sf(aes(fill = res_cat), color = “white”, size = 0.2) + scale_fill_brewer(name = “Residual (z)”, palette = “RdBu”, direction = -1, drop = FALSE) + labs(title = “Peta Residual Standar (SAR) – Klasifikasi Diskrit”) + theme_minimal(base_size = 12) + theme(axis.title = element_blank(), axis.text = element_blank(), panel.grid = element_blank()) p_resid_disc
labels <- sprintf( “%s
Prediksi: %.2f
Aktual:
%.2f
Residual (z): %.2f”, jabar_pred\(kabupaten, jabar_pred\)y_hat,
jabar_pred\(KEMISKINAN, jabar_pred\)z_resid ) |>
lapply(HTML) # <- list of HTML snippets, panjang =
nrow(jabar_pred)
pal_pred <- colorNumeric(“viridis”, domain = jabar_pred$y_hat)
leaflet(jabar_pred) %>% addProviderTiles(providers$CartoDB.Positron) %>% addPolygons( fillColor = ~pal_pred(y_hat), color = “white”, weight = 1, fillOpacity = 0.85, label = labels, # <- pakai vektor ‘labels’, bukan ~HTML(…) labelOptions = labelOptions(direction = “auto”), highlightOptions = highlightOptions(weight = 2, bringToFront = TRUE) ) %>% addLegend(“bottomright”, pal = pal_pred, values = ~y_hat, title = “Prediksi (%)”)
nb_queen <- poly2nb(jabar, queen = TRUE) lw_queen <- nb2listw(nb_queen, style = “W”)
nb_rook <- poly2nb(jabar, queen = FALSE) lw_rook <- nb2listw(nb_rook, style = “W”)
coords <- sf::st_coordinates(sf::st_centroid(jabar)) nb_knn4 <- knn2nb(knearneigh(coords, k = 4)) lw_knn4 <- nb2listw(nb_knn4, style = “W”)
sar_queen <- lagsarlm(KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_queen) sar_rook <- lagsarlm(KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_rook) sar_knn4 <- lagsarlm(KEMISKINAN ~ TPT + PDRB + IPM + GINI, data = jabar, listw = lw_knn4)
data.frame( Bobot = c(“Queen”, “Rook”, “KNN-4”), Rho = c(sar_queen\(rho, sar_rook\)rho, sar_knn4$rho), AIC = c(AIC(sar_queen), AIC(sar_rook), AIC(sar_knn4)), LogLik = c(logLik(sar_queen), logLik(sar_rook), logLik(sar_knn4)) )
res_sar <- residuals(sar_jabar) moran.test(res_sar, lw_queen, zero.policy = TRUE) # terpenuhi p > 0.05
shapiro.test(res_sar) # p > 0.05 ~ normal
library(lmtest)
res_sar <- residuals(sar_jabar)
x_vars <- sar_jabar$X[, -1, drop = FALSE] # buang intercept
df_bp <- data.frame(res_sar = res_sar, as.data.frame(x_vars))
bptest(res_sar ~ ., data = df_bp)
plot(fitted(sar_jabar), residuals(sar_jabar), xlab = “Fitted Values”, ylab = “Residuals”, main = “Uji Visual Homoskedastisitas (SAR)”) abline(h = 0, col = “red”)