Analisis Kasus Kusta di Kabupaten/Kota
Provinsi Jawa Barat Tahun
2021
Ditujukan Guna Memenuhi UAS Mata Kuliah Spasial
Disusun oleh:
Muhammad Imamul Caesar
(140610230019)
Dosen Pengampu:
I Gede Nyoman Mindra Jaya,
Ph.D
PROGRAM STUDI S-1 STATISTIKA
FAKULTAS MATEMATIKA DAN ILMU
PENGETAHUAN ALAM
UNIVERSITAS PADJADJARAN
JATINANGOR
2025
I. Pendahuluan
Penyakit kusta (lepra) merupakan salah satu penyakit menular kronis yang hingga saat ini masih menjadi tantangan kesehatan masyarakat, khususnya di negara-negara berkembang. Kusta disebabkan oleh Mycobacterium leprae dan Mycobacterium lepromatosis yang menyerang kulit, saraf perifer, saluran pernapasan atas, serta mata. Penyakit ini berkembang secara perlahan dengan masa inkubasi yang panjang, sehingga sering kali terlambat terdeteksi dan ditangani. Apabila tidak mendapatkan pengobatan yang tepat, kusta dapat menyebabkan kecacatan permanen, gangguan fungsi organ, serta dampak sosial yang serius berupa stigma dan diskriminasi terhadap penderita [1].
Meskipun kusta dapat disembuhkan melalui terapi multidrug therapy (MDT), penyakit ini belum sepenuhnya dapat dieliminasi. Organisasi Kesehatan Dunia (WHO) melaporkan bahwa lebih dari 200.000 kasus baru kusta masih ditemukan setiap tahunnya secara global, dengan konsentrasi terbesar berada di Asia Selatan, Asia Tenggara, dan Amerika Latin [2]. Indonesia secara konsisten menempati peringkat ketiga dunia dengan jumlah kasus baru kusta tertinggi setelah India dan Brasil. Fakta ini menunjukkan bahwa kusta masih menjadi masalah kesehatan masyarakat yang signifikan dan membutuhkan perhatian serius, baik dari sisi medis maupun kebijakan kesehatan [3].
Di Indonesia, kusta tidak hanya menjadi persoalan klinis, tetapi juga masalah sosial dan struktural. Data Kementerian Kesehatan Republik Indonesia menunjukkan bahwa jumlah kasus baru kusta masih berada di atas 10.000 kasus per tahun, dengan penyebaran yang tidak merata antarwilayah [4]. Beberapa provinsi, termasuk Jawa Barat, Jawa Timur, dan Jawa Tengah, tercatat sebagai wilayah dengan kontribusi kasus yang relatif tinggi. Kondisi ini mengindikasikan adanya faktor-faktor lokal yang berperan dalam mempertahankan transmisi penyakit, seperti kepadatan penduduk, kondisi sanitasi, tingkat pendidikan, serta kemiskinan [5].
Pola penyebaran penyakit kusta yang tidak homogen secara geografis menunjukkan bahwa pendekatan analisis spasial menjadi sangat relevan. Analisis spasial memungkinkan peneliti untuk mengidentifikasi klaster wilayah dengan risiko tinggi, mendeteksi adanya autokorelasi spasial, serta memahami hubungan antara faktor sosial-ekonomi dan kejadian penyakit secara lebih komprehensif [6]. Tanpa pendekatan spasial, analisis statistik konvensional berpotensi mengabaikan ketergantungan antarwilayah yang secara nyata terjadi dalam fenomena epidemiologi penyakit menular [7].
Pentingnya analisis spasial dalam kajian penyakit kusta juga berkaitan dengan upaya perencanaan dan evaluasi kebijakan kesehatan yang lebih tepat sasaran. Dengan mengetahui wilayah-wilayah yang memiliki tingkat kejadian tinggi atau membentuk klaster tertentu, pemerintah dan pemangku kepentingan dapat memfokuskan intervensi kesehatan, seperti peningkatan skrining aktif, perbaikan sanitasi lingkungan, serta edukasi masyarakat, pada daerah yang paling membutuhkan [8]. Pendekatan ini sejalan dengan prinsip efisiensi dan efektivitas dalam sistem kesehatan masyarakat.
Selain itu, penelitian mengenai kusta memiliki keterkaitan erat dengan Tujuan Pembangunan Berkelanjutan (Sustainable Development Goals / SDGs). Upaya pengendalian dan eliminasi kusta berkontribusi langsung terhadap SDG 3, yaitu menjamin kehidupan yang sehat dan mendorong kesejahteraan bagi semua orang di segala usia. Lebih jauh, kusta juga berkaitan dengan SDG 1 (pengentasan kemiskinan), SDG 4 (pendidikan berkualitas), dan SDG 10 (pengurangan ketimpangan), mengingat penyakit ini sering kali ditemukan pada kelompok masyarakat yang rentan secara sosial dan ekonomi [9]. Oleh karena itu, kajian kusta tidak dapat dipisahkan dari konteks pembangunan berkelanjutan dan keadilan sosial.
Berdasarkan latar belakang tersebut, penelitian ini bertujuan untuk menganalisis pola spasial penyakit kusta di Provinsi Jawa Barat menggunakan pendekatan statistika spasial. Melalui pemetaan kasus, analisis autokorelasi spasial, serta pemodelan spasial, diharapkan penelitian ini dapat memberikan gambaran yang lebih jelas mengenai distribusi penyakit kusta dan faktor-faktor yang memengaruhinya. Hasil penelitian ini diharapkan tidak hanya bermanfaat secara akademik sebagai penerapan metode statistika spasial, tetapi juga dapat menjadi masukan awal bagi perencanaan kebijakan kesehatan berbasis wilayah.
II. Tinjauan Pustaka
2.1. Spatial Dependence
Dalam data spasial, pengamatan pada suatu wilayah sering kali tidak independen terhadap wilayah lain. Fenomena ini dikenal sebagai spatial dependence (ketergantungan spasial), yaitu kondisi ketika nilai variabel di suatu lokasi dipengaruhi oleh nilai variabel yang sama atau variabel lain di lokasi sekitarnya. Ketergantungan ini muncul karena adanya proses difusi, keterhubungan aktivitas manusia, kedekatan geografis, maupun kesamaan karakteristik lingkungan antarwilayah [11] [12]. Pada konteks penyakit menular seperti kusta, ketergantungan spasial dapat mencerminkan adanya pola penyebaran yang mengikuti kedekatan wilayah (misalnya mobilitas penduduk antar kabupaten/kota) atau kesamaan faktor risiko lokal.
Spatial dependence biasanya dimodelkan melalui struktur ketetanggaan (neighbor structure) yang direpresentasikan dalam spatial weights matrix (matriks bobot spasial) W, misalnya berbasis kontiguitas (queen/rook) atau berbasis jarak [13]. Matriks W menjadi inti dari banyak metode analisis spasial karena mendefinisikan “siapa bertetangga dengan siapa” dan seberapa kuat pengaruhnya. Jika spatial dependence diabaikan, hasil analisis statistik klasik berisiko bias atau tidak efisien karena asumsi independensi residual dilanggar [11] [12].
2.2. Autokorelasi Spasial
Autokorelasi spasial adalah ukuran sejauh mana nilai suatu variabel di suatu lokasi berkorelasi dengan nilai variabel yang sama di lokasi-lokasi sekitar. Autokorelasi spasial dapat bersifat positif (wilayah bernilai tinggi dikelilingi tinggi, atau rendah dikelilingi rendah) atau negatif (wilayah tinggi dikelilingi rendah dan sebaliknya) [13] [14]. Dalam kajian kesehatan berbasis wilayah, autokorelasi spasial sering menjadi indikasi awal adanya klaster penyakit atau ketimpangan faktor risiko antarwilayah.
Pengujian autokorelasi spasial global yang umum digunakan adalah Moran’s I dan Geary’s C. Moran’s I lebih sensitif terhadap kesamaan nilai antar tetangga (clustering), sedangkan Geary’s C lebih sensitif terhadap perbedaan lokal (dissimilarity) [13]. Untuk mengidentifikasi pola lokal, digunakan LISA (Local Indicators of Spatial Association), khususnya local Moran’s I, yang dapat memetakan tipe klaster seperti High-High, Low-Low, High-Low, dan Low-High [14]. Output LISA membantu peneliti menentukan wilayah prioritas intervensi karena mengungkap lokasi-lokasi yang menjadi kantong kasus tinggi sekaligus berdekatan dengan wilayah tinggi lainnya.
2.3 Model Spatial Econometrics
Ketika ditemukan spatial dependence atau autokorelasi spasial, analisis regresi biasa (OLS) sering tidak memadai karena residual cenderung berkorelasi antarwilayah. Spatial econometrics menyediakan keluarga model regresi yang secara eksplisit memasukkan struktur spasial melalui matriks bobot W, sehingga estimasi parameter dan inferensi menjadi lebih valid [11] [12] [15].
Model yang sering dipakai antara lain:
Spatial Autoregressive Model (SAR / Spatial Lag Model): memasukkan pengaruh nilai variabel dependen di wilayah tetangga (Wy) ke dalam model. Model ini sesuai ketika proses penyebaran/penularan atau efek “spillover” antarwilayah dominan terjadi pada outcome [11] [15].
Spatial Error Model (SEM): menangkap autokorelasi pada komponen error (misalnya ada faktor tak terukur yang bersifat spasial seperti akses layanan kesehatan atau kualitas surveilans) [11] [15].
Spatial Durbin Model (SDM) dan SDEM: memperluas model dengan memasukkan lag spasial dari kovariat (WX) sehingga mampu menangkap spillover dari faktor risiko, bukan hanya outcome [12] [15].
SLX: fokus pada lag spasial kovariat (WX) tanpa memasukkan Wy; sering dipakai sebagai model spillover yang lebih sederhana [12].
SAC: menggabungkan lag dependen (SAR) dan error spasial (SEM) sekaligus, sehingga lebih fleksibel jika ada indikasi dua mekanisme ketergantungan berjalan bersamaan [11] [15].
Pemilihan model terbaik umumnya dilakukan melalui metrik seperti AIC/BIC, uji diagnostik (misal uji LM untuk spesifikasi), serta evaluasi kinerja residual [11] [12]. Dalam konteks penelitian kusta berbasis kabupaten/kota, model spatial econometrics berguna untuk menilai apakah variasi kasus lebih dipengaruhi oleh spillover kasus di wilayah tetangga (SAR/SDM) atau oleh faktor laten spasial yang tidak terukur (SEM/SDEM), sehingga interpretasi kebijakan menjadi lebih tepat.
2.4 Model Spatial Interpolasi
Selain pemodelan hubungan sebab-akibat, analisis spasial juga sering membutuhkan pemetaan permukaan nilai di lokasi yang tidak terobservasi. Proses ini disebut spatial interpolation. Dua pendekatan yang paling umum adalah IDW dan Kriging [16] [17].
Inverse Distance Weighting (IDW) adalah metode deterministik yang mengestimasi nilai pada suatu titik sebagai rata-rata berbobot dari titik-titik observasi terdekat, dengan bobot berbanding terbalik dengan jarak (dipangkatkan oleh parameter power p) [16]. Semakin besar p, pengaruh titik dekat semakin dominan, sehingga permukaan menjadi lebih “tajam”. IDW relatif mudah dan cepat, tetapi tidak memodelkan struktur ketergantungan spasial secara probabilistik.
Kriging adalah metode geostatistik berbasis model yang memanfaatkan variogram untuk menggambarkan struktur korelasi spasial dan menghasilkan estimator tak bias dengan varians minimum (best linear unbiased predictor) [16] [17]. Kriging tidak hanya memberi prediksi, tetapi juga ketidakpastian prediksi. Namun, kriging mensyaratkan spesifikasi variogram dan umumnya bekerja lebih stabil pada data dengan sistem koordinat proyeksi (meter) dibanding lon-lat (derajat), terutama untuk domain yang cukup luas [17].
Dalam studi kabupaten/kota, interpolasi biasanya dilakukan pada titik representatif wilayah seperti centroid. Hasil interpolasi dapat digunakan untuk menggambarkan pola spasial kasar dan sebagai pelengkap interpretasi peta administrasi (choropleth). Namun, perlu dicatat bahwa interpolasi berbasis centroid wilayah administrasi memiliki keterbatasan karena data awal bukan titik sampling kontinu, melainkan agregat wilayah; sehingga interpretasinya lebih cocok sebagai visualisasi pola daripada estimasi “nilai sebenarnya” pada setiap titik [16].
III. Metodologi Penelitian
3.1 Data Penelitian
Data jumlah kasus kusta tahun 2021 diperoleh dari data sekunder tingkat kabupaten/kota di Provinsi Jawa Barat. Variabel utama penelitian adalah:
y : jumlah kasus kusta
x1 : jumlah penduduk
x2 : persentase sanitasi layak
x3 : kepadatan penduduk
x4 : rata-rata lama sekolah
x5 : persentase penduduk miskin
3.2 Matriks Bobot Spasial
Untuk merepresentasikan hubungan spasial antarwilayah, penelitian ini membangun matriks bobot spasial \(\mathbf{W}\) menggunakan pendekatan queen contiguity, di mana dua wilayah dianggap bertetangga apabila memiliki sisi atau titik sudut yang bersinggungan [18]. Matriks bobot kemudian dinormalisasi baris (row-standardized) sehingga jumlah bobot pada setiap baris bernilai satu.
Matriks bobot spasial ini digunakan dalam seluruh tahapan analisis autokorelasi dan pemodelan spasial.
3.3 Analisis Autokorelasi Spasial
3.3.1 Autokorelasi Spasial Global (Moran’s I)
Autokorelasi spasial global diuji menggunakan indeks Moran’s I untuk mengetahui apakah distribusi kasus kusta menunjukkan pola mengelompok, menyebar, atau acak secara spasial. Statistik Moran’s I dirumuskan sebagai berikut [19][20]:
\[ I = \frac{n}{\sum_{i}\sum_{j} w_{ij}} \frac{\sum_{i}\sum_{j} w_{ij}(y_i - \bar{y})(y_j - \bar{y})} {\sum_{i}(y_i - \bar{y})^2} \] Nilai Moran’s I yang signifikan secara statistik menunjukkan adanya ketergantungan spasial global pada data kasus kusta.
3.3.2 Autokorelasi Spasial Lokal (LISA)
Untuk mengidentifikasi pola spasial pada tingkat lokal, digunakan Local Indicators of Spatial Association (LISA). Analisis ini memungkinkan identifikasi klaster spasial seperti high–high, low–low, high–low, dan low–high [20]. Statistik lokal Moran dirumuskan sebagai:
\[ I_i = (y_i - \bar{y}) \sum_{j} w_{ij}(y_j - \bar{y}) \] Wilayah dengan nilai LISA signifikan pada tingkat signifikansi \(\alpha = 0{,}05\) diinterpretasikan sebagai bagian dari klaster spasial yang bermakna.
3.4 Pemodelan Ekonometrika Spasial
3.4.1 Model Regresi Linear Klasik (OLS)
Sebagai model pembanding awal, digunakan regresi linear klasik:
\[ y = \beta_0 + \sum_{k=1}^{5} \beta_k x_k + \varepsilon \]
Namun, model ini tidak memperhitungkan efek ketergantungan spasial antarwilayah.
3.4.2 Model Spasial
Untuk mengakomodasi ketergantungan spasial, penelitian ini mengestimasi beberapa model ekonometrika spasial [21][22]:
Spatial Autoregressive Model (SAR)
\[ y = \rho Wy + X\beta + \varepsilon \]Spatial Error Model (SEM)
\[ y = X\beta + u, \quad u = \lambda Wu + \varepsilon \]Spatial Lag of X Model (SLX)
\[ y = X\beta + WX\theta + \varepsilon \]Spatial Durbin Model (SDM)
\[ y = \rho Wy + X\beta + WX\theta + \varepsilon \]Spatial Durbin Error Model (SDEM) dan
Spatial Autoregressive Combined Model (SAC)
Pemilihan model terbaik dilakukan berdasarkan nilai Akaike Information Criterion (AIC), Root Mean Square Error (RMSE), dan Mean Absolute Error (MAE).
3.5 Analisis Interpolasi Spasial
Untuk menggambarkan variasi spasial kasus kusta secara kontinu, digunakan metode interpolasi berbasis titik yang diperoleh dari centroid masing-masing kabupaten/kota.
3.5.1 Inverse Distance Weighting (IDW)
Metode IDW mengasumsikan bahwa pengaruh suatu titik pengamatan berkurang seiring dengan bertambahnya jarak. Estimator IDW dirumuskan sebagai [23]:
\[ \hat{z}(s_0) = \frac{\sum_{i=1}^{n} z(s_i) d_{i0}^{-p}} {\sum_{i=1}^{n} d_{i0}^{-p}} \]
dengan \(p\) sebagai parameter daya (power).
3.5.2 Kriging
Sebagai pembanding, digunakan metode Kriging, yang memanfaatkan struktur ketergantungan spasial melalui fungsi variogram [24]. Kriging menghasilkan estimator Best Linear Unbiased Estimator (BLUE) serta menyediakan informasi ketidakpastian prediksi.
IV. Hasil dan Pembahasan
4.1 Peta Deskriptif Persebaran Kasus Kusta di Jawa Barat Tahun 2021
Peta deskriptif persebaran kasus kusta di Provinsi Jawa Barat tahun 2021 menunjukkan bahwa kasus kusta tidak tersebar secara merata di seluruh wilayah, melainkan membentuk pola spasial yang cenderung mengelompok. Wilayah selatan Jawa Barat didominasi oleh warna dengan intensitas lebih rendah, yang mengindikasikan jumlah kasus kusta relatif kecil dibandingkan wilayah lainnya. Sebaliknya, beberapa kabupaten/kota di bagian utara hingga timur Jawa Barat memperlihatkan warna yang lebih terang, menandakan jumlah kasus yang lebih tinggi dan mengindikasikan adanya konsentrasi kasus pada wilayah tertentu.
Pola ini mengisyaratkan bahwa faktor kewilayahan berperan dalam variasi kasus kusta, seperti perbedaan kepadatan penduduk, tingkat mobilitas masyarakat, serta variasi akses dan kinerja layanan kesehatan dalam mendeteksi dan melaporkan kasus. Wilayah dengan jumlah kasus lebih tinggi tidak selalu mencerminkan tingkat risiko yang lebih besar, tetapi dapat pula menunjukkan sistem penemuan kasus yang lebih aktif. Oleh karena itu, peta ini memberikan gambaran awal yang penting untuk penentuan wilayah prioritas dalam upaya pengendalian kusta, sekaligus menjadi dasar bagi analisis spasial lanjutan yang lebih mendalam.
4.2 Uji Autokorelasi Spasial Global dan Lokal (Moran’s I dan LISA)
4.2.1 Uji Autokorelasi Spasial Global (Moran’s I)
Hasil uji autokorelasi spasial global menunjukkan bahwa persebaran kasus kusta di Jawa Barat tahun 2021 memiliki pola spasial yang tidak acak dan cenderung mengelompok. Nilai Moran’s I sebesar 0,495 dengan p-value yang sangat kecil menunjukkan adanya autokorelasi spasial positif yang signifikan, di mana wilayah dengan jumlah kasus kusta tinggi cenderung berdekatan dengan wilayah lain yang juga memiliki kasus tinggi, begitu pula sebaliknya untuk wilayah dengan kasus rendah. Temuan ini diperkuat oleh nilai Geary’s C sebesar 0,474 yang lebih kecil dari nilai harapannya (1) dan signifikan secara statistik, yang menandakan adanya kemiripan nilai kasus antarwilayah bertetangga. Secara keseluruhan, hasil ini menegaskan bahwa distribusi kasus kusta di Jawa Barat tahun 2021 membentuk pola pengelompokan spasial, sehingga pendekatan analisis dan penanganan berbasis wilayah menjadi relevan.
Peta spatial lag kasus kusta di Jawa Barat tahun 2021 menunjukkan bahwa nilai kasus kusta di suatu wilayah sangat dipengaruhi oleh kondisi wilayah-wilayah di sekitarnya. Wilayah selatan Jawa Barat didominasi oleh warna dengan intensitas lebih tinggi pada peta spatial lag, yang mengindikasikan bahwa daerah-daerah tersebut dikelilingi oleh kabupaten/kota dengan jumlah kasus kusta yang relatif tinggi, sehingga membentuk zona pengaruh spasial yang kuat. Sebaliknya, wilayah utara Jawa Barat cenderung menampilkan nilai spatial lag yang lebih rendah, menandakan lingkungan sekitar dengan jumlah kasus yang relatif kecil. Pola ini menegaskan adanya keterkaitan antarwilayah secara spasial dalam penyebaran kasus kusta, di mana kondisi suatu daerah tidak berdiri sendiri, melainkan berkaitan erat dengan situasi epidemiologis wilayah tetangganya.
4.2.2 Uji Autokorelasi Spasial Lokal (LISA)
Hasil analisis Local Indicators of Spatial Association (LISA) dengan taraf signifikansi 5% menunjukkan bahwa dari seluruh wilayah yang dianalisis, hanya sebagian kecil wilayah yang memiliki autokorelasi spasial lokal yang signifikan, yaitu 2 wilayah, sementara 24 wilayah lainnya tidak menunjukkan keterkaitan spasial lokal yang signifikan. Ringkasan nilai p-value yang didominasi oleh nilai di atas 0,05 mengindikasikan bahwa pada sebagian besar wilayah, pola keterkaitan antara nilai kasus kusta dengan wilayah tetangganya masih bersifat lemah atau tidak cukup kuat secara statistik. Temuan ini menunjukkan bahwa meskipun secara global kasus kusta membentuk pola pengelompokan spasial, pada skala lokal hanya terdapat beberapa wilayah tertentu yang berperan sebagai klaster atau titik khusus, sehingga identifikasi wilayah signifikan tersebut menjadi penting sebagai fokus analisis lanjutan dan prioritas intervensi.
Peta LISA cluster menunjukkan bahwa pada tingkat signifikansi 5%, hanya terdapat satu klaster High–High yang signifikan secara spasial dalam persebaran kasus kusta di Jawa Barat tahun 2021, yang terletak di wilayah utara provinsi. Klaster High–High ini menandakan bahwa wilayah tersebut memiliki jumlah kasus kusta yang tinggi dan dikelilingi oleh wilayah tetangga dengan jumlah kasus yang juga relatif tinggi, sehingga membentuk hotspot kusta yang nyata secara lokal. Sementara itu, sebagian besar wilayah lainnya tergolong tidak signifikan, yang berarti tidak menunjukkan keterkaitan spasial lokal yang cukup kuat antara nilai kasus kusta dengan wilayah sekitarnya. Temuan ini mengindikasikan bahwa meskipun secara global terdapat pengelompokan kasus kusta, pada skala lokal hanya wilayah tertentu yang benar-benar berperan sebagai pusat konsentrasi kasus dan layak menjadi prioritas utama dalam upaya pengendalian dan intervensi berbasis wilayah.
4.3 Estimasi Model Spasial Econometrics
Model OLS menghasilkan pola prediksi kasus kusta yang masih terfragmentasi dan kurang mencerminkan keterkaitan antarwilayah. Beberapa wilayah bahkan memperlihatkan nilai prediksi yang ekstrem atau tidak realistis, yang mengindikasikan bahwa model ini belum mampu menangkap ketergantungan spasial yang memang terdapat pada data kasus kusta. Dengan kata lain, OLS cenderung melihat setiap wilayah secara terpisah tanpa mempertimbangkan pengaruh wilayah di sekitarnya.
Ketika komponen spasial mulai dimasukkan melalui model SLX, SAR, dan SEM, pola prediksi menjadi lebih halus dan menunjukkan kesinambungan antarwilayah. Model SLX mulai memperhitungkan pengaruh variabel dari wilayah tetangga, sementara SAR dan SEM secara lebih eksplisit menangkap ketergantungan spasial baik pada nilai kasus maupun pada komponen galat. Hasilnya, prediksi yang dihasilkan lebih stabil dan pola pengelompokan wilayah dengan nilai kasus tinggi atau rendah menjadi lebih jelas.
Model yang lebih komprehensif seperti SDM, SDEM, dan SAC menghasilkan pola prediksi yang paling konsisten secara spasial. Wilayah dengan prediksi kasus kusta tinggi cenderung membentuk klaster yang selaras dengan kondisi wilayah sekitarnya, sesuai dengan temuan autokorelasi spasial global dan lokal sebelumnya. Model-model ini mampu mengurangi distorsi lokal dan menghasilkan transisi nilai prediksi yang lebih realistis secara geografis.
Secara keseluruhan, perbandingan ini menegaskan bahwa model regresi spasial memberikan gambaran persebaran kasus kusta yang lebih representatif dibandingkan OLS. Keberadaan pengaruh spasial antarwilayah menjadikan model spasial lebih relevan untuk analisis epidemiologi wilayah dan sebagai dasar perumusan kebijakan pengendalian kusta berbasis area prioritas.
4.4 Perbandingan Kinerja Model
Tabel perbandingan kinerja model menunjukkan bahwa model SEM memiliki nilai AIC terendah (275,39), yang menandakan kecocokan model terbaik secara keseluruhan dengan mempertimbangkan kompleksitas dan goodness of fit. Sementara itu, berdasarkan ukuran ketepatan prediksi, model SDM dan SDEM menunjukkan performa yang lebih baik dengan nilai RMSE dan MAE yang lebih kecil dibandingkan model lainnya, sehingga menghasilkan prediksi yang lebih akurat. Sebaliknya, model OLS dan SLX memiliki nilai AIC, RMSE, dan MAE yang paling tinggi, yang mengindikasikan kinerja paling lemah dan menegaskan bahwa model tanpa atau dengan keterbatasan komponen spasial kurang mampu merepresentasikan pola kasus kusta. Secara keseluruhan, hasil ini menunjukkan bahwa pemilihan model terbaik bergantung pada tujuan analisis, di mana SEM lebih unggul untuk pemodelan struktural, sedangkan SDM atau SDEM lebih sesuai untuk keperluan prediksi spasial kasus kusta.
4.5 Peta Residual Model Spasial Econometrics
Gambar peta residual memperlihatkan perbedaan yang cukup jelas antara model nonspasial dan model spasial dalam menangkap variasi kasus kusta di Jawa Barat tahun 2021. Pada model OLS, residual masih menunjukkan pola spasial yang kuat, terutama di wilayah selatan dan timur Jawa Barat, yang ditandai oleh pengelompokan warna serupa pada wilayah-wilayah yang berdekatan. Kondisi ini mengindikasikan bahwa masih terdapat informasi spasial yang belum terakomodasi oleh model, sehingga kesalahan prediksi tidak tersebar secara acak dan cenderung mengelompok.
Model SLX dan SAR mulai memperlihatkan perbaikan dengan berkurangnya intensitas dan luasan residual ekstrem di beberapa wilayah. Meskipun demikian, pola pengelompokan residual masih tampak, terutama pada wilayah tertentu yang konsisten menunjukkan nilai residual tinggi atau rendah. Hal ini menunjukkan bahwa meskipun pengaruh spasial telah mulai dimasukkan, ketergantungan spasial dalam data belum sepenuhnya tertangkap oleh kedua model tersebut.
Perbaikan yang lebih nyata terlihat pada model SEM, SDM, dan SAC, di mana residual tampak lebih menyebar dan tidak lagi membentuk klaster spasial yang kuat. Pola warna yang lebih acak mengindikasikan bahwa struktur ketergantungan spasial, baik pada komponen galat maupun pada variabel dependen dan independen, telah terakomodasi dengan lebih baik. Dengan demikian, model-model ini mampu mengurangi bias spasial yang sebelumnya muncul pada OLS dan model spasial yang lebih sederhana.
Secara keseluruhan, analisis peta residual menegaskan bahwa model regresi spasial yang lebih komprehensif menghasilkan kinerja yang lebih baik dalam memodelkan kasus kusta. Residual yang cenderung acak dan tidak berpola secara spasial menunjukkan bahwa model telah sesuai dengan asumsi independensi galat, sehingga hasil estimasi dan prediksi yang dihasilkan lebih reliabel untuk digunakan sebagai dasar analisis dan perumusan kebijakan pengendalian kusta berbasis wilayah.
4.6 Spasial Interpolasi
4.6.1 Metode Inverse Distance Weighting (IDW)
Peta hasil interpolasi spasial menggunakan metode Inverse Distance Weighting (IDW) dengan parameter p = 2 menggambarkan pola sebaran kontinu intensitas kasus kusta di Jawa Barat berdasarkan kedekatan spasial antarwilayah. Warna yang lebih gelap hingga kemerahan menunjukkan nilai interpolasi yang lebih tinggi, yang tampak terkonsentrasi di wilayah utara dan timur Jawa Barat, khususnya di sekitar kawasan penyangga Jakarta dan wilayah pantai utara. Pola ini menunjukkan bahwa area dengan titik observasi kasus kusta yang tinggi memberikan pengaruh yang kuat terhadap wilayah sekitarnya, sehingga membentuk zona nilai interpolasi tinggi yang menyebar secara gradual seiring bertambahnya jarak dari pusat kasus.
Sebaliknya, wilayah selatan Jawa Barat didominasi oleh warna yang lebih terang, yang mengindikasikan nilai interpolasi kasus kusta yang relatif lebih rendah. Transisi warna yang halus dari utara ke selatan mencerminkan asumsi utama metode IDW, yaitu bahwa lokasi yang berdekatan memiliki karakteristik yang lebih mirip dibandingkan lokasi yang berjauhan. Meskipun peta ini efektif dalam menampilkan gradien spasial dan potensi area prioritas secara visual, hasil interpolasi IDW tetap bersifat deskriptif dan sangat bergantung pada distribusi titik data, sehingga tidak secara langsung merepresentasikan proses epidemiologis yang mendasari penyebaran penyakit.
4.6.2 Metode Kriging
Peta hasil interpolasi spasial menggunakan metode Kriging menunjukkan pola sebaran kasus kusta di Jawa Barat tahun 2021 yang lebih halus dan terstruktur dibandingkan metode IDW. Nilai interpolasi yang tinggi tampak terkonsentrasi di wilayah utara hingga timur Jawa Barat, ditandai dengan warna merah hingga keunguan, yang mengindikasikan adanya area dengan intensitas kasus kusta yang relatif tinggi dan saling berhubungan secara spasial. Pola ini mencerminkan hasil estimasi Kriging yang mempertimbangkan struktur autokorelasi spasial antar lokasi, sehingga nilai di suatu titik tidak hanya dipengaruhi oleh jarak, tetapi juga oleh pola variabilitas spasial yang ditangkap melalui semivariogram.
Sementara itu, wilayah selatan Jawa Barat didominasi oleh warna yang lebih terang, menunjukkan nilai interpolasi kasus kusta yang relatif rendah dan lebih homogen. Transisi nilai dari wilayah utara ke selatan terlihat lebih gradual dan realistis, mencerminkan keunggulan Kriging dalam menghasilkan permukaan prediksi yang stabil serta mampu meminimalkan fluktuasi ekstrem. Dengan karakteristik tersebut, hasil interpolasi Kriging tidak hanya memberikan gambaran visual persebaran kasus kusta, tetapi juga menawarkan pendekatan yang lebih kuat secara statistik untuk mengidentifikasi zona prioritas dan memahami pola spasial penyakit secara lebih mendalam dibandingkan metode interpolasi berbasis jarak semata.
V. Kesimpulan dan Saran
Berdasarkan seluruh analisis spasial yang telah dilakukan, dapat
disimpulkan bahwa persebaran kasus kusta di Provinsi Jawa Barat tahun
2021 menunjukkan pola tidak acak dan cenderung mengelompok secara
spasial. Hal ini dibuktikan melalui peta deskriptif, uji autokorelasi
spasial global (Moran’s I dan Geary’s C) yang signifikan, serta analisis
LISA yang mengidentifikasi adanya klaster lokal bertipe
High–High pada wilayah tertentu. Pemodelan regresi spasial
menunjukkan bahwa model yang memasukkan komponen ketergantungan spasial
(SEM, SAR, SDM, SDEM, dan SAC) mampu merepresentasikan pola data dengan
lebih baik dibandingkan OLS dan SLX, baik dari sisi kesesuaian model
maupun akurasi prediksi. Temuan ini diperkuat oleh analisis residual
yang menunjukkan bahwa model spasial yang lebih komprehensif
menghasilkan residual yang lebih acak dan bebas dari pola spasial.
Selain itu, hasil interpolasi IDW dan Kriging memperlihatkan gradien
spasial yang konsisten, dengan konsentrasi nilai tinggi di wilayah utara
dan timur Jawa Barat, serta nilai lebih rendah dan homogen di wilayah
selatan.
Berdasarkan temuan tersebut, disarankan agar upaya pengendalian dan pencegahan kusta di Jawa Barat dilakukan dengan pendekatan berbasis wilayah (spatially targeted intervention), khususnya pada daerah yang teridentifikasi sebagai klaster kasus tinggi. Pemerintah daerah dan pemangku kepentingan kesehatan perlu memprioritaskan wilayah hotspot untuk peningkatan skrining aktif, edukasi masyarakat, dan penguatan sistem surveilans. Dari sisi metodologis, penelitian selanjutnya disarankan untuk mengembangkan model spasial dengan memasukkan variabel penjelas yang lebih beragam, seperti faktor sosial-ekonomi, kepadatan penduduk, dan akses layanan kesehatan, serta memanfaatkan data time-series untuk menganalisis dinamika spasial-temporal kasus kusta. Pendekatan ini diharapkan dapat menghasilkan pemahaman yang lebih komprehensif dan mendukung perumusan kebijakan kesehatan masyarakat yang lebih efektif dan tepat sasaran.
Daftar Pustaka
[1] World Health Organization. Leprosy (Hansen’s disease): Fact
Sheet. WHO, 2023.
[2] World Health Organization. Global Leprosy Update 2022. Weekly
Epidemiological Record, 2023.
[3] WHO. Global Leprosy Strategy 2021–2030. Geneva: World Health
Organization, 2021.
[4] Kementerian Kesehatan Republik Indonesia. Profil Kesehatan Indonesia
2022. Jakarta: Kemenkes RI, 2023.
[5] Kementerian Kesehatan Republik Indonesia. Situasi Kusta di
Indonesia. Jakarta: Direktorat Jenderal Pencegahan dan Pengendalian
Penyakit, 2021.
[6] Lawson, A. B. Statistical Methods in Spatial Epidemiology. 2nd
ed. Chichester: Wiley, 2013.
[7] Anselin, L. Spatial Econometrics: Methods and Models. Dordrecht:
Kluwer Academic Publishers, 1988. [8] Waller, L. A., & Gotway, C. A.
Applied Spatial Statistics for Public Health Data. Hoboken: Wiley, 2004.
[9] United Nations. Transforming our World: The 2030 Agenda for
Sustainable Development. New York: United Nations, 2015.
[10] Smith, W. C. S., & Aerts, A. Role of socioeconomic factors in
leprosy control. Leprosy Review, 2014. [11] Anselin, L. Spatial
Econometrics: Methods and Models. Kluwer Academic Publishers,
1988.
[12] LeSage, J., & Pace, R. K. Introduction to Spatial Econometrics.
CRC Press, 2009.
[13] Cliff, A. D., & Ord, J. K. Spatial Processes: Models &
Applications. Pion, 1981.
[14] Anselin, L. “Local Indicators of Spatial Association—LISA.”
Geographical Analysis, 1995.
[15] Elhorst, J. P. Spatial Econometrics: From Cross-Sectional Data to
Spatial Panels. Springer, 2014.
[16] Burrough, P. A., & McDonnell, R. A. Principles of Geographical
Information Systems. Oxford University Press, 1998.
[17] Cressie, N. Statistics for Spatial Data. Wiley, 1993.
[18] Anselin, L. (1988). Spatial Econometrics: Methods and Models.
Dordrecht: Kluwer Academic Publishers.
[19] Cliff, A. D., & Ord, J. K. (1981). Spatial Processes: Models
& Applications. London: Pion.
[20] Anselin, L. (1995). Local Indicators of Spatial Association—LISA.
Geographical Analysis, 27(2), 93–115.
[21] LeSage, J. P., & Pace, R. K. (2009). Introduction to Spatial
Econometrics. CRC Press.
[22] Elhorst, J. P. (2014). Spatial Econometrics. Springer.
[23] Cressie, N. (1993). Statistics for Spatial Data. Wiley.
[24] Bivand, R. S., Pebesma, E., & Gómez-Rubio, V. (2013). Applied
Spatial Data Analysis with R. Springer.
Lampiran (Dashboard) Dashboard Kusta Jawa Barat
Lampiran (Video Presentasi)
Video YouTube
Jika link video tidak bisa dibuka, mohon untuk buka pada tab baru
Lampiran (Syntax R)
# ============================================================
# DASHBOARD INTERAKTIF ANALISIS SPASIAL KUSTA JAWA BARAT 2021
# FINAL UPGRADE — INTERAKTIF + MODEL LENGKAP
# ============================================================
options(shiny.maxRequestSize = 500*1024^2)
library(shiny)
library(bslib)
library(sf)
library(dplyr)
library(stringr)
library(readxl)
library(leaflet)
library(spdep)
library(spatialreg)
library(DT)
library(sp)
library(gstat)
library(terra)
library(plotly)
library(shinycssloaders)
# ============================================================
# THEME
# ============================================================
theme_app <- bs_theme(
version = 5,
bootswatch = "flatly",
primary = "#1565C0",
base_font = font_google("Inter")
)
# ============================================================
# HELPER FUNCTIONS
# ============================================================
clean_kabkota <- function(x){
x <- tolower(as.character(x))
x <- str_replace_all(x, "[^a-z\\s]", " ")
x <- str_squish(x)
x <- str_replace(x, "^kabupaten\\s+", "")
x <- str_replace(x, "^kota\\s+", "") # ⬅️ INI KUNCI CIMAHI
x
}
safe_num <- function(x){
suppressWarnings(as.numeric(x))
}
make_listw <- function(sfobj){
nb <- poly2nb(sfobj, queen = TRUE)
nb2listw(nb, style="W", zero.policy=TRUE)
}
add_slx_terms <- function(df, lw){
Wx <- as.data.frame(lapply(df[,c("x1","x2","x3","x4","x5")],
function(v) lag.listw(lw, v, zero.policy=TRUE)))
names(Wx) <- paste0("W_", names(Wx))
cbind(df, Wx)
}
safe_fit <- function(expr){
tryCatch(expr, error = function(e) e)
}
is_err <- function(x) inherits(x, "error")
rmse <- function(y, yhat) {
sqrt(mean((y - yhat)^2, na.rm = TRUE))
}
mae <- function(y, yhat) {
mean(abs(y - yhat), na.rm = TRUE)
}
# ============================================================
# UI
# ============================================================
ui <- fluidPage(
theme = theme_app,
titlePanel("🧭 Dashboard Interaktif Analisis Spasial KUSTA Jawa Barat 2021"),
tags$head(
tags$style(HTML("
#box_model pre {
white-space: pre !important;
overflow-x: auto !important;
font-family: ui-monospace, SFMono-Regular, Menlo, Monaco, Consolas,
'Liberation Mono','Courier New', monospace;
font-size: 12px;
line-height: 1.25;
margin: 0;
}
"))
),
sidebarLayout(
sidebarPanel(
fileInput("shp", "Upload Shapefile (gadm41_IDN_2.zip)", accept=".zip"),
fileInput("xls", "Upload Excel Kusta 2021", accept=".xlsx"),
hr(),
selectInput(
"menu", "Pilih Analisis",
choices = c(
"Eksplorasi – Peta Kasus",
"Eksplorasi – Statistik Deskriptif",
"Autokorelasi – Moran & Geary",
"Autokorelasi – LISA",
"Interpolasi – IDW"
)
),
numericInput("alpha", "Alpha LISA", 0.05),
sliderInput("idw_p", "Power IDW", 1, 5, 2, step = 0.5),
selectInput(
"interp_method",
"Metode Interpolasi",
choices = c(
"IDW" = "idw",
"Kriging" = "kriging"
)
),
actionButton("run", "Jalankan Analisis", class="btn-primary")
),
mainPanel(
tabsetPanel(
tabPanel(
"🗺️ Peta Interaktif",
leafletOutput("map", height = 650)
),
tabPanel(
"📊 Tabel & Statistik",
DTOutput("table")
),
tabPanel(
"🧠 Output Analisis",
div(
style = "
height: 500px;
overflow-y: auto;
background-color: #f8f9fa;
padding: 12px;
border-radius: 6px;
border: 1px solid #ddd;
font-size: 13px;
white-space: pre;
",
verbatimTextOutput("model_txt")
)
),
tabPanel(
"📊 Ukuran Epidemiologi",
fluidRow(
column(
12,
h4("📌 Incidence Rate Kusta per 100.000 Penduduk"),
p("Incidence rate dihitung sebagai jumlah kasus kusta dibagi jumlah penduduk, dikalikan 100.000.")
)
),
hr(),
fluidRow(
column(
6,
leafletOutput("map_ir", height = 450)
),
column(
6,
DTOutput("table_ir")
)
),
hr(),
fluidRow(
column(
12,
h4("📈 Ringkasan Statistik Incidence Rate"),
verbatimTextOutput("ir_summary")
)
)
),
tabPanel(
"📐 Model Spasial",
## =========================
## 1️⃣ KONTROL
## =========================
fluidRow(
column(
6,
selectInput(
"model_choice",
"Pilih Model Spasial",
choices = c("OLS","SAR","SEM","SLX","SDM","SDEM","SAC")
)
),
column(
6,
radioButtons(
"map_type",
"Jenis Peta",
choices = c(
"Prediksi Model" = "fitted",
"Residual Model" = "residual"
),
inline = TRUE
)
)
),
br(),
## =========================
## 2️⃣ PETA (KIRI) + MODEL (KANAN)
## =========================
fluidRow(
column(
width = 6,
leafletOutput("map_model", height = 520)
),
column(
width = 6,
div(
id = "box_model",
style = "
height:520px;
overflow:auto;
background:#f8f9fa;
padding:12px;
border:1px solid #ddd;
border-radius:6px;
white-space: pre;
font-family: ui-monospace, SFMono-Regular, Menlo, Monaco, Consolas,
'Liberation Mono','Courier New', monospace;
font-size:12px;
",
verbatimTextOutput("model_spatial_txt")
)
)
)
),
tabPanel(
"📊 Perbandingan Model",
h4("Perbandingan Model Spasial"),
hr(),
# =========================
# 1️⃣ TABEL PERBANDINGAN
# =========================
div(
style = "
height:300px;
overflow:auto;
background:#f8f9fa;
padding:12px;
border:1px solid #ddd;
border-radius:6px;
white-space: pre;
",
verbatimTextOutput("model_compare_txt")
),
br(),
h4("🏆 Model Spasial Terbaik"),
hr(),
# =========================
# 2️⃣ PETA + RINGKASAN MODEL TERBAIK
# =========================
fluidRow(
column(
6,
leafletOutput("map_best_model", height = 450)
),
column(
6,
div(
style = "
height:450px;
overflow:auto;
background:#f8f9fa;
padding:12px;
border:1px solid #ddd;
border-radius:6px;
white-space: pre;
",
verbatimTextOutput("best_model_txt")
)
)
)
)
)
)
)
)
# ============================================================
# SERVER
# ============================================================
server <- function(input, output){
# ---------------- SHAPEFILE ----------------
shp_sf <- eventReactive(input$run, {
req(input$shp)
tmp <- tempdir()
unzip(input$shp$datapath, exdir = tmp)
shpfile <- list.files(tmp, "\\.shp$", full.names = TRUE)[1]
sf <- st_read(shpfile, quiet = TRUE)
sf <- sf %>%
filter(
NAME_1 == "Jawa Barat",
!NAME_2 %in% c("Waduk Cirata")
)
sf$key <- clean_kabkota(sf$NAME_2)
sf
})
# ---------------- EXCEL ----------------
excel_df <- eventReactive(input$run, {
req(input$xls)
df <- read_excel(input$xls$datapath)
df <- df %>%
rename(
Kabupaten_Kota = nama_kabupaten_kota,
y = jumlah_kasus,
x1 = jumlah_penduduk,
x2 = persentase_sanitasi_layak,
x3 = kepadatan_penduduk,
x4 = rata_rata_lama_sekolah,
x5 = persentase_penduduk_miskin
)
for(v in c("y","x1","x2","x3","x4","x5"))
df[[v]] <- safe_num(df[[v]])
df$key <- clean_kabkota(df$Kabupaten_Kota)
df
})
# ---------------- JOIN ----------------
data_sf <- reactive({
req(shp_sf(), excel_df())
sf_join <- left_join(shp_sf(), excel_df(), by = "key")
sf_join$y[is.na(sf_join$y)] <- 0
# 🔴 DISSOLVE GEOMETRY
sf_final <- sf_join %>%
group_by(key, NAME_2) %>%
summarise(
y = mean(y, na.rm = TRUE),
x1 = mean(x1, na.rm = TRUE),
x2 = mean(x2, na.rm = TRUE),
x3 = mean(x3, na.rm = TRUE),
x4 = mean(x4, na.rm = TRUE),
x5 = mean(x5, na.rm = TRUE),
incidence_rate = ifelse(
x1 > 0,
(y / x1) * 100000,
NA_real_
),
geometry = st_union(geometry),
.groups = "drop"
)
sf_final
})
# ============================================================
# INTERPOLASI (IDW + KRIGING) - WAJIB CRS METER (UTM)
# TARUH: tepat setelah data_sf <- reactive({ ... })
# ============================================================
utm_crs <- 32748 # UTM 48S (cukup aman untuk Jabar)
pts_utm <- reactive({
req(data_sf())
st_centroid(data_sf()) |>
st_transform(utm_crs) |>
mutate(z = as.numeric(y))
})
grid_utm <- reactive({
req(data_sf())
poly_utm <- st_transform(st_union(data_sf()), utm_crs)
# Grid lebih ringan: 5km. Kalau mau halus: 2km (lebih berat)
grid <- st_make_grid(poly_utm, cellsize = 5000, what = "centers")
grid_sf <- st_sf(geometry = grid)
# MASK: hanya titik grid yang berada di dalam polygon Jabar
st_intersection(grid_sf, poly_utm)
})
idw_utm <- reactive({
req(pts_utm(), grid_utm())
gstat::idw(
z ~ 1,
locations = as(pts_utm(), "Spatial"),
newdata = as(grid_utm(), "Spatial"),
idp = input$idw_p
)
})
kriging_utm <- reactive({
req(pts_utm(), grid_utm())
pts_sp <- as(pts_utm(), "Spatial")
grid_sp <- as(grid_utm(), "Spatial")
vgm_emp <- variogram(z ~ 1, pts_sp)
vgm_fit <- fit.variogram(vgm_emp, vgm("Sph"))
krige(z ~ 1, locations = pts_sp, newdata = grid_sp, model = vgm_fit)
})
data_desc <- reactive({
data_sf() %>%
st_drop_geometry() %>%
select(
NAME_2, # nama kab/kota
y, x1, x2, x3, x4, x5
)
})
stat_desc <- reactive({
df <- data_desc()
vars <- c("y", "x1", "x2", "x3", "x4", "x5")
var_names <- c(
"Kasus Kusta",
"Jumlah Penduduk",
"Sanitasi Layak (%)",
"Kepadatan Penduduk",
"Rata-rata Lama Sekolah (tahun)",
"Penduduk Miskin (%)"
)
result <- lapply(seq_along(vars), function(i){
v <- vars[i]
x <- df[[v]]
min_val <- min(x, na.rm = TRUE)
max_val <- max(x, na.rm = TRUE)
min_loc <- df$NAME_2[which.min(x)][1]
max_loc <- df$NAME_2[which.max(x)][1]
data.frame(
Variabel = var_names[i],
Minimum = round(min_val, 2),
Lokasi_Minimum = min_loc,
Maksimum = round(max_val, 2),
Lokasi_Maksimum = max_loc,
Rata_rata = round(mean(x, na.rm = TRUE), 2),
Median = round(median(x, na.rm = TRUE), 2),
Simpangan_Baku = round(sd(x, na.rm = TRUE), 2)
)
})
do.call(rbind, result)
})
# ---------------- LISTW ----------------
listw <- reactive({
make_listw(data_sf())
})
moran_geary <- reactive({
req(data_sf(), listw())
df <- st_drop_geometry(data_sf())
lw <- listw()
list(
moran = moran.test(df$y, lw, zero.policy = TRUE),
geary = geary.test(df$y, lw, zero.policy = TRUE)
)
})
lisa_result <- reactive({
req(data_sf(), listw())
df <- st_drop_geometry(data_sf())
lw <- listw()
lm <- localmoran(df$y, lw, zero.policy = TRUE)
# localmoran biasanya punya 5 kolom:
# Ii, E.Ii, Var.Ii, Z.Ii, Pr(z > 0)
# Kita ambil by index agar tidak tergantung nama kolom
out <- data.frame(
NAME_2 = df$NAME_2,
Ii = lm[, 1],
E_Ii = lm[, 2],
Var_Ii = lm[, 3],
Z_Ii = lm[, 4],
P_Ii = lm[, 5]
)
out
})
# ============================================================
# DATA TITIK UNTUK INTERPOLASI (CENTROID)
# ============================================================
data_point <- reactive({
req(data_sf())
st_centroid(data_sf()) |>
st_transform(4326)
})
grid_pixel <- reactive({
req(data_point())
# buat grid pixel (bukan titik)
grd <- st_make_grid(
data_point(),
cellsize = 0.15, # ukuran pixel (atur nanti)
what = "polygons"
)
st_sf(geometry = grd)
})
grid_sp <- reactive({
req(grid_pixel())
as(grid_pixel(), "Spatial")
})
points_sp <- reactive({
req(data_point())
as(data_point(), "Spatial")
})
idw_surface <- reactive({
req(points_sp(), grid_sp())
gstat::idw(
formula = y ~ 1,
locations = points_sp(),
newdata = grid_sp(),
idp = input$idw_p
)
})
kriging_surface <- reactive({
req(points_sp(), grid_sp())
vgm_emp <- variogram(y ~ 1, points_sp())
vgm_fit <- fit.variogram(vgm_emp, vgm("Sph"))
krige(
y ~ 1,
locations = points_sp(),
newdata = grid_sp(),
model = vgm_fit
)
})
interp_raster <- reactive({
req(input$interp_method)
surf <- if (input$interp_method == "idw") {
idw_surface()
} else {
kriging_surface()
}
raster::raster(surf["var1.pred"])
})
interp_masked <- reactive({
req(interp_raster(), data_sf())
poly_sp <- as(st_union(data_sf()), "Spatial")
raster::mask(interp_raster(), poly_sp)
})
# ============================================================
# GRID INTERPOLASI
# ============================================================
grid_interp <- reactive({
req(data_point())
grid <- st_make_grid(
data_point(),
n = c(60, 60),
what = "centers"
)
st_sf(geometry = grid)
})
# ============================================================
# GRID INTERPOLASI (DI-POTONG POLYGON JAWA BARAT)
# ============================================================
grid_masked <- reactive({
req(grid_interp(), data_sf())
st_intersection(
grid_interp(),
st_union(data_sf())
)
})
# ============================================================
# INTERPOLASI IDW (SF + GSTAT)
# ============================================================
idw_result <- reactive({
req(data_point(), grid_interp())
gstat::idw(
y ~ 1,
locations = data_point(),
newdata = grid_interp(),
idp = input$idw_p
)
})
# ============================================================
# INTERPOLASI KRIGING
# ============================================================
kriging_result <- reactive({
req(data_point(), grid_interp())
vgm_emp <- variogram(y ~ 1, data_point())
vgm_fit <- fit.variogram(vgm_emp, vgm("Sph"))
krige(
y ~ 1,
locations = data_point(),
newdata = grid_interp(),
model = vgm_fit
)
})
# ============================================================
# MODEL FIT (SEMUA SEKALIGUS, SEKALI HITUNG)
# ============================================================
model_pack <- reactive({
df <- st_drop_geometry(data_sf())
lw <- listw()
f <- y ~ x1 + x2 + x3 + x4 + x5
models <- list(
OLS = lm(f, df),
SAR = lagsarlm(f, df, lw),
SEM = errorsarlm(f, df, lw),
SLX = lm(y ~ x1+x2+x3+x4+x5 + W_x1+W_x2+W_x3+W_x4+W_x5, add_slx_terms(df,lw)),
SDM = lagsarlm(f, df, lw, type="mixed"),
SDEM = errorsarlm(f, df, lw, Durbin=TRUE),
SAC = sacsarlm(f, df, lw)
)
metrics <- lapply(models, function(m){
res <- residuals(m)
data.frame(
AIC = AIC(m),
RMSE = sqrt(mean(res^2)),
MAE = mean(abs(res))
)
})
comp <- do.call(rbind, metrics)
comp$Model <- rownames(comp)
comp <- comp[order(comp$AIC), ]
list(models=models, comp=comp, best=comp$Model[1])
})
# ============================================================
# OUTPUT MODEL TERBAIK (TEKS)
# ============================================================
output$model_compare_txt <- renderPrint({
comp <- model_pack()$comp
print(comp, row.names = FALSE)
})
output$best_model_txt <- renderPrint({
best <- model_pack()$best
cat("MODEL TERBAIK:", best, "\n")
cat("========================\n\n")
print(summary(model_pack()$models[[best]]))
})
# ============================================================
# MAP MODEL: PREDIKSI & RESIDUAL
# ============================================================
output$map_model <- renderLeaflet({
req(data_sf(), model_pack(), input$model_choice, input$map_type)
sfm <- st_transform(data_sf(), 4326)
model <- model_pack()$models[[input$model_choice]]
if (inherits(model, "error")) {
return(leaflet() %>% addProviderTiles("CartoDB.Positron"))
}
# =========================
# AMBIL FITTED / RESIDUAL
# =========================
if (input$map_type == "fitted") {
sfm$value <- as.numeric(fitted(model))
title_legend <- "Prediksi Kasus Kusta"
pal <- colorNumeric("viridis", sfm$value, na.color = "#f0f0f0")
} else {
sfm$value <- as.numeric(residuals(model))
title_legend <- "Residual Model"
pal <- colorNumeric("RdBu", sfm$value, reverse = TRUE, na.color = "#f0f0f0")
}
# =========================
# LEAFLET MAP
# =========================
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(value),
fillOpacity = 0.85,
color = "#333333",
weight = 1,
popup = ~paste0(
"<b>", NAME_2, "</b><br/>",
title_legend, ": ", round(value, 2)
),
highlightOptions = highlightOptions(
weight = 3,
color = "#1565C0",
bringToFront = TRUE
)
) %>%
addLegend(
pal = pal,
values = ~value,
title = title_legend,
position = "bottomright"
)
})
# ============================================================
# MAP INCIDENCE RATE
# ============================================================
output$map_ir <- renderLeaflet({
req(data_sf())
sfm <- st_transform(data_sf(), 4326)
pal <- colorNumeric(
"YlOrRd",
sfm$incidence_rate,
na.color = "#f0f0f0"
)
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(incidence_rate),
fillOpacity = 0.85,
color = "#333",
popup = ~paste0(
"<b>", NAME_2, "</b><br/>",
"IR: ", round(incidence_rate, 2),
" per 100.000"
)
) %>%
addLegend(
pal = pal,
values = ~incidence_rate,
title = "Incidence Rate /100.000",
position = "bottomright"
)
})
# ============================================================
# OUTPUT MODEL SPASIAL (TAB 📐 MODEL SPASIAL)
# ============================================================
output$model_spatial_txt <- renderPrint({
req(input$model_choice, model_pack())
models <- model_pack()$models
mname <- input$model_choice
model <- models[[mname]]
cat("MODEL SPASIAL:", mname, "\n")
cat("====================================\n\n")
if (inherits(model, "error")) {
cat("⚠️ Model gagal diestimasi\n\n")
cat("Pesan error:\n")
cat(model$message)
} else {
print(summary(model))
}
})
# ============================================================
# LEAFLET MAP
# ============================================================
output$map <- renderLeaflet({
req(data_sf(), input$menu)
sfm <- st_transform(data_sf(), 4326)
# ============================================================
# PETA MODEL TERBAIK
# ============================================================
output$map_best_model <- renderLeaflet({
best <- model_pack()$best
model <- model_pack()$models[[best]]
sfm <- data_sf()
sfm$val <- fitted(model)
pal <- colorNumeric("viridis", sfm$val)
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(val),
fillOpacity = 0.85,
color = "#333",
popup = ~paste0(NAME_2, "<br>Prediksi: ", round(val,2))
) %>%
addLegend(pal = pal, values = ~val, title = "Prediksi Model Terbaik")
})
# =====================================================
# 1️⃣ EKSPLORASI – PETA KASUS (KODE KAMU, TETAP)
# =====================================================
if (input$menu == "Eksplorasi – Peta Kasus") {
pal <- colorNumeric(
palette = "viridis",
domain = sfm$y,
na.color = "#f0f0f0"
)
popup_info <- sprintf(
"<div style='font-size:14px; line-height:1.6'>
<h4 style='margin-bottom:6px; color:#1565C0'>%s</h4>
<b>🦠 Kasus Kusta:</b> %s<br/>
<hr style='margin:6px 0'/>
<b>👥 Jumlah Penduduk:</b> %s<br/>
<b>🚿 Sanitasi Layak (%%):</b> %s<br/>
<b>🏘️ Kepadatan Penduduk:</b> %s<br/>
<b>🎓 Rata-rata Lama Sekolah (tahun):</b> %s<br/>
<b>💸 Penduduk Miskin (%%):</b> %s
</div>",
sfm$NAME_2,
formatC(sfm$y, format="d", big.mark="."),
formatC(sfm$x1, format="f", digits=2, big.mark=".", decimal.mark=","),
formatC(sfm$x2, format="f", digits=2, decimal.mark=","),
formatC(sfm$x3, format="f", digits=0, big.mark="."),
formatC(sfm$x4, format="f", digits=2, decimal.mark=","),
formatC(sfm$x5, format="f", digits=2, decimal.mark=",")
)
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(y),
fillOpacity = 0.85,
color = "#333333",
weight = 1,
label = ~paste0("<b>", NAME_2, "</b><br/>Kasus: ", y) %>%
lapply(htmltools::HTML),
popup = lapply(popup_info, htmltools::HTML),
highlightOptions = highlightOptions(
weight = 3,
color = "#1565C0",
bringToFront = TRUE
)
) %>%
addLegend(
pal = pal,
values = ~y,
title = "Kasus Kusta",
position = "bottomright"
)
}
# =====================================================
# 2️⃣ AUTOKORELASI – MORAN & GEARY (BARU)
# =====================================================
else if (input$menu == "Autokorelasi – Moran & Geary") {
lw <- listw()
sfm$lag_y <- lag.listw(lw, sfm$y, zero.policy = TRUE)
pal <- colorNumeric(
palette = "RdYlBu",
domain = sfm$lag_y,
na.color = "#f0f0f0"
)
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(lag_y),
fillOpacity = 0.85,
color = "#333333",
weight = 1,
popup = ~paste0(
"<b>", NAME_2, "</b><br/>",
"Kasus: ", y, "<br/>",
"Spatial Lag: ", round(lag_y, 2)
),
highlightOptions = highlightOptions(
weight = 3,
color = "#D32F2F",
bringToFront = TRUE
)
) %>%
addLegend(
pal = pal,
values = ~lag_y,
title = "Spatial Lag Kasus Kusta",
position = "bottomright"
)
}
# =====================================================
# 3️⃣ AUTOKORELASI – LISA (LOCAL MORAN) → PETA SAJA
# =====================================================
else if (input$menu == "Autokorelasi – LISA") {
lw <- listw()
df <- st_drop_geometry(data_sf())
z_y <- as.numeric(scale(df$y))
lag_z <- lag.listw(lw, z_y, zero.policy = TRUE)
lisa <- localmoran(df$y, lw, zero.policy = TRUE)
sfm$cluster <- "Not Significant"
alpha <- input$alpha
sfm$cluster[z_y > 0 & lag_z > 0 & lisa[,5] < alpha] <- "High-High"
sfm$cluster[z_y < 0 & lag_z < 0 & lisa[,5] < alpha] <- "Low-Low"
sfm$cluster[z_y > 0 & lag_z < 0 & lisa[,5] < alpha] <- "High-Low"
sfm$cluster[z_y < 0 & lag_z > 0 & lisa[,5] < alpha] <- "Low-High"
pal <- colorFactor(
palette = c(
"High-High" = "#D7191C",
"Low-Low" = "#2C7BB6",
"High-Low" = "#FDAE61",
"Low-High" = "#ABD9E9",
"Not Significant" = "#EEEEEE"
),
levels = c(
"High-High","Low-Low","High-Low","Low-High","Not Significant"
)
)
leaflet(sfm) %>%
addProviderTiles("CartoDB.Positron") %>%
addPolygons(
fillColor = ~pal(cluster),
fillOpacity = 0.85,
color = "#333333",
weight = 1,
popup = ~paste0("<b>", NAME_2, "</b><br/>LISA: ", cluster),
highlightOptions = highlightOptions(
weight = 3,
bringToFront = TRUE
)
) %>%
addLegend(
pal = pal,
values = ~cluster,
title = "LISA Cluster",
position = "bottomright"
)
}
else if (input$menu == "Interpolasi – IDW") {
req(input$interp_method)
# --- pilih hasil ---
res <- NULL
title_map <- NULL
pred <- NULL
if (input$interp_method == "idw") {
res <- idw_utm()
pred <- res$var1.pred
title_map <- paste0("Interpolasi IDW (p=", input$idw_p, ")")
} else {
res <- kriging_utm()
pred <- res$var1.pred
title_map <- "Interpolasi Kriging"
}
# --- ubah hasil (sp) -> sf -> lonlat untuk leaflet ---
res_df <- as.data.frame(res)
xy <- sp::coordinates(res)
res_df$X <- xy[,1]
res_df$Y <- xy[,2]
res_sf <- st_as_sf(res_df, coords = c("X","Y"), crs = utm_crs) |>
st_transform(4326)
lonlat <- st_coordinates(res_sf)
res_sf$lon <- lonlat[,1]
res_sf$lat <- lonlat[,2]
pal <- colorNumeric("YlOrRd", domain = pred, na.color = "#f0f0f0")
leaflet(res_sf) %>%
addProviderTiles("CartoDB.Positron") %>%
addCircleMarkers(
lng = ~lon,
lat = ~lat,
radius = 4,
stroke = FALSE,
fillOpacity = 0.75,
fillColor = ~pal(pred),
popup = ~paste0(title_map, "<br>Pred: ", round(pred, 2))
) %>%
addLegend(
pal = pal,
values = ~pred,
title = title_map,
position = "bottomright"
)
}
})
# ============================================================
# TABLE
# ============================================================
output$table <- renderDT({
req(stat_desc())
datatable(
stat_desc(),
rownames = FALSE,
options = list(
pageLength = 6,
dom = "tip",
autoWidth = TRUE
)
)
})
# ============================================================
# TABLE INCIDENCE RATE
# ============================================================
output$table_ir <- renderDT({
req(data_sf())
df <- data_sf() |>
st_drop_geometry() |>
select(
Kabupaten_Kota = NAME_2,
Kasus = y,
Penduduk = x1,
Incidence_Rate = incidence_rate
) |>
mutate(
Incidence_Rate = round(Incidence_Rate, 2)
)
datatable(
df,
rownames = FALSE,
options = list(
pageLength = 10,
autoWidth = TRUE
)
)
})
# ============================================================
# PLOTLY GRAPH
# ============================================================
output$plot <- renderPlotly({
req(data_sf())
df <- st_drop_geometry(data_sf())
p <- ggplot(df, aes(x = x5, y = y)) +
geom_point(color = "#1565C0", size = 3) +
geom_smooth(method = "lm", se = FALSE, color = "red") +
theme_minimal() +
labs(x = "% Penduduk Miskin", y = "Kasus Kusta",
title = "Hubungan Kemiskinan dan Kasus Kusta")
ggplotly(p)
})
# ============================================================
# MODEL OUTPUT
# ============================================================
output$model_txt <- renderPrint({
req(input$menu)
if (input$menu == "Autokorelasi – Moran & Geary") {
mg <- moran_geary()
cat("AUTOKORELASI SPASIAL GLOBAL\n")
cat("----------------------------------\n")
cat("Moran's I:\n")
print(mg$moran)
cat("\nGeary's C:\n")
print(mg$geary)
} else if (input$menu == "Autokorelasi – LISA") {
lisa <- lisa_result()
alpha <- input$alpha
cat("LOCAL INDICATORS OF SPATIAL ASSOCIATION (LISA)\n")
cat("===============================================\n\n")
cat("Alpha =", alpha, "\n\n")
sig <- lisa$P_Ii < alpha # ✅ INI YANG BENAR
cat("Jumlah wilayah signifikan:\n")
print(table(Signifikan = sig))
cat("\nRingkasan p-value:\n")
print(summary(lisa$P_Ii))
}
})
# ============================================================
# SUMMARY INCIDENCE RATE
# ============================================================
output$ir_summary <- renderPrint({
req(data_sf())
ir <- data_sf()$incidence_rate
cat("INCIDENCE RATE KUSTA\n")
cat("=========================\n\n")
print(summary(ir))
cat("\nWilayah IR tertinggi:\n")
idx_max <- which.max(ir)
cat(data_sf()$NAME_2[idx_max],
"(", round(ir[idx_max], 2), ")\n")
cat("\nWilayah IR terendah:\n")
idx_min <- which.min(ir)
cat(data_sf()$NAME_2[idx_min],
"(", round(ir[idx_min], 2), ")\n")
})
# ============================================================
# DEBUG
# ============================================================
output$debug_txt <- renderPrint({
cat("DEBUG INFO\n")
cat("===========================\n\n")
cat("Jumlah baris (row):\n")
print(nrow(data_sf()))
cat("\nJumlah kab/kota unik:\n")
print(length(unique(data_sf()$NAME_2)))
cat("\nDaftar kab/kota:\n")
print(sort(unique(data_sf()$NAME_2)))
cat("\nKolom data:\n")
print(names(data_sf()))
})
}
shinyApp(ui, server)