QUANTILE REGRESSION BELANJA KESEHATAN OUT OF POCKET (OOP) RUMAH TANGGA TAHUN 2023
Sapta Hastho Ponco
2024-04-14
Regresi Kuantil merupakan suatu pendekatan dalam analisis regresi yang dikenalkan oleh Koenker dan Bassett (1978). Pendekatan ini menduga berbagai fungsi kuantil dari suatu distribusi Y sebagai fungsi dari X. Regresi Kuantil sangat berguna jika distribusi data tidak homogen (heterogenous) dan tidak berbentuk standard seperti tidak simetris, terdapat kemencengan, atau truncated distribution. Seperti halnya dengan metode OLS, Regresi Kuantil juga meminimumkan jumlah kuadrat residual untuk mencari nilai dugaan bagi β.
Dalam tulisan ini, ingin dilihat bagaimana belanja kesehatan OOP rumah tangga di Provinsi Jawa Timur dipengaruhi dari beberapa determinan seperti pengeluaran per kapita per bulan, pendidikan tertinggi dari Kepala Rumah Tangga (KRT), Kepemilikan Jaminan Kesehatan, Jenis Kelamin KRT dan Pekerjaan KRT yang diperoleh dari Susenas Maret 2023.
Variabel Pekerjaan KRT terdiri atas: 0. Tidak bekerja 1. Berusaha sendiri 2. Pengusaha 3. Karyawan 4. Pekerja bebas, Pekerja keluarga
Variabel Pendidikan Tertinggi KRT terdiri atas: 1. Tdk/blm pernah sekolah+Tidak Tamat SD 2. SD Sederajat 3. SMP Sederajat 4. SM Sederajat 5. Perguruan Tinggi
Variabel Kepemilikan Jamkes terdiri atas: 1. Seluruh ART memiliki JKN 2. Seluruh ART memiliki jamkes selain JKN 3. Sebagian ART memiliki JKN 4. Tidak memiliki jaminan kesehatan apapun
1 Load Packages
pacman::p_load(readxl, base, lmtest, quantreg, car, stats, nortest, ggpubr)Setting Working Directory
setwd("D:\\Materi\\")2 Load Dataset
Dataset yang digunakan adalah Susenas Maret 2023.
dataku <- read_excel("Ruta_OOP 2023 - 35.xls")
dataku$LnRutaOOP <- log(dataku$Rata_set_OOP_Kesehatan)
dataku$LnKapita <- log(dataku$KAPITA)
dataku$JK_KRT <- ifelse(dataku$JK_KRT==1, 1, 0)
dataku <- subset(dataku, select=-c(KABU,KAPITA, MISKIN,Rata_set_OOP_Kesehatan))
head(dataku)## # A tibble: 6 × 6
## TAMAT_KRT JK_KRT PEKERJAAN_KRT JAMKES LnRutaOOP LnKapita
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 2 0 2 4 12.0 14.4
## 2 5 1 0 1 12.8 14.1
## 3 5 1 3 1 12.5 14.2
## 4 4 1 3 3 13.2 13.8
## 5 2 1 3 4 10.2 13.4
## 6 5 1 1 3 12.3 13.3#Data Table
#DT::datatable(dataku)3 Statistik Deskriptif
summary(dataku)## TAMAT_KRT JK_KRT PEKERJAAN_KRT JAMKES
## Min. :1.000 Min. :0.0000 Min. :0.000 Min. :1.000
## 1st Qu.:2.000 1st Qu.:1.0000 1st Qu.:1.000 1st Qu.:1.000
## Median :3.000 Median :1.0000 Median :2.000 Median :1.000
## Mean :2.734 Mean :0.8319 Mean :2.125 Mean :1.953
## 3rd Qu.:4.000 3rd Qu.:1.0000 3rd Qu.:3.000 3rd Qu.:3.000
## Max. :5.000 Max. :1.0000 Max. :4.000 Max. :4.000
## LnRutaOOP LnKapita
## Min. : 6.502 Min. :12.20
## 1st Qu.:10.779 1st Qu.:13.47
## Median :11.538 Median :13.85
## Mean :11.601 Mean :13.95
## 3rd Qu.:12.357 3rd Qu.:14.32
## Max. :18.331 Max. :17.74Variabel LnRutaOOP dan LnKapita tersebar dengan median yang tidak sama dengan rata-rata, menunjukkan distribusi yang tidak simetris.
Lebih lanjut diberikan visualisasi melalui histogram seperti dibawah ini:
par(mfrow=c(3,3),oma=c(1,1,1,1))
hist(dataku$TAMAT_KRT,col = "lightblue",xlab="TAMAT_KRT", ylab="Frequency", main = "Tingkat Pendidikan KRT")
hist(dataku$JK_KRT,col = "lightblue",xlab="JK_KRT", ylab="Frequency", main = "Jenis Kelamin KRT")
hist(dataku$JAMKES,col = "lightblue",xlab="JAMKES", ylab="Frequency", main = "Kepemilikan Jamkes")
hist(dataku$LnKapita,col = "lightblue",xlab="LnKapita", ylab="Frequency", main = "Ln(Kapita)")
hist(dataku$PEKERJAAN_KRT,col = "lightblue",xlab="PEKERJAAN_KRT", ylab="Frequency", main = "Status Pekerjaan KRT")
hist(dataku$LnRutaOOP,col = "lightblue",xlab="LnRutaOOP", ylab="Frequency", main = "Ln(Ruta OOP)")
Korelasi Pearson
datatable=data.frame(dataku$LnRutaOOP, dataku$LnKapita)
cor(datatable)## dataku.LnRutaOOP dataku.LnKapita
## dataku.LnRutaOOP 1.0000000 0.5049773
## dataku.LnKapita 0.5049773 1.0000000pairs(datatable, col="blue", main="Scatterplots")
Berdasarkan korelasi Pearson, LnKapita erat berhubungan dengan LnRutaOOP dengan arah yang sama.
Korelasi Rank-Spearman
korelasi_spearman <- cor(data.frame(dataku$TAMAT_KRT,dataku$JK_KRT,dataku$JAMKES,dataku$LnRutaOOP),method = "spearman")
cor(korelasi_spearman)## dataku.TAMAT_KRT dataku.JK_KRT dataku.JAMKES dataku.LnRutaOOP
## dataku.TAMAT_KRT 1.00000000 -0.01062186 -0.5792738 -0.1123224
## dataku.JK_KRT -0.01062186 1.00000000 -0.2613085 -0.5195133
## dataku.JAMKES -0.57927378 -0.26130845 1.0000000 -0.4472974
## dataku.LnRutaOOP -0.11232238 -0.51951335 -0.4472974 1.0000000Dari korelasi Rank-Spearman, Jenis Kelamin KRT dan Kepemilikan Jaminan Kesehatan cukup erat berhubungan dengan LnRutaOOP namun arahnya berbanding terbalik.
4 Regression Ordinary Least Square (OLS)
Metode OLS adalah salah satu metode yang terdapat dalam analisis regresi berganda. Metode ini digunakan dengan tujuan untuk meminimalisir jumlah kuadrat residual (Sum of Squares Error) dengan mengestimasi suatu garis regresi. Metode OLS termasuk jenis metode ekonometrik dengan 2 variabel yaitu variabel independen dan variabel dependen.
Hasil akhir yang dimiliki oleh metode OLS adalah fungsi regresi populasi yang akan digunakan untuk estimasi data. Untuk menghasilkan estimasi menggunakan metode OLS, diperlukan empat asumsi dasar yang bersifat BLUE (Best, Linear, Unbiased, dan Estimator).
Metode OLS juga memiliki kriteria operasi yaitu line of best fit dengan jumlah kuadrat deviasi antara titik observasi dengan garis regresi adalah minimum. Hasil estimasi dalam metode OLS memiliki sifat yang BLUE dan cenderung lebih efisien, konsisten.
Selain itu, hasil estimasi pada metode OLS juga cenderung memiliki koefisien regresi dengan distribusi yang normal. Intersep nilai dependen dari estimasi metode OLS juga cenderung memiliki distribusi yang normal jika nilai independen-nya nol.
Dalam metode OLS terdapat beberapa uji asumsi seperti Normalitas, Homoskedastisitas, Non Multikolinearitas.
Berikut dilakukan pemodelan OLS untuk memeriksa apakah ada pelanggaran pada asumsi-asumsi di dalamnya.
reg <- lm(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, data=dataku)
summary(reg)##
## Call:
## lm(formula = LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT +
## JK_KRT, data = dataku)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.8913 -0.6635 0.0043 0.6343 4.4814
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.861379 0.150812 -18.973 < 2e-16 ***
## JAMKES 0.020391 0.005188 3.931 8.49e-05 ***
## LnKapita 1.054845 0.010756 98.068 < 2e-16 ***
## TAMAT_KRT -0.062874 0.005528 -11.373 < 2e-16 ***
## PEKERJAAN_KRT -0.034450 0.005536 -6.223 4.94e-10 ***
## JK_KRT -0.060630 0.018359 -3.302 0.000959 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.144 on 30425 degrees of freedom
## Multiple R-squared: 0.2611, Adjusted R-squared: 0.2609
## F-statistic: 2150 on 5 and 30425 DF, p-value: < 2.2e-16Dari model OLS di atas, tampak seluruh variabel bebas signifikan terhadap variabel respon LnRutaOOP. Semua variabel bebas memiliki nilai p-value < α = 0.05.
Pemeriksaan Asumsi Normalitas
Pemeriksaan asumsi normalitas dilakukan dengan menggunakan Uji Anderson-Darling.
ad.test(reg$residuals)##
## Anderson-Darling normality test
##
## data: reg$residuals
## A = 100.74, p-value < 2.2e-16Hasil dari pengujian diatas memberikan kesimpulan bahwa terjadi pelanggaran asumsi normalitas. Nilai p-value < α = 0.05 (tolak \(H_{0}\)).
ggqqplot(dataku$LnKapita, main = "Ln(Ruta OOP)")
Berdasarkan hasil plot diatas tampak bahwa variabel Pengeluaran Kesehatan OOP rumah tangga tidak berada pada garis diagonal. Ini menunjukkan bahwa data tersebut tidak mengikuti distribusi normal.
Pemeriksaan Asumsi Non-Multikolinearitas
Multikolinearitas adalah suatu kondisi adanya hubungan yang linear diantara variabel-variabel bebas dalam model regresi. Multikolinearitas dapat dideteksi menggunakan Variance Inflation Factors (VIF). Nilai VIF dapat dihitung dengan rumus: \(VIF_{j}=\frac{1}{1-R_{j}^{2}}\) dengan \(R_{j}\) merupakan koefisien determinasi ke-j, j=1,2,…,k.
vif(reg)## JAMKES LnKapita TAMAT_KRT PEKERJAAN_KRT JK_KRT
## 1.035642 1.181214 1.179285 1.070945 1.096625# Tolerance
#1/vif(reg)Nilai VIF dari seluruh variabel bebas < 10 sehingga model sudah terbebas dari multikolinearitas.
Pemeriksaan Asumsi Non-Heteroskedastisitas (Homoskedastisitas)
Heteroskedastisitas terjadi jika varians dari residual suatu pengamatan ke pengamatan lain terjadi ketidaksamaan (tidak konstan). Bila dalam suatu pengamatan terdapat kasus heteroskedastisitas sedangkan asumsi lain dipenuhi, maka parameter regresi yang diduga dengan metode OLS tidak lagi memiliki sifat efisien (varians minimum) (Draper and Smith, 1998).
Hipotesis
\(H_{0}\): \(Var[ε]\) = \(σ^{2}I\) (varians residual bersifat homogen)
\(H_{1}\): \(Var[ε]\) ≠ \(σ^{2}I\) (varians residual tidak homogen)
Pemeriksaan homoskedastisitas dilakukan dengan uji Breusch-Pagan.
bptest(reg)##
## studentized Breusch-Pagan test
##
## data: reg
## BP = 418.29, df = 5, p-value < 2.2e-16Hasil dari pengujian diatas memberikan kesimpulan bahwa terjadi pelanggaran asumsi homoskedastisitas. Ini ditunjukkan dari nilai p-value < α = 0.05 (tolak \(H_{0}\)).
5 Quantile Regression
Hasil pengujian menunjukkan asumsi homoskedastisitas terlanggar. Dengan kondisi ini, apabila analisis regresi dilakukan dengan menggunakan metode OLS, maka akan mengakibatkan varians dari penduga parameter besar sehingga penduga parameternya menjadi tidak efisien. Untuk mengatasi hal tersebut dilakukan analisis dengan metode regresi kuantil.
Regresi kuantil dilakukan pada Quantile 10% (τ=0.10), Quantile 50% (τ=0.50) dan Quantile 90% (τ=0.90)
quantreg10 <- rq(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, data=dataku, tau=0.10)
summary(quantreg10)##
## Call: rq(formula = LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT +
## JK_KRT, tau = 0.1, data = dataku)
##
## tau: [1] 0.1
##
## Coefficients:
## Value Std. Error t value Pr(>|t|)
## (Intercept) -2.99148 0.30982 -9.65540 0.00000
## JAMKES 0.00757 0.01050 0.72049 0.47123
## LnKapita 0.95595 0.02201 43.42280 0.00000
## TAMAT_KRT -0.08451 0.01156 -7.30909 0.00000
## PEKERJAAN_KRT 0.01344 0.01158 1.16034 0.24592
## JK_KRT 0.05784 0.03754 1.54069 0.12340Persamaan regresi untuk kuantil 10% adalah
Y10% = -2.9915 + 0.0076(JAMKES) + 0.9559(LnKapita) - 0.0845(TAMAT_KRT) + 0.0134(PEKERJAAN_KRT) + 0.0578(JK_KRT)
quantreg50 <- rq(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, data=dataku, tau=0.50)
summary(quantreg50)##
## Call: rq(formula = LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT +
## JK_KRT, tau = 0.5, data = dataku)
##
## tau: [1] 0.5
##
## Coefficients:
## Value Std. Error t value Pr(>|t|)
## (Intercept) -2.34588 0.16178 -14.50052 0.00000
## JAMKES 0.01157 0.00508 2.27640 0.02283
## LnKapita 1.01713 0.01160 87.67423 0.00000
## TAMAT_KRT -0.05485 0.00555 -9.88706 0.00000
## PEKERJAAN_KRT -0.02032 0.00552 -3.68147 0.00023
## JK_KRT -0.08656 0.02064 -4.19418 0.00003Persamaan regresi untuk kuantil 50% adalah
Y50% = -2.3459 + 0.0116(JAMKES) + 1.0171(LnKapita) - 0.0548(TAMAT_KRT) - 0.0203(PEKERJAAN_KRT) - 0.0866(JK_KRT)
quantreg90 <- rq(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, data=dataku, tau=0.90)
summary(quantreg90)##
## Call: rq(formula = LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT +
## JK_KRT, tau = 0.9, data = dataku)
##
## tau: [1] 0.9
##
## Coefficients:
## Value Std. Error t value Pr(>|t|)
## (Intercept) -3.78452 0.34463 -10.98153 0.00000
## JAMKES 0.06645 0.01142 5.81789 0.00000
## LnKapita 1.21900 0.02481 49.12491 0.00000
## TAMAT_KRT -0.02682 0.01193 -2.24844 0.02456
## PEKERJAAN_KRT -0.10734 0.01172 -9.15791 0.00000
## JK_KRT -0.07909 0.03751 -2.10842 0.03500Persamaan regresi untuk kuantil 90% adalah
Y90% = -3.7845 + 0.0664(JAMKES) + 1.2190(LnKapita) - 0.0268(TAMAT_KRT) - 0.1073(PEKERJAAN_KRT) - 0.0791(JK_KRT)
Contoh interpretasi untuk variabel LnKapita adalah sebagai berikut:
rumah tangga yang memiliki pengeluaran per kapita per bulan pada kuantil ke-1 akan menaikkan 0.9559 pengeluaran kesehatan OOP rumah tangga mereka, sedangkan rumah tangga pada kuantil ke-5 lebih tinggi menaikkan pengeluaran kesehatan OOP rumah tangga mereka (sebesar 1.2190). Sederhananya, pengeluaran rumah tangga yang memiliki OOP akan meningkat seiring dengan meningkatnya pengeluaran per kapita per bulan.
95% CI untuk Variabel LnKapita
Estimasi β LnKapita pada OLS = 1.0548
lower10 <- 0.95595 - 1.96*0.02201
upper10 <- 0.95595 + 1.96*0.02201
lower50 <- 1.01713 - 1.96*0.01160
upper50 <- 1.01713 + 1.96*0.01160
lower90 <- 1.21900 - 1.96*0.02481
upper90 <- 1.21900 + 1.96*0.02481
BatasBawah <- c(lower10, lower50, lower90)
BatasAtas <- c(upper10, upper50, upper90)
CI <- data.frame(BatasBawah, BatasAtas)
CI## BatasBawah BatasAtas
## 1 0.9128104 0.9990896
## 2 0.9943940 1.0398660
## 3 1.1703724 1.2676276Terlihat bahwa estimasi β LnKapita dari output regresi OLS berada di luar CI untuk kuantil ke-10, ke-50, dan ke-90. Oleh karena itu dari kuantil-kuantil yang ada, estimasi β variabel LnKapita berbeda secara signifikan dari model OLS.
Quantile Regression Plots
quantregall<-rq(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, tau=seq(0.10, 0.90, by=0.05), data=dataku)
quantregplot<-summary(quantregall)
plot(quantregplot)
Dapat dilihat secara visual bahwa kesimpulan yang dibuat dari interval kepercayaan untuk Variabel LnKapita mencerminkan plot kuantil.
60th Quantile for LnKapita
Estimasi β LnKapita pada OLS = 1.0548
quantreg60 <- rq(LnRutaOOP~ JAMKES+LnKapita+TAMAT_KRT+PEKERJAAN_KRT+JK_KRT, data=dataku, tau=0.60)
summary(quantreg60)##
## Call: rq(formula = LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT +
## JK_KRT, tau = 0.6, data = dataku)
##
## tau: [1] 0.6
##
## Coefficients:
## Value Std. Error t value Pr(>|t|)
## (Intercept) -2.44089 0.15853 -15.39702 0.00000
## JAMKES 0.01165 0.00511 2.28153 0.02252
## LnKapita 1.04426 0.01132 92.20902 0.00000
## TAMAT_KRT -0.05726 0.00559 -10.24102 0.00000
## PEKERJAAN_KRT -0.03167 0.00559 -5.66702 0.00000
## JK_KRT -0.11982 0.02075 -5.77497 0.00000lower60 <- 1.04426 - 1.96*0.01132
upper60 <- 1.04426 + 1.96*0.01132
BatasBawah <- c(lower60)
BatasAtas <- c(upper60)
CI <- data.frame(BatasBawah, BatasAtas)
CI## BatasBawah BatasAtas
## 1 1.022073 1.066447Dapat dilihat bahwa koefisien β pada OLS untuk LnKapita sebesar 1.0548 berada di dalam Interval Kepercayaan kuantil ke-60, yang berarti tidak ada perbedaan signifikan dalam β antara τ=0,60 dan OLS.
Testing for difference in coefficients across different quantiles
anova(quantreg10, quantreg50, quantreg90)## Quantile Regression Analysis of Deviance Table
##
## Model: LnRutaOOP ~ JAMKES + LnKapita + TAMAT_KRT + PEKERJAAN_KRT + JK_KRT
## Joint Test of Equality of Slopes: tau in { 0.1 0.5 0.9 }
##
## Df Resid Df F value Pr(>F)
## 1 10 91283 24.669 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Lebih lanjut penguatan penggunaan regresi kuantil adalah dengan menguji perbedaan signifikan antara koefisien dari kuantil yang berbeda. Tabel Anova di atas adalah uji perbedaan signifikan antara kemiringan untuk kuantil ke-10, ke-50, dan ke-90. Dari tabel tersebut, terlihat p<0,0001, yang berarti menolak hipotesis nol bahwa kemiringan sama dan disimpulkan bahwa setidaknya ada yang satu berbeda.
6 Kesimpulan
Metode Regresi Kuantil merupakan cara estimasi dengan menggunakan pendekatan fungsi kuantil dari distribusi Y sebagai fungsi dari variabel bebas X. Regresi Kuantil sangat berguna untuk data yang tidak simetris atau tidak homogen (heteroskedaktisitas). Nilai dugaan parameter dengan regresi kuantil memberikan hasil yang berbeda dengan analisis regresi yang menggunakan metode OLS.
7 References
- Draper, N.R. and Smith, H. 1998. Applied Regression Analysis. New York: John Wiley and sons.
- Koenker, R and G. Basset. 1978. Regression Quantile Econometrica, 46, 33-50.
Direktorat Statistik Kesejahteraan Rakyat, BPS, saptahas@bps.go.id