Link Dashboard: https://firgiealdiansyahf.shinyapps.io/Dashboard-HIV-Bali/

Pendahuluan

Latar Belakang

HIV/AIDS tetap menjadi masalah kesehatan masyarakat yang penting di Indonesia meskipun prevalensinya lebih rendah dibandingkan banyak negara lain. Diperkirakan ratusan ribu orang hidup dengan HIV di Indonesia, dan pencapaian target global 95-95-95 masih menghadapi tantangan di tingkat nasional. [1]. Di tingkat provinsi, Bali tercatat sebagai salah satu provinsi dengan jumlah kasus HIV yang relatif tinggi dibandingkan provinsi lain di Indonesia dalam beberapa survei dan laporan provinsi, yang menjadikan Bali area penting untuk kajian epidemiologi dan intervensi. [2][3]. Penyebaran kasus HIV tidak selalu merata secara geografis; studi spasial menunjukkan pola klaster (hotspot) yang terkait dengan faktor sosial-demografis dan lingkungan. Analisis spasial (mis. pemetaan klaster, autokorelasi spasial, dan regresi spasial) telah digunakan di berbagai studi untuk mengidentifikasi area berisiko tinggi dan variabel penjelasnya. [4][5]. Faktor kemiskinan berperan ganda: selain meningkatkan kerentanan ekonomis yang mendorong perilaku berisiko (mis. migrasi, pekerja seks komersial, keterbatasan akses layanan kesehatan), kemiskinan juga membatasi akses ke pencegahan dan pengobatan sehingga mempengaruhi beban penyakit di suatu wilayah. Penelitian di konteks Indonesia menunjukkan hubungan antara determinan sosial ekonomi dan pengetahuan/akses terkait HIV. [6][8]. Kepadatan penduduk juga sering dikaitkan dengan dinamika penularan penyakit menular seksual: area dengan kepadatan tinggi cenderung memiliki lebih banyak interaksi sosial dan jaringan seksual yang kompleks, yang dapat memfasilitasi transmisi jika tidak diimbangi intervensi pencegahan. Studi-studi dengan pendekatan spasial di tingkat provinsi/kabupaten telah menemukan asosiasi antara kepadatan dan kasus HIV/penyakit terkait di beberapa konteks Indonesia. [7][9]. Menggabungkan analisis spasial dengan variabel kemiskinan dan kepadatan penduduk di tingkat kabupaten/kota di Bali akan membantu: (1) memetakan distribusi kasus dan klaster, (2) menguji hubungan statistik antara kasus HIV dengan kemiskinan dan kepadatan, serta (3) mengidentifikasi prioritas area untuk intervensi program kesehatan publik. Metode yang sesuai meliputi pemetaan tematik, indeks autokorelasi spasial (mis. Moran’s I), dan model regresi spasial (mis. GWR, spatial lag/error). Pendekatan ini telah dipakai sukses pada studi-studi provinsi lain di Indonesia untuk menginformasikan kebijakan lokal. [5][8][10].

Rumusan Masalah

Penyebaran kasus HIV/AIDS di Provinsi Bali menunjukkan variasi antar wilayah kabupaten/kota yang diduga tidak terjadi secara acak. Perbedaan kondisi sosial dan kependudukan antar wilayah berpotensi memengaruhi besarnya kasus HIV/AIDS yang terjadi. Namun, hingga saat ini belum diketahui secara jelas bagaimana pola distribusi spasial kasus HIV/AIDS di Provinsi Bali serta apakah terdapat kecenderungan pengelompokan kasus di wilayah tertentu.

Selain itu, faktor kemiskinan dan kepadatan penduduk sering disebut sebagai determinan sosial yang berperan dalam dinamika penularan HIV/AIDS. Akan tetapi, hubungan antara kedua faktor tersebut dengan kasus HIV/AIDS di Provinsi Bali belum banyak dikaji dengan pendekatan spasial yang mampu menangkap keterkaitan antar wilayah. Perbedaan karakteristik kabupaten/kota juga memungkinkan adanya variasi pengaruh kemiskinan dan kepadatan penduduk terhadap kasus HIV/AIDS di masing-masing wilayah.

Oleh karena itu, permasalahan utama dalam penelitian ini adalah bagaimana pola distribusi spasial kasus HIV/AIDS di Provinsi Bali serta bagaimana pengaruh faktor kemiskinan dan kepadatan penduduk terhadap kasus HIV/AIDS antar kabupaten/kota di Provinsi Bali jika ditinjau secara spasial.

Tujuan

Penelitian ini bertujuan untuk menganalisis pola distribusi spasial kasus HIV/AIDS di Provinsi Bali pada tingkat kabupaten/kota serta mengidentifikasi adanya kecenderungan pengelompokan kasus antar wilayah. Selain itu, penelitian ini bertujuan untuk mengkaji pengaruh faktor kemiskinan dan kepadatan penduduk terhadap kasus HIV/AIDS dengan pendekatan spasial, sehingga dapat diketahui perbedaan pengaruh kedua faktor tersebut antar wilayah di Provinsi Bali sebagai dasar penentuan wilayah prioritas penanggulangan HIV/AIDS. Secara khusus, tujuan penelitian ini adalah: 1. Mengetahui pola distribusi spasial kasus HIV/AIDS di Provinsi Bali pada tingkat kabupaten/kota. 2. Mengidentifikasi adanya autokorelasi spasial atau klaster kasus HIV/AIDS antar kabupaten/kota di Provinsi Bali. 3. Menganalisis pengaruh tingkat kemiskinan terhadap kasus HIV/AIDS di Provinsi Bali secara spasial. 4. Menganalisis pengaruh kepadatan penduduk terhadap kasus HIV/AIDS di Provinsi Bali secara spasial. 5. Mengkaji perbedaan pengaruh faktor kemiskinan dan kepadatan penduduk terhadap kasus HIV/AIDS antar wilayah kabupaten/kota di Provinsi Bali.

Tinjauan Pustaka

Agent-Host-Environment

Analisis epidemiologi HIV/AIDS dapat dijelaskan melalui model agent–host–environment, yang merupakan kerangka dasar dalam memahami proses terjadinya penyakit menular.

  • Agent pada kasus HIV/AIDS adalah virus HIV yang ditularkan melalui kontak darah, hubungan seksual tidak aman, penggunaan jarum suntik tidak steril, serta penularan dari ibu ke anak.

  • Host meliputi individu atau kelompok masyarakat dengan karakteristik tertentu seperti usia produktif, perilaku seksual berisiko, tingkat pendidikan, serta kondisi sosial ekonomi.

  • Environment mencakup kondisi lingkungan sosial dan struktural, seperti tingkat kemiskinan, kepadatan penduduk, mobilitas penduduk, akses layanan kesehatan, serta norma sosial.

Dalam konteks wilayah dengan kepadatan penduduk tinggi dan tingkat kemiskinan tertentu, interaksi antara agent, host, dan environment dapat mempercepat penyebaran HIV/AIDS.

Ukuran Asosiasi

Ukuran epidemiologi digunakan untuk menggambarkan besarnya masalah kesehatan serta hubungan antara faktor risiko dan kejadian penyakit. Ukuran frekuensi yang umum digunakan dalam epidemiologi meliputi prevalensi dan insidensi. Insidensi menggambarkan jumlah kasus baru yang terjadi dalam suatu periode tertentu pada populasi berisiko, sedangkan prevalensi menunjukkan proporsi kasus yang ada pada suatu waktu tertentu.

Selain ukuran frekuensi, epidemiologi juga menggunakan ukuran asosiasi untuk menilai hubungan antara faktor risiko dan kejadian penyakit. Ukuran asosiasi seperti risk ratio, odds ratio, atau incidence rate ratio (IRR) digunakan untuk mengetahui besarnya pengaruh suatu faktor terhadap kejadian penyakit. Dalam analisis berbasis data agregat wilayah, ukuran asosiasi sering diperoleh melalui pendekatan pemodelan regresi yang disesuaikan dengan karakteristik data.

Dalam konteks penelitian HIV/AIDS, ukuran frekuensi dan asosiasi digunakan untuk membandingkan tingkat kejadian antar wilayah serta menilai hubungan antara faktor sosial-demografis, seperti kemiskinan dan kepadatan penduduk, dengan jumlah atau tingkat kejadian kasus HIV/AIDS.

Desain Studi

Desain studi yang digunakan dalam penelitian ini adalah cross-sectional. Desain studi epidemiologi merupakan kerangka penelitian yang digunakan untuk menjawab pertanyaan penelitian terkait distribusi dan determinan penyakit. Beberapa desain studi yang umum digunakan dalam epidemiologi antara lain studi potong lintang (cross-sectional), case-control, cohort, dan eksperimen. Pemilihan desain studi disesuaikan dengan tujuan penelitian, jenis data, serta sumber daya yang tersedia.

Penelitian ini menggunakan desain studi potong lintang (cross-sectional) dengan unit analisis kabupaten/kota di Provinsi Bali. Desain ini digunakan untuk menggambarkan distribusi kasus HIV/AIDS serta hubungannya dengan faktor kemiskinan dan kepadatan penduduk pada satu periode waktu tertentu. Studi potong lintang sesuai digunakan pada penelitian dengan data sekunder agregat wilayah dan memungkinkan analisis deskriptif maupun analitik secara spasial.

Metodologi

Sumber Data

Penelitian ini menggunakan data sekunder yang bersumber dari Badan Pusat Statistik Bali. Data yang digunakan meliputi jumlah kasus HIV/AIDS, jumlah penduduk, persentase penduduk miskin, dan kepadatan penduduk pada tingkat kabupaten/kota di Provinsi Bali. Data spasial wilayah administrasi kabupaten/kota diperoleh dalam bentuk shapefile yang digunakan untuk analisis dan visualisasi spasial.

Variabel

Variabel independen (Y) dan variabel dependen (X) yang digunakan adalah sebagai berikut:

Variabel Keterangan
Y Jumlah Kasus HIV/AIDS
X1 Persentase Penduduk Miskin
X2 Kepadatan Penduduk

Metode Analisis

Penelitian ini menggunakan desain studi Cross-Sectional dengan pendekatan Regresi Binomial Negatif. Penelitian ini bertujuan untuk melihat pengaruh variabel prediktor yang memengaruhi variabel respon pada tingkat kabupaten/kota di Provinsi Bali. Tahapan analisis yang dilakukan adalah sebagai berikut:

  1. Eksplorasi Data Dilakukan eksplorasi data melalui analisis deskriptif dan pembuatan peta sebaran kasus HIV/AIDS untuk melihat pola distribusi antar wilayah.

  2. Uji Autokorelasi Spasial Mengukur adanya ketergantungan spasial antar wilayah menggunakan statistik Moran’s I serta menggunakan Local Moran’s I (LISA) untuk mendeteksi klaster wilayah dengan tingkat kasus HIV/AIDS yang relatif tinggi atau rendah. Statistik Moran’s I dihitung dengan rumus:

\[ I = \frac{n}{\sum_{i}\sum_{j} w_{ij}} \cdot \frac{\sum_{i}\sum_{j} w_{ij}(x_i - \bar{x})(x_j - \bar{x})} {\sum_{i}(x_i - \bar{x})^2} \]

  1. Menghitung Insidensi Insidensi HIV/AIDS dihitung untuk mengetahui tingkat kejadian kasus baru pada periode tertentu. Insidensi per 100.000 penduduk dihitung menggunakan rumus:

\[ \text{Insidensi} = \frac{\text{Jumlah Kasus HIV/AIDS}}{\text{Jumlah Penduduk}} \times 100{,}000 \]

  1. Pemodelan Regresi Binomial Negatif Regresi Binomial Negatif digunakan sebagai model regresi terbaik untuk data jumlah kasus HIV/AIDS, terutama ketika terjadi overdispersi pada model Poisson. Model yang digunakan adalah:

\[ \log(\mu_i) = \beta_0 + \beta_1 X_{1i} + \beta_2 X_{2i} + \log(\text{Populasi}_i) \]

  1. Menghitung Incidence Rate Ratio (IRR) Incidence Rate Ratio digunakan sebagai ukuran asosiasi untuk melihat pengaruh variabel prediktor terhadap laju kejadian HIV/AIDS. IRR diperoleh dari:

\[ \text{IRR} = e^{\beta} \]

Alur Kerja

  1. Pengumpulan Data Data diperoleh dari Badan dengan variabel yang digunakan yaitu jumlah kasus HIV/AIDS, jumlah penduduk, persentase penduduk miskin, dan kepadatan penduduk pada tingkat kabupaten/kota di Provinsi Bali.

  2. Eksplorasi Data Dilakukan analisis deskriptif serta visualisasi spasial untuk melihat sebaran kasus HIV/AIDS antar wilayah.

  3. Uji Autokorelasi Spasial Menguji adanya pola ketergantungan spasial secara global menggunakan Moran’s I dan secara lokal menggunakan LISA.

  4. Menghitung Insidensi Menghitung insidensi HIV/AIDS untuk melihat tingkat kejadian kasus baru pada masing-masing wilayah.

  5. Memilih Model Terbaik Pemilihan model dilakukan dengan melihat tingkat dispersi dan nilai Akaike Information Criterion (AIC) untuk menentukan model regresi yang paling sesuai.

  6. Menghitung Incidence Rate Ratio (IRR) Menghitung IRR untuk menilai hubungan antara faktor paparan (kemiskinan dan kepadatan penduduk) dengan kejadian HIV/AIDS.

Hasil dan Pembahasan

library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.3.3
## Warning: package 'ggplot2' was built under R version 4.3.3
## Warning: package 'tidyr' was built under R version 4.3.3
## Warning: package 'readr' was built under R version 4.3.3
## Warning: package 'purrr' was built under R version 4.3.3
## Warning: package 'dplyr' was built under R version 4.3.3
## Warning: package 'stringr' was built under R version 4.3.3
## Warning: package 'forcats' was built under R version 4.3.3
## Warning: package 'lubridate' was built under R version 4.3.3
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.0.4     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(MASS)
## 
## Attaching package: 'MASS'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
library(broom)
## Warning: package 'broom' was built under R version 4.3.3
library(dplyr)
library(sf)
## Warning: package 'sf' was built under R version 4.3.3
## Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUE
library(spdep)
## Warning: package 'spdep' was built under R version 4.3.3
## Loading required package: spData
## Warning: package 'spData' was built under R version 4.3.3
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
library(tmap)
library(geodata)
## Loading required package: terra
## Warning: package 'terra' was built under R version 4.3.3
## terra 1.8.29
## 
## Attaching package: 'terra'
## 
## The following object is masked from 'package:MASS':
## 
##     area
## 
## The following object is masked from 'package:tidyr':
## 
##     extract
library(stringr)
# Set working directory (sesuaikan jika pindah folder)

setwd("C:/Users/Fizky Firdzansyah/Downloads/Semester 5/Epidemiologi")

# Baca data

data <- read.csv(
"UAS.csv",
sep = ";",
dec = ".",",",
header = TRUE
)

# Struktur awal

str(data)
## 'data.frame':    9 obs. of  7 variables:
##  $ KABKOT    : chr  "Jembrana" "Tabanan" "Badung" "Gianyar" ...
##  $ Laki.laki : num  163 233 285 262 105 ...
##  $ Perempuan : num  163 235 284 265 104 ...
##  $ Populasi  : num  326 468 569 527 209 ...
##  $ Kasus     : int  24 66 132 61 30 5 25 107 317
##  $ Persentase: num  4.51 4.4 2.23 4 5.3 5.06 6.52 5.39 2.59
##  $ Kepadatan : int  383 551 1426 1447 667 498 599 616 6003
names(data)
## [1] "KABKOT"     "Laki.laki"  "Perempuan"  "Populasi"   "Kasus"     
## [6] "Persentase" "Kepadatan"
data <- data %>%
rename(
kabkot      = KABKOT,
male        = `Laki.laki`,
female      = Perempuan,
population  = Populasi,
cases       = Kasus,
poverty     = Persentase,   # X1: persentase kemiskinan
density     = Kepadatan     # X2: kepadatan penduduk
) %>%
mutate(
kabkot = str_to_upper(kabkot)
)

data

Tabel Hasil

data_deskriptif <- data %>%
  dplyr::select(cases, population, density, poverty) %>%
  dplyr::summarise(
    across(everything(), 
      list(
        mean = ~mean(., na.rm = TRUE),
        median = ~median(., na.rm = TRUE),
        sd = ~sd(., na.rm = TRUE),
        min = ~min(., na.rm = TRUE),
        max = ~max(., na.rm = TRUE)
      ),
      .names = "{.col}_{.fn}"
    )
  )

print(data_deskriptif)
##   cases_mean cases_median cases_sd cases_min cases_max population_mean
## 1   85.22222           61 96.30651         5       317        492.5856
##   population_median population_sd population_min population_max density_mean
## 1            502.33      206.7021         209.31         814.79     1354.444
##   density_median density_sd density_min density_max poverty_mean poverty_median
## 1            616   1786.549         383        6003     4.444444           4.51
##   poverty_sd poverty_min poverty_max
## 1   1.362875        2.23        6.52

Ukuran Frekuensi Epidemiologi

# Insidensi 
data <- data %>%
mutate(
incidence_per100k = (cases / population) * 100000
)

summary(data$incidence_per100k)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1906    7370   13132   14731   14333   41955

Peta Persebaran Kasus HIV

indo_shp <- geodata::gadm(
country = "IDN",
level = 2,
path = "C:/Users/Fizky Firdzansyah/Downloads/Semester 5/Epidemiologi"
)

indo_map <- st_as_sf(indo_shp)
bali <- indo_map %>%
filter(NAME_1 == "Bali") %>%
mutate(kabkot = str_to_upper(NAME_2))

data_spatial <- bali %>%
left_join(data, by = "kabkot") %>%
filter(!is.na(cases))

glimpse(data_spatial)
## Rows: 7
## Columns: 22
## $ GID_2             <chr> "IDN.2.1_1", "IDN.2.2_1", "IDN.2.3_1", "IDN.2.5_1", …
## $ GID_0             <chr> "IDN", "IDN", "IDN", "IDN", "IDN", "IDN", "IDN"
## $ COUNTRY           <chr> "Indonesia", "Indonesia", "Indonesia", "Indonesia", …
## $ GID_1             <chr> "IDN.2_1", "IDN.2_1", "IDN.2_1", "IDN.2_1", "IDN.2_1…
## $ NAME_1            <chr> "Bali", "Bali", "Bali", "Bali", "Bali", "Bali", "Bal…
## $ NL_NAME_1         <chr> NA, NA, NA, NA, NA, NA, NA
## $ NAME_2            <chr> "Badung", "Bangli", "Buleleng", "Gianyar", "Jembrana…
## $ VARNAME_2         <chr> NA, NA, NA, NA, NA, NA, NA
## $ NL_NAME_2         <chr> NA, NA, NA, NA, NA, NA, NA
## $ TYPE_2            <chr> "Kabupaten", "Kabupaten", "Kabupaten", "Kabupaten", …
## $ ENGTYPE_2         <chr> "Regency", "Regency", "Regency", "Regency", "Regency…
## $ CC_2              <chr> "5103", "5106", "5108", "5104", "5101", "5105", "510…
## $ HASC_2            <chr> "ID.BA.BD", "ID.BA.BN", "ID.BA.BL", "ID.BA.GI", "ID.…
## $ kabkot            <chr> "BADUNG", "BANGLI", "BULELENG", "GIANYAR", "JEMBRANA…
## $ male              <dbl> 285.00, 131.95, 408.95, 262.38, 163.07, 104.99, 232.…
## $ female            <dbl> 283.54, 130.35, 405.83, 264.74, 162.56, 104.32, 234.…
## $ population        <dbl> 568.55, 262.30, 814.79, 527.12, 325.63, 209.31, 467.…
## $ cases             <int> 132, 5, 107, 61, 24, 30, 66
## $ poverty           <dbl> 2.23, 5.06, 5.39, 4.00, 4.51, 5.30, 4.40
## $ density           <int> 1426, 498, 616, 1447, 383, 667, 551
## $ incidence_per100k <dbl> 23216.955, 1906.214, 13132.218, 11572.317, 7370.328,…
## $ geometry          <GEOMETRY [°]> MULTIPOLYGON (((115.2039 -8..., POLYGON ((115.347 -8…

Peta Insidensi HIV

tmap_mode("plot")
## ℹ tmap modes "plot" - "view"
## ℹ toggle with `tmap::ttm()`
tm_shape(data_spatial) +
tm_polygons(
col = "incidence_per100k",
palette = "Reds",
title = "Insidensi HIV per 100.000 Penduduk"
) +
tm_borders() +
tm_layout(
title = "Sebaran Insidensi HIV di Provinsi Bali",
legend.outside = TRUE,
frame = FALSE
)
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_tm_polygons()`: migrate the argument(s) related to the scale of
## the visual variable `fill` namely 'palette' (rename to 'values') to fill.scale
## = tm_scale(<HERE>).[v3->v4] `tm_polygons()`: use 'fill' for the fill color of polygons/symbols
## (instead of 'col'), and 'col' for the outlines (instead of 'border.col').[v3->v4] `tm_polygons()`: migrate the argument(s) related to the legend of the
## visual variable `fill` namely 'title' to 'fill.legend = tm_legend(<HERE>)'[v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`[cols4all] color palettes: use palettes from the R package cols4all. Run
## `cols4all::c4a_gui()` to explore them. The old palette name "Reds" is named
## "brewer.reds"Multiple palettes called "reds" found: "brewer.reds", "matplotlib.reds". The first one, "brewer.reds", is returned.

Autokorelasi Spasial Global (Moran’s I)

map_sp <- as_Spatial(data_spatial)

nb <- poly2nb(map_sp)
lw <- nb2listw(nb, style = "W")

moran.test(
data_spatial$incidence_per100k,
lw
)
## 
##  Moran I test under randomisation
## 
## data:  data_spatial$incidence_per100k  
## weights: lw    
## 
## Moran I statistic standard deviate = -0.55979, p-value = 0.7122
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.27663328       -0.16666667        0.03858906

LISA & Klaster Spasial

y <- st_drop_geometry(data_spatial)$incidence_per100k
lag_y <- lag.listw(lw, y)

mean_y <- mean(y)

data_spatial$lisa_cluster <- case_when(
y > mean_y & lag_y > mean(lag_y) ~ "High-High",
y < mean_y & lag_y < mean(lag_y) ~ "Low-Low",
y > mean_y & lag_y < mean(lag_y) ~ "High-Low",
y < mean_y & lag_y > mean(lag_y) ~ "Low-High",
TRUE ~ "Undefined"
)

table(data_spatial$lisa_cluster)
## 
## High-High  High-Low  Low-High 
##         1         3         3

Peta Klaster LISA

prev <- st_drop_geometry(data_spatial)$prevalence_per100k
tm_shape(data_spatial) +
tm_polygons(
col = "lisa_cluster",
palette = c("red", "blue", "pink", "lightblue", "gray"),
title = "Klaster Spasial HIV (LISA)"
) +
tm_borders() +
tm_layout(
legend.outside = TRUE,
frame = FALSE
)
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_tm_polygons()`: migrate the argument(s) related to the scale of
## the visual variable `fill` namely 'palette' (rename to 'values') to fill.scale
## = tm_scale(<HERE>).
## [v3->v4] `tm_polygons()`: migrate the argument(s) related to the legend of the
## visual variable `fill` namely 'title' to 'fill.legend = tm_legend(<HERE>)'

## Regresi Poisson

model_pois <- glm(
cases ~ poverty + density + offset(log(population)),
family = poisson(link = "log"),
data = data
)

summary(model_pois)
## 
## Call:
## glm(formula = cases ~ poverty + density + offset(log(population)), 
##     family = poisson(link = "log"), data = data)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.173e+00  1.794e-01  -6.539 6.21e-11 ***
## poverty     -2.424e-01  3.748e-02  -6.469 9.89e-11 ***
## density      1.548e-04  1.933e-05   8.009 1.16e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 367.711  on 8  degrees of freedom
## Residual deviance:  55.404  on 6  degrees of freedom
## AIC: 113.03
## 
## Number of Fisher Scoring iterations: 4

Overdisperssion

dispersion <- model_pois$deviance / model_pois$df.residual
dispersion
## [1] 9.233999

Terdapat overdispersi maka menggunakan negative binomial

Regresi Negative Binomial

if (dispersion > 1.5) {
model_nb <- glm.nb(
cases ~ poverty + density + offset(log(population)),
data = data
)
summary(model_nb)
}
## 
## Call:
## glm.nb(formula = cases ~ poverty + density + offset(log(population)), 
##     data = data, init.theta = 6.103213305, link = log)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)  
## (Intercept) -1.1804810  0.7394303  -1.596   0.1104  
## poverty     -0.2540977  0.1418668  -1.791   0.0733 .
## density      0.0001631  0.0001054   1.548   0.1216  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Negative Binomial(6.1032) family taken to be 1)
## 
##     Null deviance: 28.279  on 8  degrees of freedom
## Residual deviance: 10.384  on 6  degrees of freedom
## AIC: 90.691
## 
## Number of Fisher Scoring iterations: 1
## 
## 
##               Theta:  6.10 
##           Std. Err.:  3.70 
## 
##  2 x log-likelihood:  -82.691

Uji Korelasi

cor.test(
  data$poverty,
  data$incidence_per100k,
  method = "spearman"
)
## 
##  Spearman's rank correlation rho
## 
## data:  data$poverty and data$incidence_per100k
## S = 190, p-value = 0.108
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.5833333
cor.test(
  data$density,
  data$incidence_per100k,
  method = "spearman"
)
## 
##  Spearman's rank correlation rho
## 
## data:  data$density and data$incidence_per100k
## S = 36, p-value = 0.04325
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho 
## 0.7

Incidence Rate Ratio (IRR)

IRR <- tidy(model_nb, exponentiate = TRUE, conf.int = TRUE)
IRR

Tabel prevalensi tinggi

data[
  order(-data$incidence_per100k),
  c("kabkot", "cases", "incidence_per100k", "poverty", "density")
]

Interpretasi

Secara keseluruhan, hasil analisis menunjukkan bahwa distribusi kasus HIV/AIDS di Provinsi Bali bervariasi antar kabupaten/kota. Analisis spasial global tidak menemukan adanya autokorelasi spasial yang signifikan, namun analisis lokal (LISA) mengindikasikan adanya klaster pada beberapa wilayah tertentu.

Pemodelan regresi menunjukkan bahwa kepadatan penduduk cenderung berasosiasi positif dengan jumlah kasus HIV/AIDS, sedangkan kemiskinan menunjukkan hubungan yang lebih lemah. Hal ini mengindikasikan bahwa dinamika sosial dan interaksi penduduk di wilayah dengan kepadatan tinggi berperan lebih besar dalam penyebaran HIV/AIDS dibandingkan faktor ekonomi semata.

Hasil ini menegaskan pentingnya pendekatan spasial dan berbasis wilayah dalam perencanaan program pencegahan dan pengendalian HIV/AIDS, khususnya dengan memprioritaskan wilayah dengan kepadatan penduduk tinggi dan indikasi klaster kasus.

Kesimpulan dan Saran

Kesimpulan

Berdasarkan hasil analisis spasial dan pemodelan statistik terhadap kasus HIV/AIDS di Provinsi Bali pada tingkat kabupaten/kota, dapat disimpulkan bahwa distribusi kasus HIV/AIDS menunjukkan variasi yang cukup besar antar wilayah. Wilayah dengan karakteristik perkotaan dan kepadatan penduduk tinggi cenderung memiliki jumlah kasus dan tingkat insidensi yang lebih besar dibandingkan wilayah lainnya.

Hasil uji autokorelasi spasial global menggunakan statistik Moran’s I menunjukkan bahwa secara keseluruhan tidak terdapat pengelompokan spasial yang signifikan. Namun demikian, analisis autokorelasi spasial lokal (LISA) mengindikasikan adanya klaster lokal pada beberapa wilayah tertentu, yang menunjukkan bahwa pola spasial HIV/AIDS tidak sepenuhnya acak di tingkat lokal.

Pemodelan regresi diawali dengan Regresi Poisson dan dilanjutkan dengan Regresi Binomial Negatif akibat ditemukannya overdispersi pada model Poisson. Model Binomial Negatif dipilih sebagai model terbaik. Hasil pemodelan menunjukkan bahwa kepadatan penduduk memiliki kecenderungan hubungan positif dengan jumlah kasus HIV/AIDS, sedangkan variabel kemiskinan menunjukkan hubungan yang lebih lemah dan tidak signifikan secara statistik.

Secara keseluruhan, hasil penelitian ini menunjukkan bahwa faktor spasial dan karakteristik kependudukan, khususnya kepadatan penduduk, berperan dalam variasi kasus HIV/AIDS antar wilayah di Provinsi Bali.

Saran

Berdasarkan hasil penelitian yang telah diperoleh, beberapa saran yang dapat diberikan adalah sebagai berikut:

  1. Pemerintah daerah dan instansi kesehatan disarankan untuk memprioritaskan upaya pencegahan dan pengendalian HIV/AIDS pada wilayah dengan kepadatan penduduk tinggi serta wilayah yang teridentifikasi memiliki klaster kasus berdasarkan analisis spasial lokal.

  2. Program intervensi HIV/AIDS perlu mempertimbangkan pendekatan berbasis wilayah dengan memperhatikan karakteristik sosial dan kependudukan masing-masing kabupaten/kota.

  3. Penelitian selanjutnya disarankan untuk menggunakan data dengan rentang waktu yang lebih panjang (spatio-temporal) serta menambahkan variabel lain seperti mobilitas penduduk, akses layanan kesehatan, dan faktor perilaku untuk memperoleh gambaran yang lebih komprehensif.

  4. Pendekatan analisis spasial yang lebih lanjut, seperti regresi spasial atau Geographically Weighted Regression (GWR), dapat dipertimbangkan untuk menangkap variasi pengaruh faktor risiko antar wilayah secara lebih detail.

Daftar Pustaka

[1] U. Evaluation Office, “Indonesia (country case study) — Evaluation of the contribution of the UNAIDS Joint Programme to strengthening HIV and Primary Health Care outcomes,” 2023. [Online]. Available: http://www.unaids.org/en/whoweare/evaluation

[2] L. N. Armini, E. P. Setiawati, N. Arisanti, and D. Hilmanto, “Evaluation of Process Indicators and Challenges of the Elimination of Mother-to-Child Transmission of HIV, Syphilis, and Hepatitis B in Bali Province, Indonesia (2019–2022): A Mixed Methods Study,” Tropical Medicine and Infectious Disease, vol. 8, no. 11, Nov. 2023, doi: 10.3390/tropicalmed8110492.

[3] B. Dhyah Kunthi Wardhani et al., “Very high HIV prevalence and incidence among men who have sex with men and transgender women in Indonesia: a retrospective observational cohort study in Bali and Jakarta, 2017-2020,” Journal of the International AIDS Society, vol. 27, p. 26386, 2024, doi: 10.1002/jia2.26386/full.

[4] A. Edmonds et al., “Associations between population density and clinical and sociodemographic factors in women living with HIV in the Southern United States,” AIDS Care - Psychological and Socio-Medical Aspects of AIDS/HIV, vol. 33, no. 2, pp. 229–238, 2021, doi: 10.1080/09540121.2020.1769829.

[5] H. P. Liew and T. Brooks, “Spatial clusters of AIDS in Indonesia,” Health Policy and Technology, vol. 6, no. 2, pp. 208–213, Jun. 2017, doi: 10.1016/j.hlpt.2017.01.005.

[6] T. R. P. Lestari, “Kebijakan Pengendalian HIV/AIDS di Denpasar,” Kesmas: National Public Health Journal, vol. 8, no. 1, p. 45, Aug. 2013, doi: 10.21109/kesmas.v8i1.341.

[7] S. Fauziana Thahar and T. Sirait, “Spatial Analysis of Variabels Affecting the Number of New HIV/AIDS Cases…..………….…(Safira Fauziana Thahar, Timbang Sirait) Analisis Spasial Variabel-Variabel yang Memengaruhi Jumlah Kasus Baru HIV/AIDS di Provinsi Jawa Timur Tahun 2021 (Spatial Analysis of Variables Affecting the Number of New HIV/AIDS Cases in East Java Province in 2021).”

[8] F. D. Virdausi et al., “Socio-Economic and Demographic Factors Associated with Knowledge and Attitude of HIV/AIDS among Women Aged 15–49 Years Old in Indonesia,” Healthcare (Switzerland), vol. 10, no. 8, Aug. 2022, doi: 10.3390/healthcare10081545.

[9] G. A. D. Maha Yoga and G. I. S. Diputra, “Spatial Analysis of Regional Poverty Rates in Bali for the period 2016-2020,” Jurnal Ekonomi Pembangunan, vol. 13, no. 1, pp. 57–68, Jun. 2024, doi: 10.23960/jep.v13i1.901.

[10] P. Spasial et al., “SPATIAL MODELING OF ENVIRONMENTAL-BASED RISK FACTORS OF TUBERCULOSIS IN BALI PROVINCE: AN ECOLOGICAL STUDY,” vol. 8, no. 1, pp. 26–34, 2020, doi: 10.20473/jbe.v8i12020.

[11] Badan Pusat Statistik Provinsi Bali. Kasus Penyakit Menurut Kabupaten/Kota dan Jenis Penyakit di Provinsi Bali, 2024. Diakses pada 17 Desember 2025, dari https://bali.bps.go.id/id/statistics-table/3/YTA1Q1ptRmhUMEpXWTBsQmQyZzBjVzgwUzB4aVp6MDkjMw==/kasus-penyakit-menurut-kabupaten-kota-dan-jenis-penyakit-di-provinsi-bali--2020.html?year=2024

[12] Badan Pusat Statistik Provinsi Bali. Penduduk, Laju Pertumbuhan Penduduk, Distribusi Persentase Penduduk Kepadatan Penduduk, Rasio Jenis Kelamin Penduduk Menurut Kabupaten/Kota di Provinsi Bali, 2024. Diakses pada 16 Desember 2025, dari https://bali.bps.go.id/id/statistics-table/3/V1ZSbFRUY3lTbFpEYTNsVWNGcDZjek53YkhsNFFUMDkjMw==/jumlah-penduduk--laju-pertumbuhan-penduduk--distribusi-persentase-penduduk--kepadatan-penduduk--rasio-jenis-kelamin-penduduk-menurut-kabupaten-kota-di-provinsi-bali--2021.html?year=2024

[13] Badan Pusat Statistik Provinsi Bali. (9 September 2025). Persentase Penduduk Miskin Provinsi Bali Menurut Kabupaten/Kota, 2024. Diakses pada 17 Desember 2025, dari https://bali.bps.go.id/id/statistics-table/2/MTI1IzI=/persentase-penduduk-miskin-provinsi-bali-menurut-kabupaten-kota.html

Lampiran

Import Library

library(tidyverse)
library(MASS)
library(broom)
library(dplyr)
library(sf)
library(spdep)
library(tmap)
library(geodata)
library(stringr)

Data

# Set working directory (sesuaikan jika pindah folder)

setwd("C:/Users/Fizky Firdzansyah/Downloads/Semester 5/Epidemiologi")

# Baca data

data <- read.csv(
"UAS.csv",
sep = ";",
dec = ".",",",
header = TRUE
)

# Struktur awal

str(data)
## 'data.frame':    9 obs. of  7 variables:
##  $ KABKOT    : chr  "Jembrana" "Tabanan" "Badung" "Gianyar" ...
##  $ Laki.laki : num  163 233 285 262 105 ...
##  $ Perempuan : num  163 235 284 265 104 ...
##  $ Populasi  : num  326 468 569 527 209 ...
##  $ Kasus     : int  24 66 132 61 30 5 25 107 317
##  $ Persentase: num  4.51 4.4 2.23 4 5.3 5.06 6.52 5.39 2.59
##  $ Kepadatan : int  383 551 1426 1447 667 498 599 616 6003
names(data)
## [1] "KABKOT"     "Laki.laki"  "Perempuan"  "Populasi"   "Kasus"     
## [6] "Persentase" "Kepadatan"
data <- data %>%
rename(
kabkot      = KABKOT,
male        = `Laki.laki`,
female      = Perempuan,
population  = Populasi,
cases       = Kasus,
poverty     = Persentase,   # X1: persentase kemiskinan
density     = Kepadatan     # X2: kepadatan penduduk
) %>%
mutate(
kabkot = str_to_upper(kabkot)
)

data

Ukuran Frekuensi Epidemiologi

# Insidensi 
data <- data %>%
mutate(
incidence_per100k = (cases / population) * 100000
)

summary(data$incidence_per100k)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1906    7370   13132   14731   14333   41955

Peta Persebaran Kasus HIV

indo_shp <- geodata::gadm(
country = "IDN",
level = 2,
path = "C:/Users/Fizky Firdzansyah/Downloads/Semester 5/Epidemiologi"
)

indo_map <- st_as_sf(indo_shp)
bali <- indo_map %>%
filter(NAME_1 == "Bali") %>%
mutate(kabkot = str_to_upper(NAME_2))

data_spatial <- bali %>%
left_join(data, by = "kabkot") %>%
filter(!is.na(cases))

glimpse(data_spatial)
## Rows: 7
## Columns: 22
## $ GID_2             <chr> "IDN.2.1_1", "IDN.2.2_1", "IDN.2.3_1", "IDN.2.5_1", …
## $ GID_0             <chr> "IDN", "IDN", "IDN", "IDN", "IDN", "IDN", "IDN"
## $ COUNTRY           <chr> "Indonesia", "Indonesia", "Indonesia", "Indonesia", …
## $ GID_1             <chr> "IDN.2_1", "IDN.2_1", "IDN.2_1", "IDN.2_1", "IDN.2_1…
## $ NAME_1            <chr> "Bali", "Bali", "Bali", "Bali", "Bali", "Bali", "Bal…
## $ NL_NAME_1         <chr> NA, NA, NA, NA, NA, NA, NA
## $ NAME_2            <chr> "Badung", "Bangli", "Buleleng", "Gianyar", "Jembrana…
## $ VARNAME_2         <chr> NA, NA, NA, NA, NA, NA, NA
## $ NL_NAME_2         <chr> NA, NA, NA, NA, NA, NA, NA
## $ TYPE_2            <chr> "Kabupaten", "Kabupaten", "Kabupaten", "Kabupaten", …
## $ ENGTYPE_2         <chr> "Regency", "Regency", "Regency", "Regency", "Regency…
## $ CC_2              <chr> "5103", "5106", "5108", "5104", "5101", "5105", "510…
## $ HASC_2            <chr> "ID.BA.BD", "ID.BA.BN", "ID.BA.BL", "ID.BA.GI", "ID.…
## $ kabkot            <chr> "BADUNG", "BANGLI", "BULELENG", "GIANYAR", "JEMBRANA…
## $ male              <dbl> 285.00, 131.95, 408.95, 262.38, 163.07, 104.99, 232.…
## $ female            <dbl> 283.54, 130.35, 405.83, 264.74, 162.56, 104.32, 234.…
## $ population        <dbl> 568.55, 262.30, 814.79, 527.12, 325.63, 209.31, 467.…
## $ cases             <int> 132, 5, 107, 61, 24, 30, 66
## $ poverty           <dbl> 2.23, 5.06, 5.39, 4.00, 4.51, 5.30, 4.40
## $ density           <int> 1426, 498, 616, 1447, 383, 667, 551
## $ incidence_per100k <dbl> 23216.955, 1906.214, 13132.218, 11572.317, 7370.328,…
## $ geometry          <GEOMETRY [°]> MULTIPOLYGON (((115.2039 -8..., POLYGON ((115.347 -8…

Peta Insidensi HIV

tmap_mode("plot")
## ℹ tmap modes "plot" - "view"
tm_shape(data_spatial) +
tm_polygons(
col = "incidence_per100k",
palette = "Reds",
title = "Insidensi HIV per 100.000 Penduduk"
) +
tm_borders() +
tm_layout(
title = "Sebaran Insidensi HIV di Provinsi Bali",
legend.outside = TRUE,
frame = FALSE
)
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_tm_polygons()`: migrate the argument(s) related to the scale of
## the visual variable `fill` namely 'palette' (rename to 'values') to fill.scale
## = tm_scale(<HERE>).[v3->v4] `tm_polygons()`: migrate the argument(s) related to the legend of the
## visual variable `fill` namely 'title' to 'fill.legend = tm_legend(<HERE>)'[v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`[cols4all] color palettes: use palettes from the R package cols4all. Run
## `cols4all::c4a_gui()` to explore them. The old palette name "Reds" is named
## "brewer.reds"Multiple palettes called "reds" found: "brewer.reds", "matplotlib.reds". The first one, "brewer.reds", is returned.

Autokorelasi Spasial Global (Moran’s I)

map_sp <- as_Spatial(data_spatial)

nb <- poly2nb(map_sp)
lw <- nb2listw(nb, style = "W")

moran.test(
data_spatial$incidence_per100k,
lw
)
## 
##  Moran I test under randomisation
## 
## data:  data_spatial$incidence_per100k  
## weights: lw    
## 
## Moran I statistic standard deviate = -0.55979, p-value = 0.7122
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.27663328       -0.16666667        0.03858906

LISA & Klaster Spasial

y <- st_drop_geometry(data_spatial)$incidence_per100k
lag_y <- lag.listw(lw, y)

mean_y <- mean(y)

data_spatial$lisa_cluster <- case_when(
y > mean_y & lag_y > mean(lag_y) ~ "High-High",
y < mean_y & lag_y < mean(lag_y) ~ "Low-Low",
y > mean_y & lag_y < mean(lag_y) ~ "High-Low",
y < mean_y & lag_y > mean(lag_y) ~ "Low-High",
TRUE ~ "Undefined"
)

table(data_spatial$lisa_cluster)
## 
## High-High  High-Low  Low-High 
##         1         3         3

Peta Klaster LISA

prev <- st_drop_geometry(data_spatial)$prevalence_per100k
tm_shape(data_spatial) +
tm_polygons(
col = "lisa_cluster",
palette = c("red", "blue", "pink", "lightblue", "gray"),
title = "Klaster Spasial HIV (LISA)"
) +
tm_borders() +
tm_layout(
legend.outside = TRUE,
frame = FALSE
)
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_tm_polygons()`: migrate the argument(s) related to the scale of
## the visual variable `fill` namely 'palette' (rename to 'values') to fill.scale
## = tm_scale(<HERE>).
## [v3->v4] `tm_polygons()`: migrate the argument(s) related to the legend of the
## visual variable `fill` namely 'title' to 'fill.legend = tm_legend(<HERE>)'

Regresi Poisson

model_pois <- glm(
cases ~ poverty + density + offset(log(population)),
family = poisson(link = "log"),
data = data
)

summary(model_pois)
## 
## Call:
## glm(formula = cases ~ poverty + density + offset(log(population)), 
##     family = poisson(link = "log"), data = data)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.173e+00  1.794e-01  -6.539 6.21e-11 ***
## poverty     -2.424e-01  3.748e-02  -6.469 9.89e-11 ***
## density      1.548e-04  1.933e-05   8.009 1.16e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 367.711  on 8  degrees of freedom
## Residual deviance:  55.404  on 6  degrees of freedom
## AIC: 113.03
## 
## Number of Fisher Scoring iterations: 4

Overdisperssion

dispersion <- model_pois$deviance / model_pois$df.residual
dispersion
## [1] 9.233999

Regresi Negative Binomial

if (dispersion > 1.5) {
model_nb <- glm.nb(
cases ~ poverty + density + offset(log(population)),
data = data
)
summary(model_nb)
}
## 
## Call:
## glm.nb(formula = cases ~ poverty + density + offset(log(population)), 
##     data = data, init.theta = 6.103213305, link = log)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)  
## (Intercept) -1.1804810  0.7394303  -1.596   0.1104  
## poverty     -0.2540977  0.1418668  -1.791   0.0733 .
## density      0.0001631  0.0001054   1.548   0.1216  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Negative Binomial(6.1032) family taken to be 1)
## 
##     Null deviance: 28.279  on 8  degrees of freedom
## Residual deviance: 10.384  on 6  degrees of freedom
## AIC: 90.691
## 
## Number of Fisher Scoring iterations: 1
## 
## 
##               Theta:  6.10 
##           Std. Err.:  3.70 
## 
##  2 x log-likelihood:  -82.691

Uji Korelasi

cor.test(
  data$poverty,
  data$incidence_per100k,
  method = "spearman"
)
## 
##  Spearman's rank correlation rho
## 
## data:  data$poverty and data$incidence_per100k
## S = 190, p-value = 0.108
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.5833333
cor.test(
  data$density,
  data$incidence_per100k,
  method = "spearman"
)
## 
##  Spearman's rank correlation rho
## 
## data:  data$density and data$incidence_per100k
## S = 36, p-value = 0.04325
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho 
## 0.7

Incidence Rate Ratio (IRR)

IRR <- tidy(model_nb, exponentiate = TRUE, conf.int = TRUE)
IRR

Tabel prevalensi tinggi

data[
  order(-data$incidence_per100k),
  c("kabkot", "cases", "incidence_per100k", "poverty", "density")
]