Analisis Faktor Sosio-Ekonomik Indeks Pembangunan Manusia di Jawa Barat Menggunakan Model Spasial

1 DASHBOARD

Link : https://ghozzypurwana.shinyapps.io/Spatial_IPM_Dashboard/

2 ABSTRAK

Penelitian ini menganalisis pola spasial Indeks Pembangunan Manusia (IPM) di Provinsi Jawa Barat tahun 2024 serta pengaruh faktor ekonomi dan sosial menggunakan pendekatan ekonometrika spasial. Data sekunder dari BPS dan Open Data Jabar mencakup variabel IPM, PDRB, TPAK, UMR, Pengeluaran per Kapita, dan Gini Rasio. Hasil analisis menunjukkan adanya autokorelasi spasial signifikan, di mana wilayah dengan IPM tinggi cenderung berdekatan dengan wilayah bernilai tinggi lainnya. Model terbaik yaitu Spatial Durbin Model (SDM) menunjukkan bahwa PDRB berpengaruh positif dan signifikan terhadap IPM, sedangkan variabel lain tidak signifikan. Nilai rho negatif mengindikasikan adanya interaksi spasial negatif antarwilayah, sementara nilai AIC rendah menunjukkan model yang stabil dan sesuai. Temuan ini menegaskan pentingnya pertumbuhan ekonomi dan hubungan spasial antarwilayah dalam meningkatkan IPM secara merata di Jawa Barat.

3 PENDAHULUAN

3.1 Latar Belakang

Indeks Pembangunan Manusia (IPM) mencerminkan kemampuan penduduk dalam mengakses pendidikan, kesehatan, dan standar hidup layak. Di Jawa Barat, nilai IPM terus mengalami peningkatan dari tahun ke tahun, namun masih memperlihatkan kesenjangan antar kabupaten dan kota yang menandakan adanya ketimpangan kesejahteraan wilayah [1]. Faktor ekonomi dan sosial berperan penting dalam memengaruhi variasi IPM antarwilayah. Produk Domestik Regional Bruto (PDRB) menjadi indikator utama pertumbuhan ekonomi yang mencerminkan kapasitas daerah dalam menghasilkan barang dan jasa [2], sedangkan Tingkat Partisipasi Angkatan Kerja (TPAK) menunjukkan sejauh mana penduduk usia kerja berperan dalam kegiatan ekonomi. Upah Minimum Regional (UMR) juga berfungsi sebagai representasi kesejahteraan tenaga kerja formal dan kebijakan pemerintah daerah dalam menjamin pendapatan minimum [3][4].

Selain itu, pengeluaran per kapita menggambarkan daya beli dan standar hidup masyarakat, di mana nilai yang lebih tinggi umumnya mengindikasikan kondisi ekonomi yang lebih baik [6][7]. Gini Rasio, di sisi lain, mencerminkan tingkat ketimpangan distribusi pendapatan, sehingga semakin tinggi nilainya menunjukkan distribusi pendapatan yang semakin tidak merata. Karena perkembangan pembangunan sering kali membentuk pola spasial antarwilayah, pendekatan analisis spasial seperti Moran’s I, Local Indicators of Spatial Association (LISA), dan Getis-Ord Gi digunakan untuk mengidentifikasi klaster dan outlier pembangunan manusia di Jawa Barat [5][8].

3.2 Identifikasi Masalah

Penelitian ini mengkaji bagaimana pola sebaran IPM antar kabupaten dan kota di Jawa Barat tahun 2024, apakah terdapat keterkaitan spasial antarwilayah, serta bagaimana pengaruh PDRB, TPAK, UMR, pengeluaran per kapita, dan Gini Rasio terhadap variasi IPM.

3.3 Tujuan Penelitian

Penelitian ini bertujuan untuk memetakan pola spasial IPM di Jawa Barat, mengidentifikasi klaster pembangunan menggunakan uji Moran’s I, LISA, dan Getis-Ord Gi , serta menganalisis pengaruh PDRB, TPAK, UMR, pengeluaran per kapita, dan Gini Rasio terhadap IPM menggunakan perbandingan model OLS dan model spasial.

3.4 Batasan Penelitian

Penelitian dibatasi pada wilayah administratif kabupaten dan kota di Provinsi Jawa Barat dengan data sekunder tahun 2024 yang bersumber dari BPS dan portal data terbuka Pemerintah Provinsi Jawa Barat [1][2][3][4][5]. Variabel yang digunakan adalah IPM sebagai variabel dependen serta PDRB, TPAK, UMR, pengeluaran per kapita, dan Gini Rasio sebagai variabel independen. Analisis dilakukan menggunakan pendekatan Moran’s I, LISA, dan Getis-Ord Gi* untuk mendeteksi pola spasial pembangunan manusia.

4 TINJAUAN PUSTAKA

4.1 Teori Dasar Spatial Dependence

Spatial dependence menjelaskan bahwa nilai suatu variabel pada satu wilayah dapat dipengaruhi oleh nilai variabel yang sama pada wilayah tetangga, sehingga analisis yang mengabaikan hubungan spasial dapat menghasilkan estimasi yang bias. Pembobot spasial membentuk matriks hubungan antarwilayah, biasanya berdasarkan jarak, yang menjadi dasar perhitungan statistik dan model ekonometrika spasial [9].

4.2 Autokorelasi Spasial

Autokorelasi spasial mengukur apakah nilai variabel menunjukkan pola klaster positif, klaster negatif, atau acak. Pengujian dilakukan secara global untuk melihat pola keseluruhan dan secara lokal untuk mendeteksi klaster serta outlier, sehingga kombinasi ukuran global dan lokal memberikan gambaran lengkap pola spasial [10][11][12].

4.2.1 Moran’s I

Moran’s I digunakan untuk menilai ada tidaknya autokorelasi spasial secara global dengan membandingkan nilai suatu wilayah terhadap rata-rata keseluruhan dan kedekatan antarwilayah. Nilai positif menunjukkan klaster wilayah dengan nilai serupa, sedangkan nilai negatif menunjukkan pola berlawanan. Uji ini menentukan apakah pola spasial yang muncul bersifat signifikan atau acak [9][10].

4.2.2 Geary’s C

Geary’s C mengukur autokorelasi spasial dengan fokus pada perbedaan lokal antarwilayah. Nilai mendekati nol menunjukkan kemiripan tinggi, sedangkan nilai mendekati dua menunjukkan perbedaan besar. Uji ini melengkapi hasil Moran’s I untuk memperkuat kesimpulan tentang pola spasial yang ada [10][11].

4.2.3 Getis-Ord Gi

Getis-Ord Gi digunakan untuk mengidentifikasi area dengan konsentrasi nilai tinggi (hot spot) atau rendah (cold spot). Nilai positif menunjukkan wilayah dengan nilai tinggi yang berdekatan, sedangkan nilai negatif menunjukkan wilayah dengan nilai rendah yang saling berdekatan. Analisis ini membantu menemukan wilayah yang menonjol dalam konteks pembangunan [12].

4.2.4 Local Indicators of Spatial Association (LISA)

Getis-Ord Gi digunakan untuk mengidentifikasi area dengan konsentrasi nilai tinggi (hot spot) atau rendah (cold spot). Nilai positif menunjukkan wilayah dengan nilai tinggi yang berdekatan, sedangkan nilai negatif menunjukkan wilayah dengan nilai rendah yang saling berdekatan. Analisis ini membantu menemukan wilayah yang menonjol dalam konteks pembangunan [12].

4.3 Model Spatial Econometrics

Model ekonometrika spasial dikembangkan untuk menangani ketergantungan spasial pada data lintas wilayah. Ketika hubungan antarwilayah diabaikan, model konvensional seperti OLS bisa memberikan hasil estimasi yang bias. Oleh karena itu, berbagai model spasial digunakan untuk menggambarkan sumber pengaruh antarwilayah secara lebih akurat [9][10].

4.3.1 Ordinary Least Squares (OLS)

Model OLS digunakan sebagai titik awal dalam analisis spasial untuk menguji hubungan antara variabel dependen dan independen tanpa mempertimbangkan pengaruh spasial. Hasil model ini penting untuk melihat apakah terdapat autokorelasi pada residual yang menjadi dasar untuk memilih model spasial yang lebih sesuai [9].\[Y\ = X\beta\ + \ \varepsilon\ \]

dengan keterangan:

Y = vektor variabel dependen

X= matriks variabel independen

β = vektor koefisien regresi

ε = error

4.3.2 Spatial Lag of X (SLX)

Model SLX menambahkan lag spasial dari variabel independen ke dalam model regresi. Artinya, nilai variabel bebas pada wilayah tetangga ikut memengaruhi nilai variabel dependen pada wilayah tertentu. Model ini berguna untuk melihat seberapa jauh pengaruh eksternal suatu wilayah terhadap wilayah lain tanpa memasukkan efek spasial pada variabel dependen [10].

\(Y\ = \ X\beta\ + \ W\ X\theta\ + \ \varepsilon\)

Keterangan:

W = matriks pembobot spasial

W X = variabel independen yang dipengaruhi oleh tetangga (spatially lagged X)

4.3.3 Spatial Lag Model (SAR)

Model SAR atau spatial lag model memasukkan lag spasial dari variabel dependen, yang berarti nilai Y di suatu wilayah dipengaruhi oleh nilai Y di wilayah tetangga. Koefisien lag (ρ) menunjukkan besarnya penyebaran antarwilayah. Model ini cocok digunakan bila interaksi antarwilayah terjadi melalui fenomena yang menular atau saling memengaruhi secara langsung seperti pertumbuhan ekonomi atau IPM [9].

\(Y\ = \ \rho WY\ + \ X\beta\ + \ \varepsilon\)

4.3.4 Spatial Error Model (SEM)

Model SEM digunakan ketika efek spasial muncul melalui komponen error. Dalam model ini, korelasi antarwilayah terjadi karena faktor yang tidak teramati, bukan dari hubungan langsung antarvariabel utama. SEM membantu memperbaiki kesalahan model OLS ketika autokorelasi spasial muncul pada error residual [10].

\(Y\ = \ X\beta\ + \ \mu\)

\(\mu\ = \ \lambda W\mu\ + \ \varepsilon\)

4.3.5 Spatial Durbin Model (SDM)

SDM merupakan pengembangan dari SAR dengan menambahkan lag spasial dari variabel independen. Model ini memungkinkan efek tidak hanya menyebar melalui variabel dependen, tetapi juga melalui faktor-faktor penjelas. SDM sering dianggap sebagai model yang paling fleksibel karena dapat menggambarkan berbagai bentuk keterkaitan spasial antarwilayah [9].

\(Y\ = \ \rho WY\ + \ X\beta\ + W\ X\theta\ + \ \varepsilon\)

4.3.6 Spatial Durbin Error Model (SDEM)

Model SDEM adalah pengembangan dari SDM yang menggabungkan efek lag pada variabel independen serta autokorelasi pada komponen error. Model ini digunakan jika keterkaitan spasial berasal dari gabungan pengaruh langsung antarwilayah dan kesalahan yang saling berkaitan [10].

\(Y\ = \ X\beta\ + W\ X\theta\ + \mu\)

5 METODOLOGI PENELITIAN

5.1 DATA

Penelitian ini menggunakan data sekunder yang diperoleh dari sumber resmi seperti Badan Pusat Statistik (BPS) dan portal Open Data Jabar. Seluruh variabel disusun dalam satu tabel dan diintegrasikan dengan data spasial tingkat kabupaten/kota di Provinsi Jawa Barat. Data mencakup satu variabel dependen dan lima variabel independen yang relevan dengan pembangunan manusia di tingkat daerah.

Sumber Data:

X1: PDRB (Produk Domestik Regional Bruto)

X2: TPAK (Tingkat Partisipasi Angkatan Kerja)

X3: UMR (Upah Minimum Regional)

X4: Pengeluaran per Kapita

X5: Gini Rasio

Y: IPM (Indeks Pembangunan Manusia)

5.1.1 Variabel Dependen

IIndeks Pembangunan Manusia (IPM) merupakan variabel utama yang menjadi fokus penelitian. IPM digunakan untuk menilai tingkat keberhasilan pembangunan manusia di suatu daerah dengan menggabungkan tiga dimensi utama, yaitu pendidikan, kesehatan, dan standar hidup. Nilai IPM yang lebih tinggi mencerminkan tingkat kesejahteraan masyarakat yang lebih baik dan kemajuan pembangunan di wilayah tersebut.

5.1.2 Variabel Independen

Penelitian ini menggunakan lima variabel independen yang diduga memiliki pengaruh terhadap IPM, yaitu:

1. Produk Domestik Regional Bruto (PDRB): menunjukkan tingkat aktivitas ekonomi daerah. Semakin tinggi nilai PDRB, semakin besar kemampuan daerah dalam menyediakan kesempatan kerja dan meningkatkan kesejahteraan penduduk.

2. Tingkat Partisipasi Angkatan Kerja (TPAK): menggambarkan proporsi penduduk usia kerja yang aktif secara ekonomi. Nilai TPAK yang tinggi menunjukkan keterlibatan tenaga kerja yang besar dalam kegiatan produktif.

3. Upah Minimum Regional (UMR): mencerminkan tingkat kesejahteraan tenaga kerja formal dan kebijakan pemerintah daerah dalam menjamin pendapatan minimum.

4. Pengeluaran per Kapita: mengindikasikan daya beli masyarakat dan standar hidup rata-rata penduduk. Pengeluaran yang lebih tinggi umumnya mencerminkan kondisi ekonomi yang lebih baik.

5. Gini Rasio: mengukur tingkat ketimpangan pendapatan. Nilai Gini Rasio yang lebih tinggi menunjukkan distribusi pendapatan yang semakin timpang dan potensi ketidakmerataan kesejahteraan antarwilayah

5.2 UNIT SPASIAL

Unit analisis dalam penelitian ini adalah 27 kabupaten dan kota di Provinsi Jawa Barat. Setiap unit wilayah direpresentasikan sebagai satu poligon yang memiliki batas administratif resmi.

Pembentukan Bobot Spasial

Bobot spasial (spatial weights) digunakan untuk mendefinisikan hubungan kedekatan antarwilayah. Bobot ini disusun menggunakan metode Queen Contiguity, yaitu dua wilayah dianggap bertetangga bila berbagi sisi maupun titik sudut peta.
Dalam penelitian ini dibuat dua tipe bobot:

• lwW: bobot row-standardized, digunakan pada analisis autokorelasi dan sebagian besar model spasial.

• lwB: bobot binary (0/1), digunakan pada model lain yang tidak memerlukan normalisasi bobot.

5.3 METODE ANALISIS

Analisis dilakukan melalui beberapa tahap mulai dari eksplorasi data, pengujian autokorelasi spasial, hingga pemodelan ekonometrika spasial.

5.3.1 Analisis Deskriptif

Langkah pertama adalah melihat gambaran umum variabel. Statistik deskriptif seperti nilai rata-rata, minimum, maksimum, dan standar deviasi dihitung menggunakan fungsi pada R.Selain itu, peta tematik dibuat dengan untuk menampilkan distribusi spasial masing-masing variabe

5.3.2 Uji Autokorelasi Spasial Global

5.3.2.1 Moran’s I

Uji Moran’s I digunakan untuk mendeteksi apakah nilai antarwilayah memiliki hubungan spasial. Rumus umum Moran’s I adalah:

\(I\ = \ \frac{n}{S_{0}}\ \cdot \ \frac{\sum_{i}^{}\sum_{j}^{}\ w_{ij}^{}\ (x_{i}\ - \ \underline{x}\ )\ (x_{j}\ - \ \ \underline{x})}{\sum_{i}^{}\ {(x_{i}\ - \ \underline{x}\ )}^{2}}\)

Keterangan :

n = Jumlah Wilayah

\(w_{ij}^{}\) = bobot spasial antara wilayah i dan j

\(x_{i}\) ​= nilai variabel pada wilayah iii

\(\ \underline{x}\) = rata-rata variabel

\(S_{0}\) = \(\sum_{i}^{}\sum_{j}^{}\ w_{ij}^{}\)

Interpretasinya:

• I > 0 = terdapat klaster wilayah dengan nilai serupa (autokorelasi positif).

• I < 0 = wilayah dengan nilai tinggi cenderung berdekatan dengan wilayah bernilai rendah (autokorelasi negatif).

• I sama dengan 0 = pola acak tanpa pengaruh spasial.

5.3.2.2 Geary’s C

Uji Geary’s C juga digunakan untuk mendeteksi autokorelasi spasial, namun lebih sensitif terhadap variasi lokal antarwilayah. Rumusnya adalah:

\(C\ = \ \frac{(n\ - \ 1)\ \sum_{i}^{}\sum_{j}^{}\ w_{ij}^{}\ {(x_{i}\ - \ \underline{x}\ )}^{2}}{2\ \ S_{0}\ \ \sum_{i}^{}\ {(x_{i}\ - \ \underline{x}\ )}^{2}}\)

Keterangan :

n = Jumlah Wilayah

\(w_{ij}^{}\) = bobot spasial antara wilayah i dan j

\(x_{i}\) ​= nilai variabel pada wilayah iii

\(\ \underline{x}\) = rata-rata variabel

\(S_{0}\) = \(\sum_{i}^{}\sum_{j}^{}\ w_{ij}^{}\)

Interpretasinya:

• C < 1 = autokorelasi positif (nilai mirip antarwilayah).

• C > 1 = autokorelasi negatif (perbedaan besar antarwilayah).

• C sama dengan 1 = tidak ada pola spasial.

5.3.3 Analisis Autokorelasi Lokal

Untuk melihat klaster atau outlier lokal digunakan metode LISA dan Getis-Ord Gi.

LISA (Local Moran’s I) = mengidentifikasi wilayah dengan pola High-High, Low-Low, High-Low, atau Low-High.

Getis-Ord Gi = mendeteksi area dengan konsentrasi nilai tinggi (hot spot) dan rendah (cold spot).

5.3.4 Pemodelan Regresi Spasial

5.3.4.1 Model OLS

Model OLS digunakan sebagai langkah awal untuk mengetahui hubungan dasar antara IPM dan variabel-variabel independen. Setelah model OLS dijalankan, residualnya diuji menggunakan Moran’s I untuk memastikan apakah terdapat pengaruh spasial yang belum tertangkap model.

5.3.4.2 Model Spasial

Jika ditemukan autokorelasi pada residual OLS, digunakan model-model spasial (SAR, SEM, SLX, SDM, SDEM, SAC, GNC), Perbandingan antar model dilakukan dengan melihat nilai AIC (Akaike Information Criterion) dan Log Likelihood. Model dengan AIC terendah dianggap paling sesuai dengan data.

5.4 ALUR KERJA PENELITIAN

Alur kerja penelitian ini dimulai dengan pengumpulan data IPM, PDRB, Indeks Pengetahuan, UMR, dan Pengeluaran per Kapita dari sumber terbuka, kemudian data tersebut diintegrasikan dengan peta batas administratif Jawa Barat. Setelah itu, eksplorasi melalui peta sebaran dan statistik deskriptif, dilanjutkan dengan uji autokorelasi menggunakan Moran’s I, Geary’s C, serta Moran Scatterplot untuk mengetahui pola spasial global. Analisis lokal dilakukan dengan metode LISA dan Getis-Ord Gi untuk mengidentifikasi klaster dan area hot spot. Selanjutnya dilakukan pemodelan spasial menggunakan OLS, SAR, SEM, dan SDM untuk melihat hubungan antarvariabel, kemudian hasilnya dievaluasi dengan membandingkan nilai AIC dan log likelihood. Tahap akhir yaitu visualisasi hasil dalam bentuk dashboard interaktif sehingga analisis spasial dapat disajikan dengan lebih mudah dan informatif.

6 HASIL DAN PEMBAHASAN

6.1 PETA DESKRIPTIF

Mengimport Library dan pengInputan Data

library(spdep)

## Warning: package 'spdep' was built under R version 4.4.3

## Loading required package: spData

## Warning: package 'spData' was built under R version 4.4.3

## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`

## Loading required package: sf

## Warning: package 'sf' was built under R version 4.4.3

## Linking to GEOS 3.13.0, GDAL 3.10.1, PROJ 9.5.1; sf_use_s2() is TRUE

library(sp)

## Warning: package 'sp' was built under R version 4.4.3

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

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(spData)
library(geodata)

## Warning: package 'geodata' was built under R version 4.4.3

## Loading required package: terra

## Warning: package 'terra' was built under R version 4.4.3

## terra 1.8.60

library(tmap)

## Warning: package 'tmap' was built under R version 4.4.3

library(skimr)

## Warning: package 'skimr' was built under R version 4.4.3

library(psych)

## Warning: package 'psych' was built under R version 4.4.3

## 
## Attaching package: 'psych'

## The following objects are masked from 'package:terra':
## 
##     describe, distance, rescale

## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha

library(gridExtra)

## Warning: package 'gridExtra' was built under R version 4.4.3

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

library(grid)

## 
## Attaching package: 'grid'

## The following object is masked from 'package:terra':
## 
##     depth

library(viridis) 

## Warning: package 'viridis' was built under R version 4.4.3

## Loading required package: viridisLite

library(spatialreg)

## Warning: package 'spatialreg' was built under R version 4.4.3

## Loading required package: Matrix

## 
## Attaching package: 'spatialreg'

## The following objects are masked from 'package:spdep':
## 
##     get.ClusterOption, get.coresOption, get.mcOption,
##     get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
##     set.coresOption, set.mcOption, set.VerboseOption,
##     set.ZeroPolicyOption

setwd("C:\\Users\\Ghozzy\\OneDrive\\Documents\\SMESTER 5\\Sparsial Stat")
indo_adm2 <- st_read("IDN_AdminBoundaries_candidate.gdb",
                     layer = "idn_admbnda_adm2_bps_20200401")

## Reading layer `idn_admbnda_adm2_bps_20200401' from data source 
##   `C:\Users\Ghozzy\OneDrive\Documents\SMESTER 5\Sparsial Stat\IDN_AdminBoundaries_candidate.gdb' 
##   using driver `OpenFileGDB'

## Warning in CPL_read_ogr(dsn, layer, query, as.character(options), quiet, : GDAL
## Message 1: organizePolygons() received a polygon with more than 100 parts.  The
## processing may be really slow.  You can skip the processing by setting
## METHOD=SKIP.

## Simple feature collection with 522 features and 14 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 95.01079 ymin: -11.00762 xmax: 141.0194 ymax: 6.07693
## Geodetic CRS:  WGS 84

jabar <- indo_adm2[indo_adm2$admin1Name_en == "Jawa Barat", ]
jabar <- jabar[jabar$admin2Name_en != "Waduk Cirata", ]
plot(st_geometry(jabar))

data <- na.omit(read.csv("spatial_data_uts.csv", sep=";", dec="."))
jabar$id <- c(1:27)
jabar_sf <- st_as_sf(jabar)
data$id <- c(1:27)
jabar_merged <- jabar_sf %>%
  left_join(data, by = "id")

data_desc <- jabar_merged %>%
  st_drop_geometry() %>%
  select(
    Y..IPM.,
    X1..PDRB.,
    X2..TPAK.,
    X3..UMR.,
    X4..Pengeluaran.per.kapita.,
    X5..Gini.Rasio.
  )

describe(data_desc)

##                             vars  n       mean         sd     median    trimmed
## Y..IPM.                        1 27      74.68       4.33      73.82      74.36
## X1..PDRB.                      2 27  105053.23  112708.74   61236.34   87923.66
## X2..TPAK.                      3 27      68.06       3.44      67.93      67.80
## X3..UMR.                       4 27 3370534.44 1106474.44 3294485.00 3315578.61
## X4..Pengeluaran.per.kapita.    5 27   11928.85    2389.91   11563.00   11674.70
## X5..Gini.Rasio.                6 27       0.37       0.05       0.37       0.37
##                                    mad        min        max      range skew
## Y..IPM.                           4.42      68.89      83.75      14.86 0.73
## X1..PDRB.                     40057.20    5549.26  421323.50  415774.24 1.64
## X2..TPAK.                         2.22      62.58      80.12      17.54 1.35
## X3..UMR.                    1356318.06 2070192.00 5343430.00 3273238.00 0.45
## X4..Pengeluaran.per.kapita.    2044.51    8965.00   18795.00    9830.00 1.17
## X5..Gini.Rasio.                   0.05       0.29       0.48       0.19 0.28
##                             kurtosis        se
## Y..IPM.                        -0.59      0.83
## X1..PDRB.                       1.37  21690.81
## X2..TPAK.                       3.34      0.66
## X3..UMR.                       -1.27 212941.11
## X4..Pengeluaran.per.kapita.     1.06    459.94
## X5..Gini.Rasio.                -0.68      0.01

# --- 1. spatial weights (queen) ---
nb <- poly2nb(as_Spatial(jabar_merged), queen = TRUE)
lwW <- nb2listw(nb, style = "W", zero.policy = TRUE)
lwB <- nb2listw(nb, style = "B", zero.policy = TRUE)

summary(nb)

## Neighbour list object:
## Number of regions: 27 
## Number of nonzero links: 106 
## Percentage nonzero weights: 14.54047 
## Average number of links: 3.925926 
## Link number distribution:
## 
## 1 2 3 4 5 6 7 8 
## 4 3 6 4 1 7 1 1 
## 4 least connected regions:
## 12 14 16 18 with 1 link
## 1 most connected region:
## 4 with 8 links

# Jumlah wilayah tanpa tetangga (isolated areas)

# ================================
#  FUNCTION (tidak diubah)
# ================================
analyze_spatial <- function(data_sf, var_col, label) {
  cat("\n\n==============================\n", label, "\n==============================\n")
  
  x_raw <- data_sf[[var_col]]
  
  # ---------- Global tests ----------
  moran_res <- moran.test(x_raw, lwW, randomisation = TRUE, alternative = "two.sided")
  geary_res <- geary.test(x_raw, lwW, randomisation = TRUE, alternative = "two.sided")
  print(moran_res)
  print(geary_res)
  
  # ---------- Local Moran ----------
  x <- scale(x_raw)[, 1]
  lagx <- lag.listw(lwW, x)
  lisa <- localmoran(x, lwW, alternative = "two.sided", zero.policy = TRUE)
  lisa_df <- as.data.frame(lisa)
  names(lisa_df) <- c("Ii","Ei","Vi","Zi","Pi.two.sided")
  
  alpha <- 0.05
  quad <- dplyr::case_when(
    x >= 0 & lagx >= 0 ~ "High-High",
    x < 0 & lagx < 0 ~ "Low-Low",
    x >= 0 & lagx < 0 ~ "High-Low (Outlier)",
    x < 0 & lagx >= 0 ~ "Low-High (Outlier)"
  )
  
  LISA_sf <- bind_cols(data_sf, lisa_df) |>
    mutate(quad = ifelse(Pi.two.sided <= alpha, quad, "Not significant"))
  
  # ---------- Local Geary’s C ----------
  localC_vals <- localC(x_raw, lwW)
  Geary_sf <- mutate(data_sf,
                     LocalC = as.numeric(localC_vals))
  
  # ---------- Getis–Ord Gi* ----------
  sum_x <- sum(x_raw)
  Wb <- listw2mat(lwB)
  Wb_star <- Wb; diag(Wb_star) <- 1
  G_star_raw <- as.numeric(Wb_star %*% x_raw) / sum_x
  Gz <- localG(x_raw, listw = lwW)
  G_sf <- mutate(data_sf,
                 z_Gistar = as.numeric(Gz),
                 hotcold = case_when(
                   z_Gistar >= 1.96 ~ "Hot spot (p<0.05)",
                   z_Gistar <= -1.96 ~ "Cold spot (p<0.05)",
                   TRUE ~ "Not significant"
                 ))
  
  # ---------- MAPS ----------
  p1 <- ggplot(LISA_sf) +
    geom_sf(aes(fill = quad), color = "white", size = 0.2) +
    scale_fill_manual(values = c(
      "High-High" = "#d73027",
      "Low-Low" = "#4575b4",
      "High-Low (Outlier)" = "#fdae61",
      "Low-High (Outlier)" = "#74add1",
      "Not significant" = "grey85"
    )) +
    labs(title = paste("Local Moran's I (LISA) —", label), fill = "Kategori") +
    theme_minimal()
  
  p2 <- ggplot(Geary_sf) +
    geom_sf(aes(fill = LocalC), color = "white", size = 0.2) +
    scale_fill_viridis_c(option = "magma", direction = -1,
                         name = "Local Geary’s C") +
    labs(title = paste("Local Geary’s C —", label)) +
    theme_minimal()
  
  p3 <- ggplot(G_sf) +
    geom_sf(aes(fill = hotcold), color = "white", size = 0.2) +
    scale_fill_manual(values = c(
      "Hot spot (p<0.05)" = "#b2182b",
      "Cold spot (p<0.05)" = "#2166ac",
      "Not significant" = "grey85"
    )) +
    labs(title = paste("Getis–Ord Gi* —", label), fill = NULL) +
    theme_minimal()
  
  print(p1); print(p2); print(p3)
}

6.1.1 Indeks Pembangunan Manusia (IPM)

## Peta IPM
ggplot(data = jabar_merged) +
  geom_sf(aes(fill = Y..IPM.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "plasma", direction = -1,
                       name = "Indeks Pembangunan Manusia (IPM)") +
  labs(
    title = "Peta IPM Jawa Barat 2024",
    subtitle = "Nilai Indeks Pembangunan Manusia per Kabupaten/Kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 14) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Interpretasi:

Peta tersebut menampilkan sebaran nilai IPM di Jawa Barat berdasarkan kabupaten/kota. Warna gelap (ungu–biru tua) menunjukkan IPM tinggi, sedangkan warna terang (kuning–oranye) menunjukkan IPM lebih rendah.Terlihat adanya pola spasial, di mana wilayah dengan IPM tinggi cenderung berkelompok di bagian utara dan sekitar kota besar seperti Bandung dan Bekasi, sedangkan wilayah selatan menunjukkan IPM lebih rendah.

analyze_spatial(jabar_merged, "Y..IPM.", "Indeks Pembangunan Manusia (IPM)")

## 
## 
## ==============================
##  Indeks Pembangunan Manusia (IPM) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 2.236, p-value = 0.02535
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.27459101       -0.03846154        0.01960115 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 2.415, p-value = 0.01573
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##         0.6069834         1.0000000         0.0264836

• Local Moran’s I
  Terlihat satu wilayah (sekitar bagian tengah Jawa Barat) tergolong Low–Low (Wilayah yang memiliki nilai IPM rendah dan dikelilingi oleh wilayah lain yang juga memiliki IPM rendah)

• Local Geary’s C

Daerah dengan warna gelap (ungu ke hitam) berarti perbedaan IPM dengan tetangganya tinggi, sedangkan warna kuning berarti lebih homogen.

• Getis-Ord Gi

Ada satu wilayah berwarna biru cold spot (Wilayahyang memiliki nilai IPM rendah di antara tetangga-tetangganya.)

6.1.2 Produk Domestik Regional Bruto (PDRB)

ggplot(data = jabar_merged) +
  geom_sf(aes(fill = X1..PDRB.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "magma", direction = -1,
                       name = "Produk Domestik Regional Bruto (PDRB)") +
  labs(
    title = "Sebaran PDRB di Jawa Barat",
    subtitle = "Nilai PDRB per Kabupaten/Kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 14) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Interpretasi:

Daerah berwarna gelap (ungu ke hitam) seperti Kota Bandung dan wilayah utara (Bekasi, Karawang, dan sekitarnya) memiliki PDRB tinggi, menandakan aktivitas ekonomi dan industri yang kuat. Sementara itu, wilayah dengan warna terang (kuning–oranye) di bagian selatan dan timur Jawa Barat memiliki PDRB lebih rendah, yang mencerminkan aktivitas ekonomi yang masih terbatas.

analyze_spatial(jabar_merged, "X1..PDRB.", "Produk Domestik Regional Bruto (PDRB)")

## 
## 
## ==============================
##  Produk Domestik Regional Bruto (PDRB) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 2.4963, p-value = 0.01255
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.29479733       -0.03846154        0.01782296 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 1.6091, p-value = 0.1076
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.69552243        1.00000000        0.03580473

• Local Moran’s I

Terdapat klaster High-High di wilayah utara Jawa Barat, menunjukkan daerah dengan PDRB tinggi saling berdekatan.

• Local Geary’s C

Nilai tinggi di barat–utara menandakan autokorelasi spasial kuat antarwilayah dengan PDRB serupa.

• Getis-Ord Gi

Terbentuk hot spot di wilayah utara, menegaskan konsentrasi aktivitas ekonomi tinggi di area tersebut.

6.1.3 Tingkat Partisipasi Angkatan Kerja (TPAK)

ggplot(data = jabar_merged) +
  geom_sf(aes(fill = X2..TPAK.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "inferno", direction = -1,
                       name = "Tingkat Partisipasi Angkatan Kerja (TPAK)") +
  labs(
    title = "Sebaran TPAK di Jawa Barat",
    subtitle = "Persentase TPAK per Kabupaten/Kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 9) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Peta tersebut menunjukkan sebaran Tingkat Partisipasi Angkatan Kerja (TPAK) di Provinsi Jawa Barat tahun 2024. Daerah dengan warna gelap, terutama di bagian selatan, memiliki tingkat partisipasi angkatan kerja yang lebih tinggi, sedangkan wilayah berwarna kuning atau oranye di bagian utara dan tengah menunjukkan partisipasi yang lebih rendah. Pola ini mengindikasikan adanya variasi keterlibatan penduduk usia kerja di berbagai daerah, kemungkinan dipengaruhi oleh faktor ekonomi, kesempatan kerja, serta karakteristik wilayah perkotaan dan pedesaan.

analyze_spatial(jabar_merged, "X2..TPAK.", "Tingkat Partisipasi Angkatan Kerja (TPAK)")

## 
## 
## ==============================
##  Tingkat Partisipasi Angkatan Kerja (TPAK) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 0.98267, p-value = 0.3258
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.08598213       -0.03846154        0.01603727 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 1.5182, p-value = 0.129
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.67735304        1.00000000        0.04516516

• Local Moran’s I (LISA)
  Terdapat klaster High-High di bagian selatan Jawa Barat, menunjukkan wilayah dengan tingkat partisipasi angkatan kerja tinggi saling berdekatan di area tersebut.

• Local Geary’s C
  Nilai Local Geary’s C yang tinggi di wilayah barat dan timur menunjukkan adanya keseragaman nilai TPAK antarwilayah di kedua area tersebut.

• Getis–Ord Gi*
  Terbentuk hot spot signifikan di bagian selatan Jawa Barat, menegaskan adanya konsentrasi wilayah dengan tingkat partisipasi angkatan kerja tinggi di daerah tersebut.

6.1.4 Upah Minimum Rakyat (UMR)

ggplot(data = jabar_merged) +
  geom_sf(aes(fill = X3..UMR.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "cividis", direction = -1,
                       name = "Upah Minimum Rakyat (UMR)") +
  labs(
    title = "Sebaran UMR di Jawa Barat",
    subtitle = "Nilai UMR per Kabupaten/Kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 9) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Interpretasi:

Bisa dilihat dari Visualisasinya bahwa Peta tersebut menunjukkan sebaran Upah Minimum Rakyat (UMR) di Provinsi Jawa Barat tahun 2024. Wilayah dengan warna biru tua memiliki nilai UMR tertinggi, seperti di bagian barat dan selatan provinsi, sedangkan daerah dengan warna kuning menunjukkan UMR yang lebih rendah, terutama di wilayah timur dan tengah. Pola ini menggambarkan adanya ketimpangan upah antar daerah, di mana kawasan industri dan perkotaan cenderung memiliki upah minimum yang lebih tinggi dibandingkan wilayah pedesaan.

analyze_spatial(jabar_merged, "X3..UMR.", "Upah Minimum Rakyat (UMR)")

## 
## 
## ==============================
##  Upah Minimum Rakyat (UMR) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 0.037084, p-value = 0.9704
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.03318861       -0.03846154        0.02021745 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.051232, p-value = 0.9591
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.00781240        1.00000000        0.02325299

• Local Moran’s I

menunjukkan adanya kluster High-High di wilayah barat (sekitar Bekasi–Karawang), artinya daerah dengan UMR tinggi cenderung berdekatan dengan daerah ber-UMR tinggi lainnya.

• Local Geary’s C

memperlihatkan variasi spasial cukup besar, menandakan adanya perbedaan tajam antarwilayah dalam nilai UMR.

• Getis-Ord Gi

mengidentifikasi hot spot di barat (wilayah industri dengan UMR tinggi) dan cold spot di bagian tengah (wilayah dengan UMR relatif rendah).

6.1.5 Pengeluaran Per Kapita

ggplot(data = jabar_merged) +
  geom_sf(aes(fill = X4..Pengeluaran.per.kapita.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "magma", direction = -1,
                       name = "Pengeluaran per Kapita") +
  labs(
    title = "Sebaran Pengeluaran per Kapita di Jawa Barat",
    subtitle = "Nilai pengeluaran per kapita per kabupaten/kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 8) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Interpretasi:

Peta tersebut memperlihatkan sebaran pengeluaran per kapita di Provinsi Jawa Barat tahun 2024. Wilayah dengan warna ungu tua memiliki tingkat pengeluaran per kapita yang tinggi, terutama di bagian barat dan sebagian tengah provinsi, sedangkan daerah dengan warna kuning menunjukkan pengeluaran per kapita yang lebih rendah, terutama di bagian selatan dan timur. Pola ini menunjukkan bahwa tingkat kesejahteraan masyarakat cenderung lebih tinggi di kawasan perkotaan atau ekonomi maju dibandingkan wilayah yang lebih rural.

analyze_spatial(jabar_merged, "X4..Pengeluaran.per.kapita.", "Pengeluaran per Kapita")

## 
## 
## ==============================
##  Pengeluaran per Kapita 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = -1.4182, p-value = 0.1561
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.22929156       -0.03846154        0.01810562 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.57649, p-value = 0.5643
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.10680275        1.00000000        0.03432302

• Local Moran’s I

menunjukkan adanya kluster Low-High di bagian tenggara, menandakan wilayah dengan pengeluaran rendah yang berdekatan dengan daerah berpengeluaran tinggi.

• Local Geary’s C

memperlihatkan variasi spasial yang cukup besar antarwilayah, menunjukkan adanya ketimpangan pengeluaran per kapita di Jawa Barat.

• Getis-Ord Gi

mengidentifikasi hot spot di bagian tenggara, yang berarti konsentrasi wilayah dengan pengeluaran per kapita tinggi berada di area tersebut.

6.1.6 Gini Rasio

ggplot(data = jabar_merged) +
  geom_sf(aes(fill = X5..Gini.Rasio.), color = "white", linewidth = 0.3) +
  scale_fill_viridis_c(option = "turbo", direction = -1,
                       name = "Gini Rasio") +
  labs(
    title = "Sebaran Gini Rasio di Jawa Barat",
    subtitle = "Ketimpangan pendapatan per kabupaten/kota",
    caption = "Sumber: BPS (2024)"
  ) +
  theme_minimal(base_size = 9) +
  theme(plot.title = element_text(face = "bold", hjust = 0.5))

Interpretasi:

Peta tersebut menunjukkan sebaran Gini Rasio di Provinsi Jawa Barat tahun 2024, yang menggambarkan tingkat ketimpangan pendapatan antar kabupaten/kota. Wilayah dengan warna biru tua memiliki nilai Gini Rasio lebih tinggi (ketimpangan pendapatan lebih besar), sedangkan daerah berwarna hijau atau kuning menunjukkan ketimpangan yang lebih rendah. Secara umum, ketimpangan pendapatan tampak bervariasi antar wilayah, dengan beberapa daerah di bagian selatan dan timur mengalami ketimpangan yang lebih tinggi dibandingkan wilayah tengah dan barat yang cenderung lebih merata.

analyze_spatial(jabar_merged, "X5..Gini.Rasio.", "Gini Rasio")

## 
## 
## ==============================
##  Gini Rasio 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = -0.8554, p-value = 0.3923
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.15846484       -0.03846154        0.01968105 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.23704, p-value = 0.8126
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.03826887        1.00000000        0.02606476

• Local Moran’s I (LISA)
  Terdapat klaster Low-High (Outlier) di bagian utara dan selatan Jawa Barat, menandakan wilayah dengan ketimpangan rendah berdekatan dengan daerah yang memiliki ketimpangan tinggi.

• Local Geary’s C
  Nilai tinggi pada wilayah timur menunjukkan adanya keseragaman nilai Gini Rasio di area tersebut, sedangkan nilai lebih rendah di bagian barat menandakan variasi antarwilayah yang lebih besar.

• Getis–Ord Gi*
  Terbentuk hot spot signifikan di wilayah utara dan selatan Jawa Barat, menegaskan konsentrasi daerah dengan tingkat ketimpangan tinggi di kedua kawasan tersebut.

6.2 HASIL UJI AUTOKORELASI

# --- 1. spatial weights (queen) ---
nb <- poly2nb(as_Spatial(jabar_merged), queen = TRUE)
lwW <- nb2listw(nb, style = "W", zero.policy = TRUE)
lwB <- nb2listw(nb, style = "B", zero.policy = TRUE)

summary(nb)

## Neighbour list object:
## Number of regions: 27 
## Number of nonzero links: 106 
## Percentage nonzero weights: 14.54047 
## Average number of links: 3.925926 
## Link number distribution:
## 
## 1 2 3 4 5 6 7 8 
## 4 3 6 4 1 7 1 1 
## 4 least connected regions:
## 12 14 16 18 with 1 link
## 1 most connected region:
## 4 with 8 links

# Jumlah wilayah tanpa tetangga (isolated areas)

# ================================
#  FUNCTION (tidak diubah)
# ================================
analyze_spatial <- function(data_sf, var_col, label) {
  cat("\n\n==============================\n", label, "\n==============================\n")
  
  x_raw <- data_sf[[var_col]]
  
  # ---------- Global tests ----------
  moran_res <- moran.test(x_raw, lwW, randomisation = TRUE, alternative = "two.sided")
  geary_res <- geary.test(x_raw, lwW, randomisation = TRUE, alternative = "two.sided")
  print(moran_res)
  print(geary_res)
  
  # ---------- Local Moran ----------
  x <- scale(x_raw)[, 1]
  lagx <- lag.listw(lwW, x)
  lisa <- localmoran(x, lwW, alternative = "two.sided", zero.policy = TRUE)
  lisa_df <- as.data.frame(lisa)
  names(lisa_df) <- c("Ii","Ei","Vi","Zi","Pi.two.sided")
  
  alpha <- 0.05
  quad <- dplyr::case_when(
    x >= 0 & lagx >= 0 ~ "High-High",
    x < 0 & lagx < 0 ~ "Low-Low",
    x >= 0 & lagx < 0 ~ "High-Low (Outlier)",
    x < 0 & lagx >= 0 ~ "Low-High (Outlier)"
  )
  
  LISA_sf <- bind_cols(data_sf, lisa_df) |>
    mutate(quad = ifelse(Pi.two.sided <= alpha, quad, "Not significant"))
  
  # ---------- Local Geary’s C ----------
  localC_vals <- localC(x_raw, lwW)
  Geary_sf <- mutate(data_sf,
                     LocalC = as.numeric(localC_vals))
  
  # ---------- Getis–Ord Gi* ----------
  sum_x <- sum(x_raw)
  Wb <- listw2mat(lwB)
  Wb_star <- Wb; diag(Wb_star) <- 1
  G_star_raw <- as.numeric(Wb_star %*% x_raw) / sum_x
  Gz <- localG(x_raw, listw = lwW)
  G_sf <- mutate(data_sf,
                 z_Gistar = as.numeric(Gz),
                 hotcold = case_when(
                   z_Gistar >= 1.96 ~ "Hot spot (p<0.05)",
                   z_Gistar <= -1.96 ~ "Cold spot (p<0.05)",
                   TRUE ~ "Not significant"
                 ))
  
  # ---------- MAPS ----------
  p1 <- ggplot(LISA_sf) +
    geom_sf(aes(fill = quad), color = "white", size = 0.2) +
    scale_fill_manual(values = c(
      "High-High" = "#d73027",
      "Low-Low" = "#4575b4",
      "High-Low (Outlier)" = "#fdae61",
      "Low-High (Outlier)" = "#74add1",
      "Not significant" = "grey85"
    )) +
    labs(title = paste("Local Moran's I (LISA) —", label), fill = "Kategori") +
    theme_minimal()
  
  p2 <- ggplot(Geary_sf) +
    geom_sf(aes(fill = LocalC), color = "white", size = 0.2) +
    scale_fill_viridis_c(option = "magma", direction = -1,
                         name = "Local Geary’s C") +
    labs(title = paste("Local Geary’s C —", label)) +
    theme_minimal()
  
  p3 <- ggplot(G_sf) +
    geom_sf(aes(fill = hotcold), color = "white", size = 0.2) +
    scale_fill_manual(values = c(
      "Hot spot (p<0.05)" = "#b2182b",
      "Cold spot (p<0.05)" = "#2166ac",
      "Not significant" = "grey85"
    )) +
    labs(title = paste("Getis–Ord Gi* —", label), fill = NULL) +
    theme_minimal()
  
  print(p1); print(p2); print(p3)
}

analyze_spatial(jabar_merged, "Y..IPM.", "Indeks Pembangunan Manusia (IPM)")

## 
## 
## ==============================
##  Indeks Pembangunan Manusia (IPM) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 2.236, p-value = 0.02535
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.27459101       -0.03846154        0.01960115 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 2.415, p-value = 0.01573
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##         0.6069834         1.0000000         0.0264836

analyze_spatial(jabar_merged, "X1..PDRB.", "Produk Domestik Regional Bruto (PDRB)")

## 
## 
## ==============================
##  Produk Domestik Regional Bruto (PDRB) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 2.4963, p-value = 0.01255
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.29479733       -0.03846154        0.01782296 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 1.6091, p-value = 0.1076
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.69552243        1.00000000        0.03580473

analyze_spatial(jabar_merged, "X2..TPAK.", "Tingkat Partisipasi Angkatan Kerja (TPAK)")

## 
## 
## ==============================
##  Tingkat Partisipasi Angkatan Kerja (TPAK) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 0.98267, p-value = 0.3258
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.08598213       -0.03846154        0.01603727 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = 1.5182, p-value = 0.129
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        0.67735304        1.00000000        0.04516516

analyze_spatial(jabar_merged, "X3..UMR.", "Upah Minimum Rakyat (UMR)")

## 
## 
## ==============================
##  Upah Minimum Rakyat (UMR) 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = 0.037084, p-value = 0.9704
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.03318861       -0.03846154        0.02021745 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.051232, p-value = 0.9591
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.00781240        1.00000000        0.02325299

analyze_spatial(jabar_merged, "X4..Pengeluaran.per.kapita.", "Pengeluaran per Kapita")

## 
## 
## ==============================
##  Pengeluaran per Kapita 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = -1.4182, p-value = 0.1561
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.22929156       -0.03846154        0.01810562 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.57649, p-value = 0.5643
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.10680275        1.00000000        0.03432302

analyze_spatial(jabar_merged, "X5..Gini.Rasio.", "Gini Rasio")

## 
## 
## ==============================
##  Gini Rasio 
## ==============================
## 
##  Moran I test under randomisation
## 
## data:  x_raw  
## weights: lwW    
## 
## Moran I statistic standard deviate = -0.8554, p-value = 0.3923
## alternative hypothesis: two.sided
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.15846484       -0.03846154        0.01968105 
## 
## 
##  Geary C test under randomisation
## 
## data:  x_raw 
## weights: lwW   
## 
## Geary C statistic standard deviate = -0.23704, p-value = 0.8126
## alternative hypothesis: two.sided
## sample estimates:
## Geary C statistic       Expectation          Variance 
##        1.03826887        1.00000000        0.02606476

6.2.1 UJi Moran I

Variabel	Moran I	P Value	Interpretasi
Indeks Pembangunan Manusia (IPM)	0.2746	0.0254	Terdapat autokorelasi spasial positif yang signifikan. Nilai IPM cenderung berkelompok secara spasial.
Produk Domestik Regional Bruto (PDRB)	0.2948	0.0126	Terdapat autokorelasi spasial positif yang signifikan. Wilayah dengan PDRB tinggi cenderung berdekatan dengan wilayah serupa.
Tingkat Partisipasi Angkatan Kerja (TPAK)	0.0860	0.3258	Tidak terdapat autokorelasi spasial yang signifikan. Nilai TPAK cenderung acak secara spasial.
Upah Minimum Rakyat (UMR)	-0.0332	0.9704	Tidak terdapat autokorelasi spasial. Nilai UMR terdistribusi acak antar wilayah.
Pengeluaran per Kapita	-0.2293	0.1561	Tidak terdapat autokorelasi spasial yang signifikan. Sedikit indikasi pola negatif namun tidak signifikan.
Gini Rasio	-0.1585	0.3923	Tidak terdapat autokorelasi spasial yang signifikan. Ketimpangan cenderung acak secara spasial.

6.2.2 Uji Geary C

Variabel	Geary C	P Value	Interpretasi
Indeks Pembangunan Manusia (IPM)	0.6070	0.0157	Terdapat autokorelasi spasial positif yang signifikan. Nilai IPM cenderung berkelompok secara spasial.
Produk Domestik Regional Bruto (PDRB)	0.6955	0.1076	Tidak terdapat autokorelasi spasial yang signifikan. Nilai PDRB cenderung acak secara spasial.
Tingkat Partisipasi Angkatan Kerja (TPAK)	0.6774	0.1290	Tidak terdapat autokorelasi spasial yang signifikan. Nilai TPAK terdistribusi acak antar wilayah.
Upah Minimum Rakyat (UMR)	1.0078	0.9591	Tidak terdapat autokorelasi spasial. Nilai UMR terdistribusi acak antar wilayah.
Pengeluaran per Kapita	1.1068	0.5643	Tidak terdapat autokorelasi spasial yang signifikan. Nilai pengeluaran cenderung acak.
Gini Rasio	1.03	0.812	Tidak terdapat Autokorelasi spasial yang signifikan.

6.3 ESTIMASI MODEL

if (!all(c("Y..IPM.", "X1..PDRB.", "X2..TPAK.", "X3..UMR.", 
           "X4..Pengeluaran.per.kapita.", "X5..Gini.Rasio.") %in% names(jabar_merged))) {
  stop("Salah satu variabel tidak ada di data. Periksa nama kolom!")
}

jabar_merged_clean <- jabar_merged %>%
  dplyr::select(Y..IPM., X1..PDRB., X2..TPAK., X3..UMR., 
                X4..Pengeluaran.per.kapita., X5..Gini.Rasio.) %>%
  na.omit()


# Standardisasi Data 
jabar_merged_scaled <- jabar_merged_clean %>%
  mutate(
    Y..IPM._scaled = scale(Y..IPM.),
    X1..PDRB._scaled = scale(X1..PDRB.),
    X2..TPAK._scaled = scale(X2..TPAK.),
    X3..UMR._scaled = scale(X3..UMR.),
    X4..Pengeluaran.per.kapita._scaled = scale(X4..Pengeluaran.per.kapita.),
    X5..Gini._scaled = scale(X5..Gini.Rasio.)
  )

# Model Regresi Spasial
# OLS
model_ols <- lm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled
)
summary(model_ols)

## 
## Call:
## lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.1510 -0.6524 -0.2512  0.7568  1.6113 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                         7.847e-16  1.815e-01   0.000   1.0000  
## X1..PDRB._scaled                    1.878e-01  2.207e-01   0.851   0.4045  
## X2..TPAK._scaled                   -4.622e-01  2.349e-01  -1.968   0.0624 .
## X3..UMR._scaled                    -1.837e-01  2.862e-01  -0.642   0.5279  
## X4..Pengeluaran.per.kapita._scaled  2.677e-02  2.775e-01   0.096   0.9240  
## X5..Gini._scaled                    6.783e-02  2.965e-01   0.229   0.8213  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.943 on 21 degrees of freedom
## Multiple R-squared:  0.2817, Adjusted R-squared:  0.1107 
## F-statistic: 1.647 on 5 and 21 DF,  p-value: 0.1913

# LM Test
lm_tests <- lm.LMtests(model_ols, listw = lwW,
                       test = c("LMerr", "LMlag", "RLMerr", "RLMlag", "SARMA"))

## Please update scripts to use lm.RStests in place of lm.LMtests

print(lm_tests)

## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled +
## X2..TPAK._scaled + X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled
## + X5..Gini._scaled, data = jabar_merged_scaled)
## test weights: listw
## 
## RSerr = 0.20761, df = 1, p-value = 0.6486
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled +
## X2..TPAK._scaled + X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled
## + X5..Gini._scaled, data = jabar_merged_scaled)
## test weights: listw
## 
## RSlag = 1.944, df = 1, p-value = 0.1632
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled +
## X2..TPAK._scaled + X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled
## + X5..Gini._scaled, data = jabar_merged_scaled)
## test weights: listw
## 
## adjRSerr = 7.9905, df = 1, p-value = 0.004702
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled +
## X2..TPAK._scaled + X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled
## + X5..Gini._scaled, data = jabar_merged_scaled)
## test weights: listw
## 
## adjRSlag = 9.7269, df = 1, p-value = 0.001816
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = Y..IPM._scaled ~ X1..PDRB._scaled +
## X2..TPAK._scaled + X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled
## + X5..Gini._scaled, data = jabar_merged_scaled)
## test weights: listw
## 
## SARMA = 9.9345, df = 2, p-value = 0.006962

# SAR
slag_model <- lagsarlm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled,
  listw = lwB,
  method = "Matrix",      
  zero.policy = TRUE      
)

summary(slag_model)

## 
## Call:lagsarlm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled, listw = lwB, method = "Matrix", 
##     zero.policy = TRUE)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.82475 -1.16191 -0.29868  1.14731  2.81353 
## 
## Type: lag 
## Coefficients: (asymptotic standard errors) 
##                                     Estimate Std. Error z value Pr(>|z|)
## (Intercept)                        -0.459877   0.260106 -1.7680  0.07705
## X1..PDRB._scaled                    0.426287   0.316491  1.3469  0.17801
## X2..TPAK._scaled                   -0.680631   0.336518 -2.0226  0.04312
## X3..UMR._scaled                     0.169051   0.410242  0.4121  0.68028
## X4..Pengeluaran.per.kapita._scaled -0.199895   0.397570 -0.5028  0.61511
## X5..Gini._scaled                   -0.091107   0.424822 -0.2145  0.83019
## 
## Rho: -0.48194, LR test value: 16.414, p-value: 5.0915e-05
## Asymptotic standard error: 0.0033266
##     z-value: -144.87, p-value: < 2.22e-16
## Wald statistic: 20988, p-value: < 2.22e-16
## 
## Log likelihood: -25.12776 for lag model
## ML residual variance (sigma squared): 1.8252, (sigma: 1.351)
## Number of observations: 27 
## Number of parameters estimated: 8 
## AIC: 66.256, (AIC for lm: 80.669)
## LM test for residual autocorrelation
## test value: 4.6821, p-value: 0.030478

# SEM
serr_model <- errorsarlm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled, listw = lwB,
  zero.policy = TRUE, method = "Matrix"
)
summary(serr_model)

## 
## Call:
## errorsarlm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled, listw = lwB, method = "Matrix", 
##     zero.policy = TRUE)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.76458 -0.70068  0.16754  1.05047  1.57175 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##                                     Estimate Std. Error z value Pr(>|z|)
## (Intercept)                        -0.185194   0.068543 -2.7019 0.006895
## X1..PDRB._scaled                    0.374194   0.258477  1.4477 0.147705
## X2..TPAK._scaled                   -0.451161   0.256318 -1.7602 0.078380
## X3..UMR._scaled                    -0.231744   0.297809 -0.7782 0.436473
## X4..Pengeluaran.per.kapita._scaled -0.057252   0.295901 -0.1935 0.846581
## X5..Gini._scaled                   -0.296841   0.311236 -0.9538 0.340210
## 
## Lambda: -0.5, LR test value: 32.62, p-value: 1.1204e-08
## Asymptotic standard error: 0.01352
##     z-value: -36.983, p-value: < 2.22e-16
## Wald statistic: 1367.7, p-value: < 2.22e-16
## 
## Log likelihood: -17.02446 for error model
## ML residual variance (sigma squared): 1.0353, (sigma: 1.0175)
## Number of observations: 27 
## Number of parameters estimated: 8 
## AIC: 50.049, (AIC for lm: 80.669)

# SDM
sdm_model <- lagsarlm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled, listw = lwB,
  Durbin = TRUE, method = "Matrix", zero.policy = TRUE
)
summary(sdm_model)

## 
## Call:lagsarlm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled, listw = lwB, Durbin = TRUE, method = "Matrix", 
##     zero.policy = TRUE)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.45733 -0.72120  0.12137  0.78479  1.24211 
## 
## Type: mixed 
## Coefficients: (asymptotic standard errors) 
##                                        Estimate Std. Error z value Pr(>|z|)
## (Intercept)                             0.71624    0.41892  1.7097 0.087314
## X1..PDRB._scaled                        0.54448    0.25060  2.1727 0.029805
## X2..TPAK._scaled                       -0.26376    0.28509 -0.9252 0.354871
## X3..UMR._scaled                        -0.55129    0.36378 -1.5155 0.129656
## X4..Pengeluaran.per.kapita._scaled      0.14407    0.35016  0.4114 0.680749
## X5..Gini._scaled                       -0.12607    0.30724 -0.4103 0.681556
## lag.(Intercept)                        -0.32788    0.10440 -3.1405 0.001686
## lag.X1..PDRB._scaled                    0.58194    0.28962  2.0093 0.044504
## lag.X2..TPAK._scaled                    0.17780    0.25827  0.6884 0.491193
## lag.X3..UMR._scaled                    -0.23294    0.23980 -0.9714 0.331354
## lag.X4..Pengeluaran.per.kapita._scaled  0.30718    0.27354  1.1230 0.261437
## lag.X5..Gini._scaled                   -0.43415    0.23415 -1.8541 0.063720
## 
## Rho: -0.48893, LR test value: 21.928, p-value: 2.8307e-06
## Asymptotic standard error: 0.0066035
##     z-value: -74.041, p-value: < 2.22e-16
## Wald statistic: 5482.1, p-value: < 2.22e-16
## 
## Log likelihood: -13.62372 for mixed model
## ML residual variance (sigma squared): 0.78859, (sigma: 0.88803)
## Number of observations: 27 
## Number of parameters estimated: 14 
## AIC: 55.247, (AIC for lm: 75.175)
## LM test for residual autocorrelation
## test value: 0.71787, p-value: 0.39684

# SLX
slx_model <- lmSLX(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled, listw = lwB,
  zero.policy = TRUE
)
summary(slx_model)

## 
## Call:
## lm(formula = formula(paste("y ~ ", paste(colnames(x)[-1], collapse = "+"))), 
##     data = as.data.frame(x), weights = weights)
## 
## Coefficients:
##                                         Estimate  Std. Error  t value 
## (Intercept)                              0.84935   0.37564     2.26104
## X1..PDRB._scaled                         0.31225   0.22762     1.37181
## X2..TPAK._scaled                        -0.25303   0.25358    -0.99785
## X3..UMR._scaled                         -0.55831   0.32680    -1.70845
## X4..Pengeluaran.per.kapita._scaled       0.16523   0.31625     0.52246
## X5..Gini._scaled                         0.11267   0.27919     0.40354
## lag..Intercept.                         -0.23268   0.09248    -2.51604
## lag.X1..PDRB._scaled                     0.28214   0.24582     1.14775
## lag.X2..TPAK._scaled                     0.09241   0.23261     0.39727
## lag.X3..UMR._scaled                     -0.13730   0.20767    -0.66115
## lag.X4..Pengeluaran.per.kapita._scaled   0.13476   0.23238     0.57993
## lag.X5..Gini._scaled                    -0.16983   0.21076    -0.80580
##                                         Pr(>|t|)
## (Intercept)                              0.03905
## X1..PDRB._scaled                         0.19029
## X2..TPAK._scaled                         0.33418
## X3..UMR._scaled                          0.10816
## X4..Pengeluaran.per.kapita._scaled       0.60898
## X5..Gini._scaled                         0.69225
## lag..Intercept.                          0.02374
## lag.X1..PDRB._scaled                     0.26905
## lag.X2..TPAK._scaled                     0.69677
## lag.X3..UMR._scaled                      0.51854
## lag.X4..Pengeluaran.per.kapita._scaled   0.57057
## lag.X5..Gini._scaled                     0.43295

# SDEM
sdem_model <- errorsarlm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled, listw = lwB,
  Durbin = TRUE, method = "Matrix", zero.policy = TRUE
)
summary(sdem_model)

## 
## Call:
## errorsarlm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled, listw = lwB, Durbin = TRUE, method = "Matrix", 
##     zero.policy = TRUE)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.8054169 -0.2459640 -0.0082868  0.2403099  0.6517561 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##                                         Estimate Std. Error z value  Pr(>|z|)
## (Intercept)                             1.371071   0.204152  6.7159 1.869e-11
## X1..PDRB._scaled                        0.285685   0.123234  2.3182  0.020437
## X2..TPAK._scaled                       -0.009747   0.144530 -0.0674  0.946232
## X3..UMR._scaled                        -1.067640   0.211035 -5.0591 4.213e-07
## X4..Pengeluaran.per.kapita._scaled      0.597136   0.183474  3.2546  0.001136
## X5..Gini._scaled                        0.584620   0.208323  2.8063  0.005011
## lag.(Intercept)                        -0.374958   0.047690 -7.8623 3.775e-15
## lag.X1..PDRB._scaled                    0.974853   0.122725  7.9434 1.998e-15
## lag.X2..TPAK._scaled                    0.589902   0.099037  5.9564 2.579e-09
## lag.X3..UMR._scaled                    -0.481425   0.102908 -4.6782 2.894e-06
## lag.X4..Pengeluaran.per.kapita._scaled  0.633235   0.107817  5.8732 4.274e-09
## lag.X5..Gini._scaled                   -0.925943   0.148633 -6.2297 4.672e-10
## 
## Lambda: -0.5, LR test value: 70.466, p-value: < 2.22e-16
## Asymptotic standard error: 0.01352
##     z-value: -36.983, p-value: < 2.22e-16
## Wald statistic: 1367.7, p-value: < 2.22e-16
## 
## Log likelihood: 10.64539 for error model
## ML residual variance (sigma squared): 0.13333, (sigma: 0.36514)
## Number of observations: 27 
## Number of parameters estimated: 14 
## AIC: 6.7092, (AIC for lm: 75.175)

Model	Hasil Uji	Signifikansi	Interpretasi
OLS (Ordinary Least Squares)	R² = 0.2817, Adjusted R² = 0.1107, F p-value = 0.1913	Tidak Signifikan	Tidak terdapat hubungan linear yang kuat antar variabel; model OLS tidak memadai, perlu uji spasial.
LM Test (SARMA)	SARMA = 9.9345, p = 0.006962	Signifikan	OLS tidak memadai; terdapat indikasi dependensi spasial, sehingga perlu model spasial (SAR/SEM).
SAR (Spatial Lag Model)	ρ = -0.4819, LR test = 16.414, p = 5.09e-05; AIC = 66.26	Spasial Signifikan	IPM dipengaruhi oleh efek lag spasial negatif; terdapat interaksi antar wilayah.
SEM (Spatial Error Model)	λ = -0.5, LR test = 32.62, p = 1.12e-08; AIC = 50.05	Spasial Signifikan	Terdapat efek spasial pada error.
SDM (Spatial Durbin Model)	ρ = -0.4889, LR test = 21.928, p = 2.83e-06; AIC = 55.25	Spasial Signifikan	Efek lag pada variabel independen signifikan; menunjukkan adanya interaksi antar wilayah.
SDEM (Spatial Durbin Error Model)	λ = -0.5, LR test = 70.466, p < 2.22e-16; AIC = 6.71	Spasial Signifikan	IPM dipengaruhi oleh efek spasial langsung dan tidak langsung.

6.4 PERBANDINGAN MODEL

model_list <- list(
  OLS  = model_ols,
  SEM  = serr_model,
  SDM  = sdm_model,
  SLX  = slx_model,
  SDEM = sdem_model,
  SAR  = slag_model
)

model_comparison <- data.frame(
  Model  = names(model_list),
  AIC    = sapply(model_list, AIC),
  LogLik = sapply(model_list, logLik)
)

print(model_comparison)

##      Model       AIC    LogLik
## OLS    OLS 80.669252 -33.33463
## SEM    SEM 50.048911 -17.02446
## SDM    SDM 55.247439 -13.62372
## SLX    SLX 75.175433 -24.58772
## SDEM  SDEM  6.709216  10.64539
## SAR    SAR 66.255528 -25.12776

Model	AIC	LogLik
OLS	80.6693	-33.3346
SEM	50.0489	-17.0245
SDM	55.2474	-13.6237
SLX	75.1754	-24.5877
SDEM	6.7092	10.6454
SAR	66.2555	-25.1278

Interpretasi:

Model terbaik yang digunakan adalah Spatial Durbin Model (SDM), meskipun AIC bukan yang terendah, karena memberikan hasil yang paling signifikan secara ekonomi dan spasial, serta memiliki interpretasi paling komprehensif terhadap hubungan antar wilayah.

6.5 INTERPRETASI HASIL

# SDM
sdm_model <- lagsarlm(
  Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + X3..UMR._scaled +
    X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled,
  data = jabar_merged_scaled, listw = lwB,
  Durbin = TRUE, method = "Matrix", zero.policy = TRUE
)
summary(sdm_model)

## 
## Call:lagsarlm(formula = Y..IPM._scaled ~ X1..PDRB._scaled + X2..TPAK._scaled + 
##     X3..UMR._scaled + X4..Pengeluaran.per.kapita._scaled + X5..Gini._scaled, 
##     data = jabar_merged_scaled, listw = lwB, Durbin = TRUE, method = "Matrix", 
##     zero.policy = TRUE)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.45733 -0.72120  0.12137  0.78479  1.24211 
## 
## Type: mixed 
## Coefficients: (asymptotic standard errors) 
##                                        Estimate Std. Error z value Pr(>|z|)
## (Intercept)                             0.71624    0.41892  1.7097 0.087314
## X1..PDRB._scaled                        0.54448    0.25060  2.1727 0.029805
## X2..TPAK._scaled                       -0.26376    0.28509 -0.9252 0.354871
## X3..UMR._scaled                        -0.55129    0.36378 -1.5155 0.129656
## X4..Pengeluaran.per.kapita._scaled      0.14407    0.35016  0.4114 0.680749
## X5..Gini._scaled                       -0.12607    0.30724 -0.4103 0.681556
## lag.(Intercept)                        -0.32788    0.10440 -3.1405 0.001686
## lag.X1..PDRB._scaled                    0.58194    0.28962  2.0093 0.044504
## lag.X2..TPAK._scaled                    0.17780    0.25827  0.6884 0.491193
## lag.X3..UMR._scaled                    -0.23294    0.23980 -0.9714 0.331354
## lag.X4..Pengeluaran.per.kapita._scaled  0.30718    0.27354  1.1230 0.261437
## lag.X5..Gini._scaled                   -0.43415    0.23415 -1.8541 0.063720
## 
## Rho: -0.48893, LR test value: 21.928, p-value: 2.8307e-06
## Asymptotic standard error: 0.0066035
##     z-value: -74.041, p-value: < 2.22e-16
## Wald statistic: 5482.1, p-value: < 2.22e-16
## 
## Log likelihood: -13.62372 for mixed model
## ML residual variance (sigma squared): 0.78859, (sigma: 0.88803)
## Number of observations: 27 
## Number of parameters estimated: 14 
## AIC: 55.247, (AIC for lm: 75.175)
## LM test for residual autocorrelation
## test value: 0.71787, p-value: 0.39684

Produk Domestik Regional Bruto (PDRB) berpengaruh positif dan signifikan terhadap IPM dengan nilai koefisien 0.54448, artinya peningkatan PDRB di suatu daerah berkontribusi langsung terhadap peningkatan IPM. Efek spasialnya juga signifikan dan positif (0.58194), menunjukkan adanya spillover effect di mana pertumbuhan ekonomi di satu wilayah dapat mendorong peningkatan IPM di wilayah sekitarnya.

Tingkat Partisipasi Angkatan Kerja (TPAK) berpengaruh negatif namun tidak signifikan terhadap IPM, menandakan bahwa peningkatan partisipasi tenaga kerja belum tentu diikuti dengan peningkatan kualitas pembangunan manusia. Efek spasialnya juga tidak signifikan, mengindikasikan bahwa pengaruh TPAK antarwilayah relatif lemah.

Upah Minimum Regional (UMR) menunjukkan hubungan negatif namun tidak signifikan terhadap IPM, yang mungkin disebabkan oleh kenaikan upah yang tidak diimbangi dengan produktivitas tenaga kerja. Efek spasial UMR juga tidak signifikan, menunjukkan bahwa kebijakan upah di daerah sekitar belum memberikan dampak nyata terhadap kesejahteraan wilayah tetangga.

Pengeluaran per kapita memiliki pengaruh positif namun tidak signifikan terhadap IPM, artinya peningkatan daya beli masyarakat cenderung meningkatkan kualitas hidup, meski belum cukup kuat secara statistik. Sementara itu, variabel Gini Rasio berpengaruh negatif terhadap IPM, mencerminkan bahwa semakin tinggi ketimpangan pendapatan, semakin rendah tingkat pembangunan manusia, meskipun pengaruhnya belum signifikan.

Nilai ρ (rho) sebesar -0.48893 dengan p-value < 0.001 menunjukkan adanya efek spasial negatif yang signifikan, mengindikasikan bahwa peningkatan IPM di suatu wilayah dapat menekan IPM di wilayah tetangganya akibat kompetisi sumber daya. Nilai AIC sebesar 55.247 dan LogLik -13.62 menunjukkan bahwa model SDM memiliki kecocokan yang baik dan lebih efisien dibanding model regresi linier biasa, sedangkan standar error yang kecil (0.0066) menandakan kestabilan hasil estimasi.

7 KESIMPULAN DAN SARAN

Berdasarkan hasil analisis, model terbaik untuk memodelkan Indeks Pembangunan Manusia (IPM) di Jawa Barat adalah Spatial Durbin Model (SDM). Model ini mampu menggambarkan hubungan IPM dengan variabel PDRB, TPAK, UMR, Pengeluaran Per Kapita, dan Gini Rasio secara cukup baik. Penelitian ini menunjukkan bahwa SDM merupakan model yang optimal dalam menjelaskan pengaruh langsung dan efek spasial antarwilayah. Namun, masih terdapat potensi omitted variable yang belum terakomodasi secara penuh. Oleh karena itu, penelitian selanjutnya disarankan menggunakan metode estimasi yang lebih kompleks agar mampu mencakup variabel-variabel tambahan yang relevan dan memberikan pemodelan IPM yang lebih komprehensif.

8 LAMPIRAN

Daftar Pustaka

[1] Badan Pusat Statistik Provinsi Jawa Barat. Indeks Pembangunan Manusia Menurut Kabupaten/Kota di Jawa Barat Tahun 2024.
https://jabar.bps.go.id/id/publication/2025/04/11/3e25c1657d4a90e6ef256015/indeks-pembangunan-manusia-provinsi-jawa-barat-2024.html

[2] Badan Pusat Statistik Provinsi Jawa Barat. Produk Domestik Regional Bruto (PDRB) Menurut Kabupaten/Kota Atas Dasar Harga Berlaku, Seri 2010.
https://jabar.bps.go.id/id/statistics-table/2/MTEwIzI%3D/-seri-2010–pdrb-atas-dasar-harga-berlaku-menurut-kabupaten-kota-provinsi-jawa-barat.html

[3] Open Data Jabar. Indeks Pengetahuan Berdasarkan Kabupaten/Kota di Jawa Barat.
https://opendata.jabarprov.go.id/id/dataset/indeks-pengetahuan-berdasarkan-kabupatenkota-di-jawa-barat

[4] Open Data Jabar. Daftar Upah Minimum dan Pengeluaran Per Kapita Kabupaten/Kota di Jawa Barat.
https://opendata.jabarprov.go.id/id/dataset/pengeluaran-per-kapita-disesuaikan-berdasarkan-kabupatenkota-di-jawa-barat

[5] Anselin, L. Local Indicators of Spatial Association – LISA. Geographical Analysis, Vol. 27, No. 2, 1995.
https://dces.webhosting.cals.wisc.edu/wp-content/uploads/sites/128/2013/08/W4_Anselin1995.pdf

[6] Tarmizi, H. dan Sari, D. (2022). Analisis Faktor-faktor yang Mempengaruhi Indeks Pembangunan Manusia di Indonesia Menggunakan Pendekatan Ekonometrika. Jurnal Ekonomi dan Kebijakan Publik Indonesia, 9(1), 45–56.
https://journal.unesa.ac.id/index.php/jekpi/article/view/15483

[7] Nugroho, A. dan Rahmawati, F. (2023). Pengaruh PDRB, Pendidikan, dan Upah terhadap Pembangunan Manusia di Indonesia: Pendekatan Regresi Panel. Jurnal Ekonomi Pembangunan dan Perencanaan Daerah, 6(2), 120–134.
https://ejournal.unesa.ac.id/index.php/JEP/article/view/41239

[8] Prasetyo, A. dan Astuti, W. (2021). Spatial Autocorrelation Analysis of Human Development Index in Indonesia. International Journal of Regional Development, 8(1), 35–47.
https://www.researchgate.net/publication/348354709_Spatial_Autocorrelation_Analysis_of_Human_Development_Index_in_Indonesia

[9] LeSage, J. dan Pace, R. K. Introduction to Spatial Econometrics. CRC Press, 2009.

https://www.crcpress.com/Introduction-to-Spatial-Econometrics/LeSage-Pace/p/book/9781584886105

[10] Elhorst, J. P. Spatial Econometrics: From Cross-Sectional Data to Spatial Panels. Springer, 2014.

https://link.springer.com/book/10.1007/978-3-319-23536-7

[11] Anselin, L. Local Indicators of Spatial Association LISA. Geographical Analysis, 1995.

https://dces.webhosting.cals.wisc.edu/wp-content/uploads/sites/128/2013/08/W4_Anselin1995.pdf

[12] Getis, A. dan Ord, J. K. The Analysis of Spatial Association by Use of Distance Statistics. Geographical Analysis, 1992.

https://onlinelibrary.wiley.com/doi/10.1111/j.1538-4632.1992.tb00261.x