1 Abstrak

Tuberkulosis (TB) merupakan penyakit menular yang disebabkan oleh bak-teri Mycobacterium tuberculosis dan masih menjadi salah satu masalah kesehatan masyarakat utama di dunia. Menurut World Health Organization (WHO), TB termasuk penyebab kematian tertinggi setelah HIV/AIDS dan menjadi salah satu fokus dalam Tujuan Pembangunan Berkelanjutan (Sustainable Development Goals). Indonesia meru-pakan salah satu negara dengan beban TB yang tinggi, sehingga diperlukan upaya pengendalian yang efektif melalui identifikasi wilayah berisiko tinggi. Analisis spasial dan temporal diperlukan untuk memahami variasi risiko TB antarwilayah dan antarwaktu guna mendukung perencanaan program pengendalian yang lebih tepat sasaran. Penelitian ini bertujuan untuk memetakan risiko relatif insidensi TB di Provinsi Jawa Barat selama periode 2020–2025 menggunakan model Bayesian spatio-temporal. Data yang digunakan berupa jumlah kasus insidensi TB dan jumlah penduduk pada 27 kabupaten/kota di Jawa Barat tahun 2020–2025. Pemodelan dilakukan menggunakan model hierarkis Bayesian dengan distribusi Poisson, yang mengakomodasi efek spasial, temporal, dan interaksi spasial-temporal. Estimasi parameter dilakukan menggunakan metode Integrated Nested Laplace Approximation (INLA) untuk memperoleh estimasi yang efisien dan akurat. Hasil penelitian menunjukkan adanya variasi risiko relatif TB yang signifikan antarwilayah dan antarperiode waktu. Pada tahun 2020, wilayah dengan risiko relatif tinggi terkonsentrasi di Kota Cirebon, Kota Tasikmalaya, dan Kota Bekasi. Risiko relatif meningkat pada tahun 2021–2022 dengan penyebaran yang lebih luas terutama di kawasan timur dan utara Jawa Barat, meliputi Kabupaten Karawang, Kabupaten Bekasi, Kabupaten Cirebon, Kota Bekasi, dan Kota Cirebon. Tahun 2022 merupakan periode dengan tingkat risiko relatif tertinggi. Selanjutnya, pada periode 2023–2025 terjadi kecenderungan penurunan risiko relatif, meskipun Kota Cirebon dan Kota Banjar tetap menunjukkan risiko tinggi secara konsisten. Sebaliknya, Kabupaten Bogor, Kabupaten Sukabumi, Kabupaten Cianjur, Kabupaten Garut, dan Kabupaten Tasikmalaya cenderung memiliki risiko relatif rendah selama periode penelitian. Hasil pemetaan juga menunjukkan adanya pola pengelompokan spasial (spatial clustering), terutama di wilayah timur dan utara Jawa Barat. Temuan ini mengindikasikan bahwa faktor spasial dan temporal berperan penting dalam persebaran kasus TB. Oleh karena itu, model Bayesian spatio-temporal mampu memberikan gambaran risiko relatif yang lebih komprehensif dan dapat menjadi dasar dalam penyusunan kebijakan pengendalian TB yang lebih efektif dan berbasis wilayah.

Keywords: Tuberkulosis, Risiko Relatif, Bayesian Spatio-Temporal, INLA

2 Pendahuluan

Tuberkulosis (TB) merupakan penyakit menular yang disebabkan oleh bakteri Mycobacterium tuberculosis yang menyebar dari penderita TB melalui udara. Sebagian besar bakteri tersebut menyerang organ paru, tetapi dapat juga menyerang organ lain (ekstra paru). Hampir seperempat penduduk dunia terinfeksi dengan bakteri tersebut. Sekitar 90% TB diderita oleh orang dewasa (55% laki-laki dan 33% perempuan) dan 12% diderita oleh anak-anak. Hingga saat ini, TB masih menjadi penyebab kematian tertinggi setelah HIV/AIDS dan termasuk 20 penyebab kematian utama dunia (WHO, 2024).

TB masih menjadi masalah kesehatan baik di tingkat nasional maupun global, serta termasuk salah satu prioritas dalam Tujuan Pembangunan Berkelanjutan (SDGs), khu-susnya tujuan tiga (3), yaitu kehidupan sehat dan sejahtera. Seluruh negara anggota World Health Organization (WHO) dan Perserikatan Bangsa-Bangsa (PBB) berkomit-men untuk mengakhiri epidemi TB melalui End TB Strategy (Kemenkes RI, 2024a). Menurut WHO, salah satu indikator utama yang digunakan untuk menilai beban pen-yakit TB adalah insidensi, yaitu jumlah kasus TB baru dan kambuh dalam 100.000 penduduk pada waktu tertentu (Kemenkes RI, 2024a). Dalam mendukung target SDGs 2030 dan End TB Strategy, WHO menetapkan target penurunan angka kematian akibat TB sebesar 75% dan penurunan insidensi TB sebesar 50% pada tahun 2025, serta penurunan angka kematian sebesar 90% dan penurunan insidensi sebesar 80% pada tahun 2030.

Berdasarkan Peta Jalan Eliminasi TB di Indonesia yang sesuai dengan target global, Indonesia berkomitmen menurunkan angka insidensi TB menjadi 163 kasus per 100.000 penduduk pada tahun 2025. Selanjutnya, pada tahun 2030 ditargetkan turun menjadi 65 kasus per 100.000 penduduk, dengan angka kematian akibat TB menjadi 6 per 100.000 penduduk (Kemenkes RI, 2023b). Di tingkat nasional, melalui Rencana Pembangunan Jangka Menengah Nasional (RPJMN) 2025–2029, Indonesia menargetkan penurunan angka insidensi TB berdasarkan baseline tahun 2024 sebesar 387 kasus per 100.000 penduduk, menjadi 329 pada tahun 2025, dan 190 pada tahun 2029 (Bappenas, 2025).

Berdasarkan Global Tuberculosis Report 2024 (WHO, 2024), Indonesia menempati peringkat kedua di dunia dalam hal beban kasus TB setelah India, dengan estimasi 1.090.000 kasus dan 125.000 kematian setiap tahunnya atau setara dengan sekitar 14 kematian per jam. Jumlah kasus TB yang ditemukan terus mengalami peningkatan. Pada tahun 2024, tercatat sekitar 885.000 kasus, naik dari 821.200 kasus pada tahun 2023. Peningkatan ini terutama terjadi di provinsi dengan jumlah penduduk besar, seperti Jawa Barat, Jawa Timur, dan Jawa Tengah. Secara khusus, Provinsi Jawa Barat mencatat jumlah kasus TB tertinggi di Indonesia. Pada tahun 2024, jumlah kasus yang dilaporkan mencapai 229.683 kasus, mengalami kenaikan sebesar 7,77% dibandingkan tahun sebe-lumnya.

Kondisi ini mendorong program pengendalian tuberkulosis nasional untuk terus mempercepat upaya eliminasi TB pada tahun 2030 sesuai target yang telah ditetapkan. Untuk menggambarkan perkembangan kasus TB, berikut disajikan tren insidensi TB di Provinsi Jawa Barat selama tahun 2020–2024.

Gambar 1.1 Insidensi TB di Provinsi Jawa Barat Tahun 2020 – 2025
Gambar 1.1 Insidensi TB di Provinsi Jawa Barat Tahun 2020 – 2025

Berdasarkan Gambar 1.1, insidensi TB di Provinsi Jawa Barat mengalami pening-katan dari tahun ke tahun, dengan lonjakan signifikan pada tahun 2022. Tren kenaikan ini diperkirakan akan terus berlanjut di tahun berikutnya. Oleh karena itu, untuk mendukung upaya percepatan eliminasi TB, diperlukan identifikasi faktor- faktor yang berpengaruh signifikan terhadap TB melalui sebuah pemodelan, guna menekan angka insidensi TB di Provinsi Jawa Barat. Penelitian menunjukkan bahwa infeksi HIV meru-pakan faktor risiko yang berkontribusi terhadap peningkatan kasus baru TB (WHO, 2024). Kepadatan penduduk dan penduduk miskin juga berperan dalam meningkatkan penyebaran TB (Aulia et al., 2024; Kustanto, 2020). Sedangkan, penerapan perilaku hidup bersih dan sehat, serta akses terhadap air minum layak dan sanitasi layak diketahui mengurangi risiko terinfeksi TB dan menurunkan kasus TB (Chang & Renaud, 2022; Hasanuddin et al., 2023; Islam et al., 2025; Sukmawati & Sri Hazanah, 2025; Widyaningsih & Budiawan, 2023).

Analisis regresi merupakan metode statistika yang bertujuan untuk mengestimasi dan memodelkan hubungan linier antara variabel yang diamati dengan sejumlah faktor yang diduga memengaruhinya (Draper & Smith, 1998). Namun pendekatan regresi konven-sional tidak dapat menangkap pola penyebaran TB yang dipengaruhi oleh kedekatan geografis antarwilayah. Penelitian menunjukkan bahwa insidensi TB di suatu wilayah dipengaruhi oleh wilayah sekitarnya, menandakan adanya korelasi spasial antarwilayah yang berdekatan (dependensi spasial) yang mengindikasikan adanya rantai penularan aktif TB secara geografis (Ng, In Chan et al., 2012).

Dependensi spasial berarti observasi pada suatu lokasi dipengaruhi oleh lokasi sekitarnya. Hal ini menyebabkan asumsi independensi antar pengamatan dalam analisis regresi klasik tidak terpenuhi (Ng, In Chan et al., 2012). Sebagi contoh, Kota Cimahi dan Kota Bandung, yang berdekatan secara geografis sama- sama menunjukkan insidensi TB yang tinggi, hal ini sejalan dengan prinsip Waldo Tobler bahwa lokasi yang berdekatan cenderung memiliki pola yang lebih mirip dibandingkan lokasi yang berjauhan. Salah satu bentuk ketergantungan spasial ini tercermin dalam mobilitas harian penduduk Cimahi menuju Kota Bandung, yang turut mempercepat penyebaran TB karena ting-ginya kemungkinan perpindahan penduduk antarwilayah yang memiliki jarak geografis dekat (Rahman et al., 2019).

Selain dependensi spasial, terdapat heterogenitas spasial pada kasus TB (Lin et al., 2023). Heterogenitas spasial menggambarkan perbedaan hubungan antara variabel penjelas dan respons di setiap wilayah (Wubuli et al., 2015). Ketika varians residual tidak seragam antar lokasi atau terjadi heteroskedastisitas, hal ini mencerminkan adanya heterogenitas spasial yang tidak dapat dijelaskan oleh model regresi global (Brunsdon et al., 1998). Akibatnya, regresi klasik yang mengasumsikan hubungan homogen di seluruh wilayah hanya mampu menghasilkan estimasi parameter yang berlaku secara global dan tidak dapat merepresentasikan variasi spasial.

Karakteristik wilayah menjadi faktor penting dalam analisis insidensi TB meng-ingat setiap daerah memiliki kondisi yang berbeda, sehingga efek yang diberikan oleh faktor-faktor penyebab insidensi TB juga bervariasi antar wilayah (Shaweno et al., 2018). Salah satu karakteristik wilayah yang memengaruhi insidensi TB adalah kepadatan penduduk (Wang, Xin et al., 2025). Hubungan antara karakteristik wilayah dan variasi insidensi TB dapat diamati pada Gambar 1.2.

Gambar 1.2 Peta Insidensi Tuberkulosis Provinsi Jawa Barat Tahun 2020−2024
Gambar 1.2 Peta Insidensi Tuberkulosis Provinsi Jawa Barat Tahun 2020−2024

Berdasarkan Gambar 1.2, Kota Cimahi dengan kepadatan penduduk tinggi mem-iliki angka insidensi TB yang juga tinggi. Semakin padat penduduk suatu wilayah, maka penyebaran penyakit TB lebih mudah dan cepat akibat tingginya intensitas kontak antarindividu. Sementara itu, Kota Depok yang meskipun memiliki kepadatan penduduk tinggi menunjukkan tingkat insidensi TB yang relatif sedang. Kondisi ini dapat dikaitkan dengan akses layanan kesehatan yang lebih baik, tingginya cakupan deteksi dini, serta kedekatan dengan fasilitas kesehatan rujukan di wilayah DKI Jakarta sehingga penanganan kasus TB dapat dilakukan lebih cepat dan efektif.

Dalam analisis epidemiologi TB, data kejadian penyakit umumnya menunjukkan adanya dependensi spasial dan temporal. Wilayah yang berdekatan cenderung memiliki karakteristik insidensi yang mirip, sementara jumlah kasus pada suatu waktu juga di-pengaruhi oleh kondisi pada periode sebelumnya. Selain itu, data kasus TB sering kali menunjukkan variasi yang tinggi antarwilayah dan antarwaktu akibat perbedaan kon-disi sosial, ekonomi, kepadatan penduduk, akses layanan kesehatan, serta mobilitas masyarakat. Oleh karena itu, pendekatan statistik konvensional yang mengasumsikan independensi antarobservasi menjadi kurang sesuai untuk digunakan dalam pemodelan insidensi TB.

Untuk mengakomodasi ketergantungan spasial dan temporal secara simultan, digunakan pendekatan Bayesian Spatio-Temporal Model. Pendekatan Bayesian memungkinkan integrasi informasi spasial dan temporal melalui distribusi prior dan struktur dependensi antarwilayah maupun antarperiode waktu. Model ini mampu mengakomodasi heterogenitas spasial, autokorelasi spasial, serta dinamika temporal secara lebih fleksibel dibandingkan metode regresi klasik. Dalam konteks insidensi TB, Bayesian spatio-temporal model dapat menghasilkan estimasi risiko relatif yang lebih stabil, terutama pada wilayah dengan jumlah kasus kecil atau data yang berfluktuasi tinggi. Selain itu, pendekatan Bayesian juga memungkinkan pemetaan risiko penyakit yang lebih akurat melalui proses smoothing spasial-temporal sehingga pola persebaran TB dapat diidentifikasi secara lebih jelas.

Insidensi TB di Jawa Barat menunjukkan adanya variasi temporal yang cukup sig-nifikan, terutama setelah pandemi COVID-19. Pada tahun 2020–2021 terjadi penurunan pelaporan kasus akibat terganggunya sistem pelayanan kesehatan dan rendahnya akses masyarakat terhadap pemeriksaan TB. Banyak penderita TB mengalami keterlambatan diagnosis maupun pengobatan karena keterbatasan mobilitas dan kekhawatiran masyarakat terhadap risiko penularan COVID-19 di fasilitas kesehatan. Selain itu, seba-gian besar sumber daya kesehatan dialihkan untuk penanganan pandemi sehingga program pengendalian TB menjadi kurang optimal. Setelah pandemi mulai terkendali, terjadi peningkatan penemuan kasus TB pada tahun 2022 hingga 2024 sebagai dampak dari tertundanya proses deteksi dan pengobatan pada periode sebelumnya (Kemenkes RI, 2024a; PPID Dinkes Provinsi Jawa Barat, 2022).

Pendekatan Bayesian spatio-temporal model dinilai lebih sesuai untuk memodelkan fenomena tersebut karena mampu mempertimbangkan interaksi antara ruang dan waktu secara bersamaan. Risiko TB pada suatu kabupaten/kota tidak hanya dipengaruhi oleh kondisi wilayah tersebut pada periode sebelumnya, tetapi juga dapat dipengaruhi oleh kondisi wilayah lain yang berdekatan pada waktu tertentu. Dengan memasukkan efek spasial dan temporal secara simultan, model Bayesian spatio-temporal mampu menggambarkan pola penyebaran TB secara lebih realistis dan komprehensif. Selain itu, pendekatan ini dapat digunakan untuk mengidentifikasi wilayah dengan risiko tinggi secara lebih akurat sehingga dapat mendukung perencanaan intervensi kesehatan yang lebih efektif dan berbasis wilayah.

3 Sumber Data dan Variabel

Data yang digunakan adalah data insidensi TB di 27 kabupaten/kota di Provinsi Jawa Barat selama kurun waktu enam tahun, yaitu periode tahun 2020 sampai dengan tahun 2025, yang diperoleh dari Dinas Kesehatan Provinsi Jawa Barat dan Open Data Jawa Barat. Variabel yang digunakan dalam penelitian ini terdiri atas satu variabel respon dan enam variabel prediktor yang mewakili faktor sosial ekonomi, lingkungan, kependudukan, dan kesehatan. Deskripsi masing-masing variabel penelitian disajikan pada Tabel 1.

Tabel 1. Variabel Penelitian
Variabel Jenis Satuan Sumber
Insidensi TB Variabel Respon (Y) Kasus TB per 100.000 penduduk Dinas Kesehatan Jawa Barat
Tingkat Kemiskinan Variabel Prediktor (X1) Persentase Open Data Jabar
Air Minum Layak Variabel Prediktor (X2) Persentase Open Data Jabar
Sanitasi Layak Variabel Prediktor (X3) Persentase Open Data Jabar
PHBS Variabel Prediktor (X4) Persentase Open Data Jabar
Kepadatan Penduduk Variabel Prediktor (X5) Jiwa/km² Open Data Jabar
Insidensi HIV Variabel Prediktor (X6) Kasus HIV per 100.000 penduduk Open Data Jabar

#Hasil Analisis

## Warning: package 'dplyr' was built under R version 4.5.2
## 
## Attaching package: 'dplyr'
## The following object is masked from 'package:kableExtra':
## 
##     group_rows
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
## Warning: package 'ggplot2' was built under R version 4.5.2
## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
## Loading required package: sp
## 
## Attaching package: 'raster'
## The following object is masked from 'package:dplyr':
## 
##     select
## Loading required package: spData
## 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: 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
## Package 'mclust' version 6.1.2
## Type 'citation("mclust")' for citing this R package in publications.
## 
## Attaching package: 'mclust'
## The following object is masked from 'package:dplyr':
## 
##     count
## This is INLA_25.06.07 built 2025-06-11 19:15:03 UTC.
##  - See www.r-inla.org/contact-us for how to get help.
##  - List available models/likelihoods/etc with inla.list.models()
##  - Use inla.doc(<NAME>) to access documentation
## Loaded glmnet 4.1-10
## Registering fonts with R

3.1 Load Shapefile

jabar_sf <- st_read("Jabar.dbf")
## Reading layer `Jabar' from data source 
##   `/Users/kikiamelia/Desktop/Amelia/Jabar.dbf' using driver `ESRI Shapefile'
## Simple feature collection with 27 features and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 106.3705 ymin: -7.823398 xmax: 108.8338 ymax: -5.91377
## CRS:           NA
# Cek data shapefile
plot(st_geometry(jabar_sf))

# Plot peta dasar dengan label ID
B <- ggplot(data = jabar_sf) +
  geom_sf(fill = "white", color = "black") +
  geom_sf_text(
    aes(label = ID),
    size       = 3,
    color      = "black",
    fontface   = "bold",
    check_overlap = TRUE
  ) +
  theme_minimal()
B

3.2 Read Data

TB <- read.csv("DataTBJabar7.csv", header = TRUE, sep = ";", dec = ".")

jabar_tb_join <- jabar_sf %>%
  left_join(TB, by = c("Kabupaten" = "Kab.Kota"))

3.3 Visualisasi Peta Distribusi Insidensi Tuberkulosis

brks   <- c(39, 148, 206, 337, 551, 1400)
labels <- c()
format(brks, scientific = FALSE)
## [1] "  39" " 148" " 206" " 337" " 551" "1400"
labels <- labels[1:(length(labels) - 1)]

jabar_tb_join$brks <- cut(
  jabar_tb_join$insidensiTB,
  breaks        = brks,
  include.lowest = TRUE,
  labels        = labels
)

brks_scale  <- levels(jabar_tb_join$brks)
labels_scale <- rev(brks_scale)

A_tb<- ggplot(jabar_tb_join) +
  geom_sf(aes(fill = brks), color = "black", size = 0.1) +
  geom_sf_text(
    aes(label = id),
    size      = 1.5,
    color     = "black",
    fontface  = "bold",
    check_overlap = TRUE
  ) +
  facet_wrap(. ~ Tahun, ncol = 3) +
  scale_fill_manual(
    labels = c("<148", "[148 – 206)", "[206 – 337)", "[337 – 551)", ">=551"),
    values = c("grey95", "#FFFDEB", "#F97316", "#DC2626", "#7F1D1D")
  ) +
  theme_bw() +
  theme(
    axis.title.x  = element_blank(),
    axis.text.x   = element_blank(),
    axis.ticks.x  = element_blank(),
    axis.title.y  = element_blank(),
    axis.text.y   = element_blank(),
    axis.ticks.y  = element_blank(),
    legend.position = "right",
    text = element_text(size = 13)
  ) +
  labs(fill = "Jumlah Kasus TB")
A_tb

Peta insidensi Tuberkulosis (TB) menunjukkan adanya perubahan distribusi spasial kejadian TB di wilayah penelitian selama periode 2020–2025. Pada tahun 2020 hingga 2022, sebagian besar kabupaten/kota masih berada pada kategori insidensi rendah (<159,62%) dan sedang (159,62–229,89%), yang ditunjukkan oleh dominasi warna putih dan krem. Hanya beberapa wilayah yang memiliki insidensi lebih tinggi, terutama pada wilayah tertentu yang ditandai dengan warna oranye hingga merah.

Pada tahun 2023 terjadi peningkatan insidensi TB yang cukup nyata. Jumlah wilayah yang masuk dalam kategori tinggi (330,03–469,39%) dan sangat tinggi (≥469,39%) bertambah dibandingkan tahun-tahun sebelumnya. Peningkatan tersebut mulai terlihat pada beberapa wilayah bagian barat, tengah, dan utara.

Pada tahun 2024, wilayah dengan insidensi tinggi semakin meluas. Beberapa kabupaten/kota yang sebelumnya berada pada kategori sedang berubah menjadi kategori tinggi, menunjukkan peningkatan jumlah kasus TB dibandingkan periode sebelumnya. Selain itu, mulai terbentuk pengelompokan wilayah dengan insidensi tinggi yang berdekatan secara geografis.

Tahun 2025 menunjukkan kondisi dengan tingkat insidensi TB tertinggi selama periode pengamatan. Sebagian besar wilayah berada pada kategori tinggi hingga sangat tinggi, yang ditunjukkan oleh dominasi warna merah dan merah tua pada peta. Wilayah bagian barat, utara, dan beberapa wilayah tengah memperlihatkan konsentrasi insidensi yang relatif lebih tinggi dibandingkan wilayah lainnya.

Secara keseluruhan, peta insidensi TB menunjukkan tren peningkatan kejadian TB dari tahun 2020 hingga 2025. Selain peningkatan jumlah wilayah dengan insidensi tinggi, terlihat pula adanya pola pengelompokan spasial pada beberapa wilayah yang secara konsisten memiliki insidensi tinggi. Temuan ini mengindikasikan bahwa beban penyakit TB tidak tersebar secara merata, melainkan terkonsentrasi pada wilayah-wilayah tertentu yang memerlukan perhatian lebih dalam perencanaan program pengendalian dan pencegahan TB.

3.4 Matriks Bobot Spasial

# Beri nama baris shapefile sesuai urutan
row.names(jabar_sf) <- as.character(seq_len(nrow(jabar_sf)))

# Bangun matriks ketetanggaan (queen)
WB_nb  <- poly2nb(jabar_sf, queen = TRUE)
lw     <- nb2listw(WB_nb, style = "W", zero.policy = TRUE)
W_matrix <- listw2mat(lw)
W1     <- as(as_dgRMatrix_listw(lw), "CsparseMatrix")

# Plot peta + jaringan ketetanggaan
CoordK <- st_coordinates(st_centroid(jabar_sf))
## Warning: st_centroid assumes attributes are constant over geometries
plot(st_geometry(jabar_sf), axes = TRUE, col = "gray90")
text(CoordK[, 1], CoordK[, 2],
     labels = row.names(jabar_sf),
     col = "black", cex = 0.8, pos = 3)
points(CoordK[, 1], CoordK[, 2], pch = 19, cex = 0.7, col = "blue")
plot(WB_nb, CoordK, add = TRUE, col = "red", lwd = 2)

3.5 Persiapan Data Pemodelan

DataL <- read.csv("DataTBJabar7.csv", header = TRUE, sep = ";", dec = ".")
head(DataL)
##   id        Kab.Kota Tahun JumlahPenduduk JumlahKasusTB        Ei       SMR
## 1  1       Kab Bogor  2020        6067956         16332 13531.512 1.2069605
## 2  2    Kab Sukabumi  2020        2526878          5828  5634.925 1.0342639
## 3  3     Kab Cianjur  2020        2329635          5068  5195.074 0.9755394
## 4  4     Kab Bandung  2020        3773706          7595  8415.346 0.9025179
## 5  5       Kab Garut  2020        2644843          3516  5897.987 0.5961356
## 6  6 Kab Tasikmalaya  2020        1802165          2620  4018.819 0.6519328
##         lnSMR AirMinumLayak SanitasiLayak  PHBS KepadatanPenduduk HIV
## 1  0.18810521         92.13          62.4 52.10              2002 417
## 2  0.03368999         79.98          82.9 54.86               657 110
## 3 -0.02476469         82.73          83.2 64.72               645 189
## 4 -0.10256676         95.87          75.5 58.90              2050 176
## 5 -0.51728708         81.86          31.1 60.10               841   5
## 6 -0.42781375         85.17          91.4 57.20               731  95
##   TingkatKemiskinan insidensiTB InsidensiHIV IT IS         ST
## 1              7.69    269.1516    6.8721659  1  1 0.01378559
## 2              7.09    230.6403    4.3531979  1  2 0.02184470
## 3             10.36    217.5448    8.1128589  1  3 0.02203697
## 4              6.91    201.2610    4.6638503  1  4 0.01579764
## 5              9.98    132.9379    0.1890471  1  5 0.03025267
## 6             10.34    145.3807    5.2714374  1  6 0.01415867
# Pisahkan komponen data
DataLY <- data.frame(Case = DataL$insidensiTB)
DataLE <- data.frame(E    = DataL$Ei)

DataLX <- data.frame(
  TingkatKemiskinan = DataL$TingkatKemiskinan,
  AirMinumLayak     = DataL$AirMinumLayak,
  SanitasiLayak     = DataL$SanitasiLayak,
  PHBS              = DataL$PHBS,
  KepadatanPenduduk = DataL$KepadatanPenduduk,
  insidensiHIV      = DataL$InsidensiHIV
)


DataLlnSMR <- DataL$lnSMR
DataLXlnSMR <- data.frame(DataLlnSMR, DataLX)

# Gabungkan menjadi data frame penuh
FullDataL <- data.frame(Tahun = TB$Tahun, DataLY, DataLE, DataLX)

# Buat versi terstandarisasi (scaled) + indeks spasial & temporal
DATAL1  <- FullDataL
DATAL1S <- cbind(DATAL1[, 1:2], scale(DATAL1[, c(-1, -2)]))

DATAL1S$IS1 <- DataL$IS
DATAL1S$IS2 <- DataL$IS
DATAL1S$IS3 <- DataL$IS
DATAL1S$IT1 <- DataL$IT
DATAL1S$IT2 <- DataL$IT
DATAL1S$IT3 <- DataL$IT

3.6 Statistik Deskriptif

## Warning: package 'tidyr' was built under R version 4.5.2
## Warning: package 'purrr' was built under R version 4.5.2
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ lubridate 1.9.4     ✔ tibble    3.3.0
## ✔ purrr     1.2.2     ✔ tidyr     1.3.2
## ✔ readr     2.1.5     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ mclust::count()     masks dplyr::count()
## ✖ tidyr::expand()     masks Matrix::expand()
## ✖ tidyr::extract()    masks raster::extract()
## ✖ dplyr::filter()     masks stats::filter()
## ✖ dplyr::group_rows() masks kableExtra::group_rows()
## ✖ dplyr::lag()        masks stats::lag()
## ✖ purrr::lift()       masks caret::lift()
## ✖ purrr::map()        masks mclust::map()
## ✖ tidyr::pack()       masks Matrix::pack()
## ✖ raster::select()    masks dplyr::select()
## ✖ tidyr::unpack()     masks Matrix::unpack()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
Tabel 2. Statistik Deskriptif Variabel Penelitian
Variabel Tahun Mean SD Min Max
Insidensi TB (Per 100.000 Penduduk) 2020 248.49 129.31 82.16 551.33
2021 183.16 79.96 89.21 382.59
2022 204.25 117.23 96.86 594.27
2023 365.50 190.96 163.76 840.40
2024 420.91 224.80 172.42 1081.92
2025 448.21 217.38 208.68 1077.82
Tingkat Kemiskinan (%) 2020 8.42 2.81 2.45 12.97
2021 8.97 2.89 2.58 13.13
2022 8.65 2.82 2.53 12.77
2023 8.18 2.68 2.38 12.13
2024 8.01 2.65 2.34 11.93
2025 7.60 2.50 2.31 11.02
Air Minum Layak (%) 2020 93.72 6.35 79.98 100.00
2021 93.36 6.76 77.35 100.00
2022 93.71 5.84 78.15 99.44
2023 94.13 5.54 81.13 99.61
2024 95.10 4.22 85.37 99.85
2025 95.78 4.59 84.79 100.00
Sanitasi Layak (%) 2020 78.81 16.59 31.10 98.80
2021 86.92 13.02 43.90 100.00
2022 88.56 12.35 55.29 100.00
2023 86.07 17.72 29.66 100.00
2024 89.45 14.89 34.56 100.00
2025 76.79 15.98 46.61 97.93
PHBS (%) 2020 61.36 9.57 41.73 80.18
2021 62.75 11.45 30.94 80.77
2022 63.63 10.78 43.88 81.59
2023 64.26 10.73 41.00 83.24
2024 67.87 10.03 45.48 84.67
2025 67.87 10.03 45.48 84.67
Kepadatan Penduduk (jiwa/km²) 2020 3860.15 4551.68 419.00 14577.00
2021 3895.78 4584.33 423.00 14630.00
2022 3823.96 4547.67 383.00 14776.00
2023 3873.26 4624.98 382.00 15047.00
2024 3910.93 4668.12 385.00 15176.00
2025 3941.33 4701.59 387.00 15300.00
Insidensi HIV 2020 16.14 22.24 0.19 107.16
2021 13.74 17.07 0.99 79.59
2022 22.38 20.67 4.75 107.88
2023 23.73 18.39 4.67 91.23
2024 23.41 13.99 4.37 60.32
2025 21.09 13.35 6.42 62.09

Berdasarkan Tabel X, rata-rata jumlah kasus tuberkulosis di Jawa Barat menunjukkan kecenderungan meningkat dari 248,49 kasus pada tahun 2020 menjadi 448,21 kasus pada tahun 2025. Pada tahun 2020, jumlah kasus TB terendah sebesar 82,16 kasus terdapat di Kabupaten Pangandaran, sedangkan jumlah kasus tertinggi sebesar 551,33 kasus terdapat di Kota Sukabumi. Pada tahun 2021, jumlah kasus TB terendah sebesar 89,21 kasus ditemukan di Kabupaten Bandung, sementara jumlah kasus tertinggi sebesar 382,59 kasus terdapat di Kota Cirebon.

Pada tahun 2022, jumlah kasus TB terendah sebesar 96,86 kasus terdapat di Kabupaten Indramayu, sedangkan jumlah kasus tertinggi sebesar 594,27 kasus terdapat di Kota Cirebon. Selanjutnya, pada tahun 2023 jumlah kasus TB terendah sebesar 163,76 kasus ditemukan di Kabupaten Tasikmalaya, sementara jumlah kasus tertinggi sebesar 840,40 kasus terdapat di Kota Cirebon.

Pada tahun 2024, jumlah kasus TB terendah sebesar 172,42 kasus masih terdapat di Kabupaten Tasikmalaya, sedangkan jumlah kasus tertinggi mencapai 1.081,92 kasus di Kota Cirebon. Nilai tersebut merupakan jumlah kasus tertinggi selama periode penelitian. Pada tahun 2025, jumlah kasus TB terendah sebesar 208,68 kasus terdapat di Kabupaten Pangandaran, sedangkan jumlah kasus tertinggi sebesar 1.077,82 kasus terdapat di Kota Banjar.

Secara umum, Kota Cirebon merupakan wilayah yang paling konsisten memiliki jumlah kasus TB tertinggi selama periode 2021–2024, sedangkan Kota Sukabumi dan Kota Banjar masing-masing menjadi wilayah dengan jumlah kasus tertinggi pada tahun 2020 dan 2025. Di sisi lain, Kabupaten Pangandaran, Kabupaten Bandung, Kabupaten Indramayu, dan Kabupaten Tasikmalaya merupakan wilayah yang beberapa kali mencatat jumlah kasus TB terendah selama periode penelitian. Temuan ini menunjukkan adanya perbedaan beban kasus TB yang cukup besar antar kabupaten/kota di Jawa Barat, yang mengindikasikan pentingnya pendekatan spasial dalam upaya pengendalian penyakit TB.

3.7 ESTIMASI SMR (Standardized Morbidity Ratio)

Ei <- DataL$Ei

SMR_20 <- DataL$JumlahKasusTB[1:27]  / Ei[1:27]
SMR_21 <- DataL$JumlahKasusTB[28:54] / Ei[28:54]
SMR_22 <- DataL$JumlahKasusTB[55:81] / Ei[55:81]
SMR_23 <- DataL$JumlahKasusTB[82:108]/ Ei[82:108]
SMR_24 <- DataL$JumlahKasusTB[109:135]  / Ei[109:135] 
SMR_25 <- DataL$JumlahKasusTB[136:162] / Ei[136:162] 


SMR <- c(SMR_20, SMR_21, SMR_22, SMR_23, SMR_24, SMR_25)
summary(SMR)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.3684  0.7242  0.9298  1.1054  1.2095  3.2648
par(mar = c(3, 3, 2, 1))
boxplot(SMR_20, SMR_21, SMR_22, SMR_23, SMR_24, SMR_25,
        main  = "Boxplot SMR",
        names = c("2020","2021","2022", "2023", "2024", "2025"))

3.8 Pemodelan Bayesian Spatio Tempral

3.8.1 DEFINISI PRIOR UNTUK INLA

# (1) Inverse Gamma Prior
prior.c5  <- c(1, 1e-5)
igprior5  <- list(theta = list(param = prior.c5))

prior.c1  <- c(1, 0.01)
igprior1  <- list(theta = list(param = prior.c1))

# (2) Penalized Complexity Prior
pcprior <- list(prec = list(prior = "pc.prec", param = c(3, 0.01)))

# (3) Half-Cauchy Prior
halfcauchy <- "expression:
  lambda    = 25;
  precision = exp(log_precision);
  logdens   = -1.5*log_precision - log(pi*lambda) - log(1 + 1/(precision*lambda^2));
  log_jacobian = log_precision;
  return(logdens + log_jacobian);"

hcprior <- list(prec = list(prior = halfcauchy))

# (4) Uniform Prior
sdunif <- "expression:
  logdens = -log_precision/2;
  return(logdens)"

uprior <- list(prec = list(prior = sdunif))

# Control list INLA
control <- list(
  predictor = list(compute = TRUE, cdf = c(log(1))),
  results   = list(return.marginals.random = TRUE, return.marginals.predictor = TRUE),
  compute   = list(
    hyperpar = TRUE, return.marginals = TRUE, dic = TRUE,
    mlik = TRUE, cpo = TRUE, po = TRUE, waic = TRUE,
    graph = TRUE, gdensity = TRUE, openmp.strategy = "huge"
  )
)

3.8.2 EKSPLORASI DATA (KORELASI & VIF)

DataLX1       <- DataLX[-c(1:60), ]
VIF           <- diag(solve(cor(DataLX1)))

DataLXlnSMR1  <- DataLXlnSMR[-c(1:60), ]
Spearman.Correlation <- round(cor(DataLXlnSMR1, method = "spearman"), 1)

r <- ggcorrplot(Spearman.Correlation, type = "lower", lab = TRUE)
## Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
## ℹ Please use tidy evaluation idioms with `aes()`.
## ℹ See also `vignette("ggplot2-in-packages")` for more information.
## ℹ The deprecated feature was likely used in the ggcorrplot package.
##   Please report the issue at <https://github.com/kassambara/ggcorrplot/issues>.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
r

3.8.3 DATA TERSTANDARISASI UNTUK INLA

FullDataLS <- data.frame(
  FullDataL[, 1:2],
  scale(FullDataL[, -c(1, 2)])
)
head(FullDataLS)
##   Tahun     Case           E TingkatKemiskinan AirMinumLayak SanitasiLayak
## 1  2020 269.1516  2.00198309        -0.2262510    -0.3872416   -1.39940415
## 2  2020 230.6403  0.08470931        -0.4466739    -2.5556225   -0.09740126
## 3  2020 217.5448 -0.02208552         0.7546310    -2.0648367   -0.07834756
## 4  2020 201.2610  0.75978922        -0.5128008     0.2802271   -0.56739255
## 5  2020 132.9379  0.14858003         0.6150298    -2.2201035   -3.38734026
## 6  2020 145.3807 -0.30767776         0.7472836    -1.6293758    0.44245359
##           PHBS KepadatanPenduduk insidensiHIV
## 1 -1.183567914        -0.4144635   -0.7323041
## 2 -0.922716600        -0.7106292   -0.8719440
## 3  0.009165267        -0.7132716   -0.6635259
## 4 -0.540890764        -0.4038940   -0.8547228
## 5 -0.427477150        -0.6701129   -1.1027850
## 6 -0.701560052        -0.6943346   -0.8210411
m <- nrow(FullDataLS)
FullDataLS$beta1 <- rep(1, m)
FullDataLS$beta2 <- rep(2, m)
FullDataLS$beta3 <- rep(3, m)
FullDataLS$beta4 <- rep(4, m)
FullDataLS$beta5 <- rep(5, m)
FullDataLS$beta6 <- rep(6, m)
FullDataLS$beta7 <- rep(7, m)

# Tambahkan indeks spasial & temporal
FullDataLS$IS1 <- DataL$IS
FullDataLS$IS2 <- DataL$IS
FullDataLS$IS3 <- DataL$IS
FullDataLS$IT1 <- DataL$IT
FullDataLS$IT2 <- DataL$IT
FullDataLS$IT3 <- DataL$IT

param.beta <- list(prec = list(param = c(1e-6, 1e-6)))

3.8.4 PEMODELAN INLA

# Helper INLA control yang dipakai semua model
ctrl_fixed     <- list(mean = 0, mean.intercept = 0, prec = 0.0001, prec.intercept = 0.0001)
ctrl_predictor <- list(compute = TRUE, cdf = c(log(1)))
ctrl_compute   <- list(dic = TRUE, cpo = TRUE, waic = TRUE)
ctrl_inla      <- list(tolerance = 1e-20, h = 1e-8, int.strategy = "eb", strategy = "simplified.laplace")

# ---- GLM Poisson (seleksi variabel awal) ----
GLM0 <- glm(
  Case ~ TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
    PHBS + KepadatanPenduduk + insidensiHIV + offset(log(E)),
  family = "gaussian", data = FullDataL
)



# ---- ModelRun1: Fixed saja ----
ModelRun1 <- Case ~ TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
  PHBS + KepadatanPenduduk + insidensiHIV

RModelRun1 <- inla(ModelRun1, family = "Gaussian", data = FullDataLS, E = E,
                   control.fixed     = ctrl_fixed,
                   control.predictor = ctrl_predictor,
                   control.compute   = ctrl_compute,
                   control.inla      = ctrl_inla)
## 
##  *** inla.core.safe:  The inla program failed, but will rerun in case better initial values may help. try=1/1 
## 
##  *** inla.core.safe:  rerun with improved initial values
# ---- ModelRun2: Fixed + Spasial BYM ----
ModelRun2 <- Case ~ f(IS1, model = "bym", graph = W1, constr = TRUE) +
  TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
  PHBS + KepadatanPenduduk + insidensiHIV

RModelRun2 <- inla(ModelRun2, family = "Gaussian", data = FullDataLS, E = E,
                   control.fixed     = ctrl_fixed,
                   control.predictor = ctrl_predictor,
                   control.compute   = ctrl_compute,
                   control.inla      = ctrl_inla)
## 
##  *** inla.core.safe:  The inla program failed, but will rerun in case better initial values may help. try=1/1 
## 
##  *** inla.core.safe:  rerun with improved initial values
# ---- ModelRun3: Fixed + Temporal RW2 ----
ModelRun3 <- Case ~ f(IT1, model = "rw2", constr = TRUE) +
  TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
  PHBS + KepadatanPenduduk + insidensiHIV

RModelRun3 <- inla(ModelRun3, family = "Gaussian", data = FullDataLS, E = E,
                   control.fixed     = ctrl_fixed,
                   control.predictor = ctrl_predictor,
                   control.compute   = ctrl_compute,
                   control.inla      = ctrl_inla)



# ---- ModelRun4: Fixed + Temporal RW2 + Spasial BYM ----
ModelRun4 <- Case ~ TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
  PHBS + KepadatanPenduduk + insidensiHIV +
  f(IT1, model = "rw2", constr = TRUE) +
  f(IS1, model = "bym", graph = W1, constr = TRUE)

RModelRun4 <- inla(ModelRun4, family = "Gaussian", data = FullDataLS, E = E,
                   control.fixed     = ctrl_fixed,
                   control.predictor = ctrl_predictor,
                   control.compute   = ctrl_compute,
                   control.inla      = ctrl_inla)


# ---- ModelRun5: Fixed + Temporal RW2 + Spasial BYM + Interaksi ----
ModelRun5 <- Case ~ TingkatKemiskinan + AirMinumLayak + SanitasiLayak +
  PHBS + KepadatanPenduduk + insidensiHIV +
  f(IT1, model = "rw2", constr = TRUE) +
  f(IS1, model = "bym", graph = W1, constr = TRUE) +
  f(IS2, model = "bym", graph = W1, initial = 1, constr = TRUE,
    group = IS1, control.group = list(model = "rw2"))

RModelRun5 <- inla(ModelRun5, family = "Gaussian", data = FullDataLS, E = E,
                   control.fixed     = ctrl_fixed,
                   control.predictor = ctrl_predictor,
                   control.compute   = ctrl_compute,
                   control.inla      = ctrl_inla)
## 
##  *** inla.core.safe:  The inla program failed, but will rerun in case better initial values may help. try=1/1 
## 
##  *** inla.core.safe:  rerun with improved initial values

3.8.5 EVALUASI MODEL

pred1 <- RModelRun1$summary.fitted.values$mean[1:108]
pred2 <- RModelRun2$summary.fitted.values$mean[1:108]
pred3 <- RModelRun3$summary.fitted.values$mean[1:108]
pred4 <- RModelRun4$summary.fitted.values$mean[1:108]
pred5 <- RModelRun5$summary.fitted.values$mean[1:108]

Case  <- FullDataL$Case

R2 <- function(actual, predicted) {
  1 - sum((actual - predicted)^2) / sum((actual - mean(actual))^2)
}

R2_1 <- R2(Case, pred1)
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
R2_2 <- R2(Case, pred2)
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
R2_3 <- R2(Case, pred3)
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
R2_4 <- R2(Case, pred4)
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
R2_5 <- R2(Case, pred5)
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
DIC  <- c(RModelRun1$dic$dic, RModelRun2$dic$dic,
          RModelRun3$dic$dic, RModelRun4$dic$dic, RModelRun5$dic$dic)
rsq  <- c(R2_1, R2_2, R2_3, R2_4,R2_5)
library(knitr)
library(kableExtra)

model_fit <- data.frame(
  No = c(1,2,3,4,5),
  Model = c(
    "ηit = β0 + Xitᵀβ",
    "ηit = β0 + Xitᵀβ + ui + si",
    "ηit = β0 + Xitᵀβ + ωt",
    "ηit = β0 + Xitᵀβ + ui + si + ωt",
    "ηit = β0 + Xitᵀβ + ui + si + ωt + δit"
  ),
  DIC = c(2057.78,2057.83,2006.95,2005.68,1992.55),
  WAIC = c(2062,2060,2009,2006,1931),
  R_Square = c(0.53967,0.53967,0.669,0.673,0.853)
)

kbl(
  model_fit,
  caption = "Tabel 5. Hasil Ukuran Kecocokan Model",
  align = c("c","l","c","c","c"),
  digits = c(0,0,2,0,5),
  booktabs = TRUE
) %>%
  kable_styling(
    bootstrap_options = c(
      "striped",
      "hover",
      "condensed"
    ),
    full_width = FALSE,
    position = "center"
  )
Tabel 5. Hasil Ukuran Kecocokan Model
No Model DIC WAIC R_Square
1 ηit = β0 + Xitᵀβ 2057.78 2062 0.53967
2 ηit = β0 + Xitᵀβ + ui + si 2057.83 2060 0.53967
3 ηit = β0 + Xitᵀβ + ωt 2006.95 2009 0.66900
4 ηit = β0 + Xitᵀβ + ui + si + ωt 2005.68 2006 0.67300
5 ηit = β0 + Xitᵀβ + ui + si + ωt + δit 1992.55 1931 0.85300

Berdasarkan Tabel 5, apabila efek spasial dan temporal tidak dimasukkan ke dalam model memperoleh nilai DIC terbesar dan R-Square terkecil yaitu pada Model 1, sedangkan jika memasukan efek spasial, temporal, maupun interaksinya memperoleh nilai DIC yang jauh lebih kecil dan nilai R-Square yang jauh lebih besar. Sehingga berdasarkan hasil analisis yang ditampilkan pada Tabel 5 model yang paling cocok yaitu model

\[ \eta_{it}=\beta_0+\mathbf{X}*{it}^{T}\boldsymbol{\beta}+u_i+s_i+\omega_t+\delta*{it}\]

dengan nilai DIC sebesar 1992,55 dan nilai R-Square sebesar 0,853 yang menunjukkan kemampuan model dalam menjelaskan Risiko Relatif penyakit TB di Jawa Barat.

3.8.6 RISIKO RELATIF SPATIO-TEMPORAL (ST)

ST_20 <- RModelRun5$summary.fitted.values$mean[1:27]  / Ei[1:27]
ST_21 <- RModelRun5$summary.fitted.values$mean[28:54] / Ei[28:54]
ST_22 <- RModelRun5$summary.fitted.values$mean[55:81] / Ei[55:81]
ST_23 <- RModelRun5$summary.fitted.values$mean[82:108]/ Ei[82:108]
ST_24 <- RModelRun5$summary.fitted.values$mean[109:135]/ Ei[109:135]
ST_25 <- RModelRun5$summary.fitted.values$mean[136:162]/ Ei[136:162]

graphics.off()
par(mar = c(3, 3, 2, 1))
boxplot(ST_20, ST_21, ST_22, ST_23, ST_24, ST_25,
        main  = "Boxplot ST",
        names = c("2020","2021","2022", "2023", "2024", "2025"))

ST  <- c(ST_20, ST_21,ST_22, ST_23, ST_24, ST_25)
ST1 <- data.frame(ST_22, ST_23, ST_24, ST_25)


boxplot(SMR, ST,
        main  = "Boxplot SMR dan Spatio-Temporal Model",
        names = c("SMR", "ST"))

Gambar 2 menunjukan perbandingan hasil taksiran Risiko Relatif menggunakan SMR dan Spatio-Temporal Model. Hasil menunjukkan bahwa boxplot SMR memiliki ukuran boxplot yang lebih lebar dibandingkan dengan boxplot Spa-tio-Temporal Model. Artinya hasil taksiran Risiko Relatif menggunakan SMR memiliki varians yang lebih besar dibandingkan dengan hasil taksiran Risiko Relatif menggunakan Spatio-Temporal Model. Hal tersebut diduga karena SMR dipengaruhi oleh ukuran populasi sehingga hasil taksiran kurang tepat. Sehingga yang lebih baik digunakan yaitu Spatio-Temporal Model. Selain untuk melihat perbandingan hasil taksiran Risiko Relatif pada kedua metode tersebut, dengan menggunakan pendekatan Bayesian dapat dipelajari struktur ketergantungan spasial dan temporal yang dapat di-jadikan rujukan untuk mengkaji faktor risiko yang tidak dimasukkan secara eksplisit dalam model.

3.8.7 LOAD DATA TB & VISUALISASI PETA RISIKO RELATIF

TB <- read.csv("DataTBJabar7.csv", header = TRUE, sep = ";", dec = ".")
head(TB)
##   id        Kab.Kota Tahun JumlahPenduduk JumlahKasusTB        Ei       SMR
## 1  1       Kab Bogor  2020        6067956         16332 13531.512 1.2069605
## 2  2    Kab Sukabumi  2020        2526878          5828  5634.925 1.0342639
## 3  3     Kab Cianjur  2020        2329635          5068  5195.074 0.9755394
## 4  4     Kab Bandung  2020        3773706          7595  8415.346 0.9025179
## 5  5       Kab Garut  2020        2644843          3516  5897.987 0.5961356
## 6  6 Kab Tasikmalaya  2020        1802165          2620  4018.819 0.6519328
##         lnSMR AirMinumLayak SanitasiLayak  PHBS KepadatanPenduduk HIV
## 1  0.18810521         92.13          62.4 52.10              2002 417
## 2  0.03368999         79.98          82.9 54.86               657 110
## 3 -0.02476469         82.73          83.2 64.72               645 189
## 4 -0.10256676         95.87          75.5 58.90              2050 176
## 5 -0.51728708         81.86          31.1 60.10               841   5
## 6 -0.42781375         85.17          91.4 57.20               731  95
##   TingkatKemiskinan insidensiTB InsidensiHIV IT IS         ST
## 1              7.69    269.1516    6.8721659  1  1 0.01378559
## 2              7.09    230.6403    4.3531979  1  2 0.02184470
## 3             10.36    217.5448    8.1128589  1  3 0.02203697
## 4              6.91    201.2610    4.6638503  1  4 0.01579764
## 5              9.98    132.9379    0.1890471  1  5 0.03025267
## 6             10.34    145.3807    5.2714374  1  6 0.01415867
ST<-TB$ST
jabar_tb2_join <- jabar_sf %>%
  left_join(TB, by = c("Kabupaten" = "Kab.Kota"))

# Interval risiko relatif
brks <- c(0.0115, 0.0319, 0.0459, 0.0637, 0.167, Inf)
labels <- c()
format(brks, scientific = FALSE)
## [1] "0.0115" "0.0319" "0.0459" "0.0637" "0.1670" "   Inf"
labels_tb <- labels[1:(length(labels) - 1)]

jabar_tb2_join$brks <- cut(
  jabar_tb2_join$ST,
  breaks        = brks,
  include.lowest = TRUE,
  labels        = labels
)

brks_scale_tb  <- levels(jabar_tb2_join$brks)
labels_scale_tb <- rev(brks_scale_tb)

A_tbST <- ggplot(jabar_tb2_join) +
  geom_sf(aes(fill = brks), color = "black", size = 0.1) +
  geom_sf_text(
    aes(label = id),
    size      = 1.5,
    color     = "black",
    fontface  = "bold",
    check_overlap = TRUE
  ) +
  facet_wrap(. ~ Tahun, ncol = 3) +
  scale_fill_manual(
    labels = c(
      "<0.0319",
      "[0.0319 – 0.0459)",
      "[0.0459 – 0.0637)",
      "[0.0637 – 0.167)",
      ">=1.67"
    ),
    values = c("grey95", "#FFFDEB", "#F97316", "#DC2626", "#7F1D1D")
  ) +
  theme_bw() +
  theme(
    axis.title.x  = element_blank(),
    axis.text.x   = element_blank(),
    axis.ticks.x  = element_blank(),
    axis.title.y  = element_blank(),
    axis.text.y   = element_blank(),
    axis.ticks.y  = element_blank(),
    legend.position = "right",
    text = element_text(size = 13)
  ) +
  labs(fill = "RR TB")

A_tbST

Berdasarkan peta risiko relatif TB tahun 2020–2025, terlihat adanya perbedaan tingkat risiko relatif antar kabupaten/kota di Jawa Barat. Warna yang semakin gelap menun-jukkan risiko relatif yang semakin tinggi, sedangkan warna yang lebih terang menun-jukkan risiko relatif yang lebih rendah.

Pada tahun 2020, wilayah dengan risiko relatif tinggi terkonsentrasi di Kota Cirebon (22), Kota Tasikmalaya (26), Kota Bekasi (23), dan beberapa wilayah di seki-tarnya. Sebagian besar kabupaten lainnya masih berada pada kategori risiko rendah hingga sedang. Memasuki tahun 2021 dan 2022, terjadi peningkatan jumlah wilayah dengan risiko relatif tinggi, terutama di kawasan timur dan utara Jawa Barat, seperti Kabupaten Karawang (15), Kabupaten Bekasi (16), Kota Bekasi (23), Kabupaten Cirebon (9), dan Kota Cirebon (22). Tahun 2022 menunjukkan penyebaran risiko tertinggi, ditandai dengan semakin banyak wilayah yang masuk kategori risiko tinggi dan sangat tinggi.

Pada tahun 2023 hingga 2025, pola risiko relatif cenderung mengalami penurunan dibandingkan tahun 2022. Namun demikian, wilayah timur Jawa Barat masih konsisten menunjukkan risiko relatif yang lebih tinggi dibandingkan wilayah lainnya. Kota Cirebon (22) dan Kota Banjar (27) tetap berada pada kategori risiko tinggi pada sebagian besar periode pengamatan. Sementara itu, Kabupaten Bogor (1), Kabu-paten Sukabumi (2), Kabupaten Cianjur (3), Kabupaten Garut (5), dan Kabupaten Tasikmalaya (6) cenderung berada pada kategori risiko rendah selama periode penelitian.

Hasil pemetaan ini sejalan dengan hasil estimasi efek spasial-temporal sebe-lumnya, di mana Kabupaten Bogor merupakan wilayah yang paling sering memiliki nilai efek spasial-temporal terendah, sedangkan Kota Cirebon secara konsisten mem-iliki nilai efek spasial-temporal tertinggi pada periode 2020–2024. Kondisi tersebut menunjukkan bahwa wilayah perkotaan dengan kepadatan penduduk yang tinggi dan mobilitas masyarakat yang besar cenderung memiliki risiko relatif TB yang lebih tinggi dibandingkan wilayah kabupaten yang didominasi kawasan pedesaan. Secara umum, peta risiko relatif menunjukkan adanya pengelompokan spasial (spatial clustering), di mana wilayah dengan risiko tinggi cenderung berdekatan secara geografis, khususnya di bagian timur dan utara Jawa Barat. Temuan ini mengindikasi-kan bahwa faktor spasial dan temporal berperan penting dalam persebaran kasus TB sehingga pendekatan pemodelan spasio-temporal tepat digunakan untuk menggam-barkan variasi risiko relatif antarwilayah dan antarwaktu.

4 DAFTAR PUSTAKA

Badan Perencanaan Pembangunan Nasional. (2025). Rencana Pembangunan Jangka Menengah Nasional (RPJMN) 2025–2029. Jakarta: Bappenas.

Brunsdon, C., Fotheringham, A. S., & Charlton, M. E. (1998). Geographically weighted regression: Modelling spatial non-stationarity. Journal of the Royal Statistical Society: Series D (The Statistician), 47(3), 431–443.

Chang, A., & Renaud, F. (2022). Environmental determinants and tuberculosis risk: A spatial epidemiological perspective. International Journal of Environmental Research and Public Health.

Draper, N. R., & Smith, H. (1998). Applied Regression Analysis (3rd ed.). New York: John Wiley & Sons.

Hasanuddin, H., dkk. (2023). Hubungan sanitasi lingkungan dengan kejadian tuberkulosis di Indonesia. Jurnal Kesehatan Lingkungan Indonesia.

Islam, M., et al. (2025). Determinants of tuberculosis incidence using spatio-temporal modelling approaches. BMC Public Health.

Kementerian Kesehatan Republik Indonesia. (2023). Peta Jalan Eliminasi Tuberkulosis Indonesia 2030. Jakarta: Kementerian Kesehatan RI.

Kementerian Kesehatan Republik Indonesia. (2024a). Laporan Program Penanggulangan Tuberkulosis Indonesia Tahun 2023. Jakarta: Kementerian Kesehatan RI.

Kustanto, A. (2020). Kemiskinan dan kejadian tuberkulosis di Indonesia. Jurnal Kependudukan Indonesia.

Lin, Y., et al. (2023). Spatial heterogeneity and tuberculosis risk: A Bayesian modelling approach. Scientific Reports.

Ng, I. C., et al. (2012). Spatial dependence of tuberculosis incidence and its epidemiological implications. BMC Infectious Diseases.

PPID Dinas Kesehatan Provinsi Jawa Barat. (2022). Profil Kesehatan Provinsi Jawa Barat Tahun 2021. Bandung: Dinas Kesehatan Provinsi Jawa Barat.

Rahman, A., et al. (2019). Spatial clustering of tuberculosis and population mobility in urban areas. BMC Public Health.

Shaweno, D., et al. (2018). Bayesian spatio-temporal modelling of tuberculosis risk. Spatial and Spatio-temporal Epidemiology.

Sukmawati, S., & Sri Hazanah, N. (2025). Perilaku hidup bersih dan sehat serta kejadian tuberkulosis. Jurnal Kesehatan Masyarakat Indonesia.

Wang, X., et al. (2025). Population density and tuberculosis incidence: Evidence from a spatio-temporal analysis. Spatial Statistics.

Widyaningsih, N., & Budiawan, A. (2023). Faktor lingkungan dan kejadian tuberkulosis di Indonesia. Jurnal Epidemiologi Kesehatan Indonesia.

World Health Organization. (2024). Global Tuberculosis Report 2024. Geneva: World Health Organization.

Wubuli, A., et al. (2015). Bayesian spatio-temporal analysis of tuberculosis incidence. International Journal of Environmental Research and Public Health, 12(2), 2055–2070.