Regresi Poisson Kasus TBC
Alfihsyahrina Hapsery*;Ratu Bunga Prawesti Arie Salim;Saputri Dyah Pratiwi;Muhammad Alim Alhakam
2023-07-09
Kasus Tuberkulosis di Provinsi Jawa Timur pada Tahun 2020
summary(cars)## speed dist
## Min. : 4.0 Min. : 2.00
## 1st Qu.:12.0 1st Qu.: 26.00
## Median :15.0 Median : 36.00
## Mean :15.4 Mean : 42.98
## 3rd Qu.:19.0 3rd Qu.: 56.00
## Max. :25.0 Max. :120.00Packages
library(readr)
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.2 ✔ purrr 1.0.1
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ── 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 errorslibrary(DT)
library(DataExplorer)
library(Kendall)
library(vroom)##
## Attaching package: 'vroom'
##
## The following objects are masked from 'package:readr':
##
## as.col_spec, col_character, col_date, col_datetime, col_double,
## col_factor, col_guess, col_integer, col_logical, col_number,
## col_skip, col_time, cols, cols_condense, cols_only, date_names,
## date_names_lang, date_names_langs, default_locale, fwf_cols,
## fwf_empty, fwf_positions, fwf_widths, locale, output_column,
## problems, speclibrary(MASS)##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:dplyr':
##
## selectlibrary(ggplot2)
library(car)## Loading required package: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:dplyr':
##
## recode
##
## The following object is masked from 'package:purrr':
##
## someInput Data
library(readxl)
data<- read_excel("data tbc fix.xlsx")
str(data)## tibble [38 × 10] (S3: tbl_df/tbl/data.frame)
## $ No : num [1:38] 1 2 3 4 5 6 7 8 9 10 ...
## $ Kabupaten/Kota : chr [1:38] "Pacitan" "Ponorogo" "Trenggalek" "Tulungagung" ...
## $ Persentase TBC : num [1:38] 0.00688 0.02156 0.00962 0.01984 0.01494 ...
## $ TBC : num [1:38] 287 899 401 827 623 ...
## $ Usia Produktif : num [1:38] 1.46 2.11 1.22 2.55 2.8 3.9 6.54 2.62 5.99 3.97 ...
## $ Tamatan Pendidikan : num [1:38] 11.9 18.2 14.3 20.8 14.4 ...
## $ Rumah Sehat : num [1:38] 69.9 85 76.8 89.5 80.4 ...
## $ Tempat Pengelolaan Makanan: num [1:38] 36.4 53.5 85.9 78.6 63.5 41.4 47.3 54.5 71.9 65.6 ...
## $ Tempat Umum : num [1:38] 62.3 59.7 77.2 81.8 79.1 52.5 44 60.5 43.9 65.5 ...
## $ Pengeluaran Non Makanan : num [1:38] 47.3 55 50.9 51.2 52.7 ...data.poisson <- data[-1:-2:-3]## Warning in -1:-2:-3: numerical expression has 2 elements: only the first usedstr(data.poisson)## tibble [38 × 7] (S3: tbl_df/tbl/data.frame)
## $ TBC : num [1:38] 287 899 401 827 623 ...
## $ Usia Produktif : num [1:38] 1.46 2.11 1.22 2.55 2.8 3.9 6.54 2.62 5.99 3.97 ...
## $ Tamatan Pendidikan : num [1:38] 11.9 18.2 14.3 20.8 14.4 ...
## $ Rumah Sehat : num [1:38] 69.9 85 76.8 89.5 80.4 ...
## $ Tempat Pengelolaan Makanan: num [1:38] 36.4 53.5 85.9 78.6 63.5 41.4 47.3 54.5 71.9 65.6 ...
## $ Tempat Umum : num [1:38] 62.3 59.7 77.2 81.8 79.1 52.5 44 60.5 43.9 65.5 ...
## $ Pengeluaran Non Makanan : num [1:38] 47.3 55 50.9 51.2 52.7 ...Statistika Deskriptif
Summary Data
summary(data.poisson)## TBC Usia Produktif Tamatan Pendidikan Rumah Sehat
## Min. : 157.0 Min. :0.340 Min. :11.92 Min. :44.07
## 1st Qu.: 574.2 1st Qu.:1.468 1st Qu.:17.16 1st Qu.:77.84
## Median : 912.0 Median :2.395 Median :19.91 Median :84.75
## Mean :1097.2 Mean :2.630 Mean :20.42 Mean :82.07
## 3rd Qu.:1404.2 3rd Qu.:3.135 3rd Qu.:23.32 3rd Qu.:91.64
## Max. :4101.0 Max. :7.950 Max. :31.09 Max. :98.71
## Tempat Pengelolaan Makanan Tempat Umum Pengeluaran Non Makanan
## Min. : 7.80 Min. :18.80 Min. :37.10
## 1st Qu.:54.95 1st Qu.:62.30 1st Qu.:46.77
## Median :65.50 Median :68.90 Median :49.45
## Mean :62.99 Mean :68.07 Mean :49.75
## 3rd Qu.:72.12 3rd Qu.:77.75 3rd Qu.:54.63
## Max. :95.80 Max. :96.00 Max. :62.96Line Chart Perkembangan Jumlah Kasus Tuberkulosis di Provinsi Jawa Timur
library(ggplot2)
Tahun <- c("2017","2018","2019","2020")
Tahun <- as.integer(Tahun)
TBC <- c(48183,54863,48516,41693)
datatbc <- data.frame(Tahun, TBC)
datatbc## Tahun TBC
## 1 2017 48183
## 2 2018 54863
## 3 2019 48516
## 4 2020 41693lineChart <- ggplot(datatbc, aes(x = Tahun)) +
geom_line(aes(y = TBC), color = "#00ba38", size=1) ## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.lineChart
Interpretasi: Jumlah kasus tuberkulosis meningkat pada tahun 2015.Pada
tahun 2019 kasus tuberkulosis mengalami penurunan hingga tahun 2020.
Histogram Persentase Kasus Tuberkulosis di Kabupaten/Kota Provinsi Jawa Timur 2020
library(RColorBrewer)
colourCount = length(unique(data$`Kabupaten/Kota`))
getPalette = colorRampPalette(brewer.pal(9, "Paired"))
ggplot(data = data, mapping = aes(x = No, y = `Persentase TBC`)) +
geom_col(aes(fill = `Kabupaten/Kota`), alpha = 0.7) +
geom_text(aes(label = paste0(round(`Persentase TBC`, 2), "%")), vjust = 1.5,
colour = "black", size = 3) +
theme(axis.text.x=element_text(angle=50, hjust=1), legend.position="bottom") +
guides(fill=guide_legend(nrow=3)) +
labs(title = "",
x = "Kabupaten/Kota",
y = "Persentase Kasus TBC")
Interpretasi: Kabupaten/kota dengan kasus TBC tertinggi adalah kota
Surabaya sebanyak 4.101 kasus. Sedangkan kabupaten/kota dengan kasus TBC
terendah adalah Kota Batu sebanyak 157 kasus.
Boxplot Jumlah Kasus Tuberkulosis dan Persentase Tingkat Pendidikan yang Ditamatkan (SMA Sederajat)
boxplot(data.poisson$`Tamatan Pendidikan`, data = data.poisson,
horizontal = TRUE,
col = c ("lightblue"),
main = ""
)
Interpretasi: Berdasarkan boxplot di atas dapat dilihat bahwa tidak
terdapat pencilan (outlier) pada boxplot sehingga asumsi dalam analisis
regresi terpenuhi dan tidak terjadi bias.
Boxplot Jumlah Kasus Tuberkulosis dan Persentase Pengeluaran Non Makanan
boxplot(data.poisson$`Pengeluaran Non Makanan`, data = data.poisson,
horizontal = TRUE,
col = c ("lightgreen"),
main = ""
)
Interpretasi: Interpretasi: Berdasarkan boxplot di atas dapat dilihat
bahwa tidak terdapat pencilan (outlier) pada boxplot sehingga asumsi
dalam analisis regresi terpenuhi dan tidak terjadi bias.
Scatter Plot Jumlah Kasus Tuberkulosis dan Persentase Penduduk Usia Produktif (15-59 Tahun)
x1 <- data.poisson$`Usia Produktif`
y1 <- data.poisson$TBC
plot(x1, y1, main = "",
xlab = "Persentase Penduduk Usia Produktif", ylab = "Jumlah Kasus TBC",
pch = 19, frame = FALSE)
abline(lm(y1 ~ x1, data = data.poisson), col ="red")
Interpretasi: Jumlah kasus TBC (y) dan persentase penduduk usia
produktif (x1) memiliki hubungan yang positif. Peningkatan yang terjadi
pada variabel y juga diikuti peningkatan pada variabel x1. Dan jika
variabel 1 mengalami penurunan, variabel 2 juga mengalami penurunan.
Create Model Regresi Poisson
model.poisson <- glm(TBC~., family = "poisson", data = data.poisson)
summary(model.poisson)##
## Call:
## glm(formula = TBC ~ ., family = "poisson", data = data.poisson)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 7.0195473 0.0448604 156.475 < 2e-16 ***
## `Usia Produktif` 0.2700105 0.0034027 79.352 < 2e-16 ***
## `Tamatan Pendidikan` 0.0377711 0.0015926 23.717 < 2e-16 ***
## `Rumah Sehat` -0.0048959 0.0005928 -8.258 < 2e-16 ***
## `Tempat Pengelolaan Makanan` 0.0012327 0.0003866 3.188 0.00143 **
## `Tempat Umum` -0.0006261 0.0005725 -1.094 0.27410
## `Pengeluaran Non Makanan` -0.0261757 0.0013737 -19.054 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 18006.0 on 37 degrees of freedom
## Residual deviance: 2217.3 on 31 degrees of freedom
## AIC: 2558.3
##
## Number of Fisher Scoring iterations: 4Interpretasi: Berdasarkan hasil output di atas menunjukkan bahwa terdapat lima variabel memiliki nilai P-Value lebih kecil dari (alpha = 5%) yang artinya berpengaruh secara signifikan terhadap model regresi poisson.
Uji Multikolinearitas
library(car)
vif(model.poisson)## `Usia Produktif` `Tamatan Pendidikan`
## 2.355014 3.121941
## `Rumah Sehat` `Tempat Pengelolaan Makanan`
## 2.405180 1.826780
## `Tempat Umum` `Pengeluaran Non Makanan`
## 2.578281 2.980914Interpretasi: Berdasarkan hasil uji multikolinearitas, semua variabel prediktor telah memenuhi asumsi non multikolinieritas karena nilai VIF dari setiap variabel prediktor < 10. Hal ini menunjukkan bahwa tidak ada variabel prediktor yang saling berkorelasi dengan variabel prediktor lainnya.
Overdispersi
Residual deviance: 2217.3 on 31 degrees of freedom
2217.3/31## [1] 71.52581Interpretasi: Diketahui bahwa nilai residual devians sebesar 2217,3 dengan derajat bebasnya 31. Nilai hasil bagi antara residual devians dan derajat bebas adalah 71,52581. Sehingga disimpulkan bahwa pada model regresi poisson untuk jumlah kasus tuberkulosis di Jawa Timur pada tahun 2020 terdapat kasus overdispersi. Karena terdapat overdispersi pada model regresi poisson, maka dilanjutkan untuk melakukan pemodelan pada regresi binomial negatif.
Alternatif Model (Regresi Binomial Negatif)
model.nb <- glm.nb(TBC~., data.poisson)
summary(model.nb)##
## Call:
## glm.nb(formula = TBC ~ ., data = data.poisson, init.theta = 13.34753705,
## link = log)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.828e+00 4.109e-01 16.619 < 2e-16 ***
## `Usia Produktif` 3.212e-01 3.046e-02 10.544 < 2e-16 ***
## `Tamatan Pendidikan` 4.178e-02 1.357e-02 3.080 0.00207 **
## `Rumah Sehat` -6.149e-03 4.976e-03 -1.236 0.21656
## `Tempat Pengelolaan Makanan` 8.211e-05 3.155e-03 0.026 0.97924
## `Tempat Umum` 1.212e-03 4.764e-03 0.254 0.79916
## `Pengeluaran Non Makanan` -2.616e-02 1.152e-02 -2.270 0.02318 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(13.3475) family taken to be 1)
##
## Null deviance: 231.31 on 37 degrees of freedom
## Residual deviance: 38.68 on 31 degrees of freedom
## AIC: 541.69
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 13.35
## Std. Err.: 3.10
##
## 2 x log-likelihood: -525.688Interpretasi: Berdasarkan hasil output di atas menunjukkan bahwa terdapat tiga variabel memiliki nilai P-Value lebih kecil dari (alpha = 5%) yang artinya berpengaruh secara signifikan terhadap model regresi binomial negatif.
Alternatif Model (Regresi Binomial Negatif dengan Seleksi Variabel)
model.nb2 <- glm.nb(TBC~`Usia Produktif`+`Tamatan Pendidikan`+`Pengeluaran Non Makanan`, data.poisson)
summary(model.nb2)##
## Call:
## glm.nb(formula = TBC ~ `Usia Produktif` + `Tamatan Pendidikan` +
## `Pengeluaran Non Makanan`, data = data.poisson, init.theta = 12.66081794,
## link = log)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.704699 0.399761 16.772 < 2e-16 ***
## `Usia Produktif` 0.322546 0.027191 11.862 < 2e-16 ***
## `Tamatan Pendidikan` 0.037099 0.012599 2.945 0.00323 **
## `Pengeluaran Non Makanan` -0.030161 0.009727 -3.101 0.00193 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(12.6608) family taken to be 1)
##
## Null deviance: 219.595 on 37 degrees of freedom
## Residual deviance: 38.694 on 34 degrees of freedom
## AIC: 537.7
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 12.66
## Std. Err.: 2.93
##
## 2 x log-likelihood: -527.697Interpretasi: Berdasarkan hasil output di atas menunjukkan bahwa terdapat tiga variabel memiliki nilai P-Value lebih kecil dari (alpha = 5%) yang artinya berpengaruh secara signifikan terhadap model regresi binomial negatif dengan seleksi variabel.
Perbandingan Nilai AIC
AIC(model.poisson) # Model regresi poisoon## [1] 2558.329AIC(model.nb) # Model regresi binomial negatif## [1] 541.6881AIC(model.nb2) # Model regresi binomial negatif dengan seleksi variabel## [1] 537.697Interpretasi: Model regresi poisson dengan seleksi variabel memiliki nilai AIC paling kecil, sehingga model regresi binomial negatif dengan seleksi variabel (AIC = 537,697) sesuai untuk memodelkan faktor-faktor yang mempengaruhi jumlah kasus tuberkulosis setiap kabupaten/kota di Provinsi Jawa Timur 2020.