LATAR BELAKANG
PACKAGES
library(readxl) #Membaca file data excel
library(plm) #untuk membuat model
library(kableExtra) #untuk tampilan tabel
library(lmtest) #uji homoskedastisitas
library(RMySQL)
library(tidyverse)DATA
data_adp <- read_excel("C:/Users/HP/Downloads/Kelompok 6 ADP/TUGAS ADP GABUNGAN/Dataset Kelompok 6 ADP.xlsx", sheet = "Data FIX")
head(data_adp)## # A tibble: 6 x 7
## `Kab/Kota` Tahun TPT RLS TPAK PDRB PPM
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Bogor 2019 9.11 8.82 65.4 5.85 6.66
## 2 Bogor 2020 14.3 8.55 62.7 -1.76 7.69
## 3 Bogor 2021 12.2 8.74 62.6 3.55 8.13
## 4 Sukabumi 2019 8.05 7.65 62.6 5.64 6.22
## 5 Sukabumi 2020 9.6 7.82 61.6 -0.91 7.09
## 6 Sukabumi 2021 9.51 7.71 64.9 3.74 7.7DATA PANEL
datapanel <- pdata.frame(data_adp)
summary(datapanel)## Kab.Kota Tahun TPT RLS TPAK
## Bandung : 3 2019:27 Min. : 3.250 Min. : 6.610 Min. :55.74
## Bandung Barat: 3 2020:27 1st Qu.: 7.700 1st Qu.: 7.959 1st Qu.:62.65
## Bekasi : 3 2021:27 Median : 9.210 Median : 8.553 Median :64.74
## Bogor : 3 Mean : 9.077 Mean : 8.914 Mean :65.07
## Ciamis : 3 3rd Qu.:10.880 3rd Qu.: 9.797 3rd Qu.:67.39
## Cianjur : 3 Max. :14.290 Max. :11.628 Max. :76.79
## (Other) :63
## PDRB PPM
## Min. :-3.80 Min. : 2.070
## 1st Qu.:-0.77 1st Qu.: 6.650
## Median : 3.56 Median : 8.260
## Mean : 2.59 Mean : 8.267
## 3rd Qu.: 5.03 3rd Qu.:10.340
## Max. : 7.85 Max. :13.130
## glimpse(datapanel)## Rows: 81
## Columns: 7
## $ Kab.Kota <fct> Bandung, Bandung, Bandung, Bandung Barat, Bandung Barat, Band~
## $ Tahun <fct> 2019, 2020, 2021, 2019, 2020, 2021, 2019, 2020, 2021, 2019, 2~
## $ TPT <pseries> 5.51, 8.58, 8.32, 8.24, 12.25, 11.65, 9.00, 11.54, 10.09,~
## $ RLS <pseries> 9.128577, 9.268625, 9.404190, 8.712663, 8.671729, 8.39515~
## $ TPAK <pseries> 65.31591, 62.19562, 65.12070, 61.96509, 59.91116, 60.7499~
## $ PDRB <pseries> 6.36, -1.80, 3.56, 5.05, -2.41, 3.46, 3.95, -3.39, 3.62, ~
## $ PPM <pseries> 5.94, 6.91, 7.15, 9.38, 10.49, 11.30, 4.01, 4.82, 5.21, 6~Check Balance
datapanel %>%
is.pbalanced()## [1] TRUECek Dimensi Waktu
datapanel %>%
is.pconsecutive()## Bandung Bandung Barat Bekasi Bogor
## TRUE TRUE TRUE TRUE
## Ciamis Cianjur Cirebon Garut
## TRUE TRUE TRUE TRUE
## Indramayu Karawang Kota Bandung Kota Banjar
## TRUE TRUE TRUE TRUE
## Kota Bekasi Kota Bogor Kota Cimahi Kota Cirebon
## TRUE TRUE TRUE TRUE
## Kota Depok Kota Sukabumi Kota Tasikmalaya Kuningan
## TRUE TRUE TRUE TRUE
## Majalengka Pangandaran Purwakarta Subang
## TRUE TRUE TRUE TRUE
## Sukabumi Sumedang Tasikmalaya
## TRUE TRUE TRUEEksplorasi Data
ggplot(data = datapanel, aes(x = Tahun, y = TPT)) +
geom_line() +
labs(x = "Tahun", y = "Tingkat Pengangguran Terbuka") +
theme(legend.position = "none")+ theme_bw()
ggplot() + geom_point(aes(datapanel$TPT, datapanel$RLS),colour = 'red') + ggtitle('Tingkat Pengangguran Terbuka Vs Rata-Rata Lama Sekolah di Atas 15 Tahun') +
xlab('Pengangguran') + ylab('Rata-Rata Lama Sekolah')
ggplot() + geom_point(aes(datapanel$TPT, datapanel$TPAK),colour = 'red') + ggtitle('Tingkat Pengangguran Terbuka Vs Tingkat Partisipasi Angkatan Kerja') +
xlab('Pengangguran') + ylab('Partisipasi Angkatan Kerja')
ggplot() + geom_point(aes(datapanel$TPT, datapanel$PDRB),colour = 'red') + ggtitle('Tingkat Pengangguran Terbuka Vs Produk Domestik Regional Bruto') +
xlab('Pengangguran') + ylab('PDRB')
ggplot() + geom_point(aes(datapanel$TPT, datapanel$PPM),colour = 'red') + ggtitle('Tingkat Pengangguran Terbuka Vs Persentase Penduduk Miskin') +
xlab('Pengangguran') + ylab('Penduduk Miskin')
library(gplots)##
## Attaching package: 'gplots'## The following object is masked from 'package:stats':
##
## lowessplotmeans(TPT ~ `Kab/Kota`, main="Keragaman Antar Individu",data_adp)
plotmeans(TPT ~ Tahun, main="Keragaman Antar Waktu",data_adp)
Common Effect Model (CEM)
Common effect model (CEM) adalah suatu model dalam statistik inferensial yang digunakan untuk menguji hubungan antara satu variabel independen dengan satu variabel dependen, ketika terdapat satu atau lebih variabel mediasi yang dapat mempengaruhi hubungan tersebut. Dalam CEM, variabel mediasi tersebut diperlakukan sebagai variabel bebas (independent variable) dan digunakan untuk menjelaskan hubungan antara variabel independen dan variabel dependen.
Sebagai contoh, misalkan seorang peneliti ingin menguji hubungan antara tingkat pendidikan seseorang dengan penghasilannya. Namun, ia juga menyadari bahwa variabel mediasi seperti pengalaman kerja dan jenis pekerjaan dapat mempengaruhi hubungan ini. Dalam hal ini, peneliti dapat menggunakan CEM untuk memeriksa apakah pengalaman kerja dan jenis pekerjaan memediasi hubungan antara tingkat pendidikan dan penghasilan.
Dalam CEM, variabel mediasi digunakan sebagai variabel bebas dalam analisis regresi dan pengujian statistik. Peneliti dapat menguji apakah variabel mediasi mempengaruhi hubungan antara variabel independen dan variabel dependen dengan menghitung nilai koefisien jalur (path coefficient) dari masing-masing variabel mediasi dalam model. Jika koefisien jalur variabel mediasi signifikan, maka hal tersebut menunjukkan bahwa variabel mediasi mempengaruhi hubungan antara variabel independen dan variabel dependen.
cem <- plm(TPT ~ RLS + TPAK + PDRB + PPM, data=data_adp, model = "pooling")
summary(cem)## Pooling Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = data_adp,
## model = "pooling")
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -3.084195 -1.377728 0.088007 0.997320 3.482821
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## (Intercept) 32.942443 4.637132 7.1041 5.588e-10 ***
## RLS 0.211543 0.199544 1.0601 0.2924
## TPAK -0.392173 0.053028 -7.3956 1.566e-10 ***
## PDRB -0.285390 0.057823 -4.9356 4.604e-06 ***
## PPM 0.061240 0.087130 0.7029 0.4843
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 440.77
## Residual Sum of Squares: 191.5
## R-Squared: 0.56552
## Adj. R-Squared: 0.54266
## F-statistic: 24.7309 on 4 and 76 DF, p-value: 3.9325e-13Berdasarkan output model CEM di atas, didapatkan peubah yang tidak signifikan terhadap model, yaitu peubah RLS dan PPM. Selain itu, didapatkan juga koefisien determinasi dari model sebesar 56,55%. Hal tersebut menandakan bahwa keragaman yang dapat dijelaskan oleh model tergolong cukup besar sehingga model CEM sudah baik digunakan untuk memodelkan pengangguran di Indonesia.
Uji Normalitas
Dalam pengujian ini, digunakan Uji Jarque Bera dengan hipotesis
sebagai berikut:
H0 = Sisaan menyebar normal
H1 = Sisaan tidak
menyebar normal
library(tseries)
res.cem <- residuals(cem)(normal <- jarque.bera.test(res.cem))##
## Jarque Bera Test
##
## data: res.cem
## X-squared = 2.2328, df = 2, p-value = 0.3275Histogram
hist(res.cem,
xlab = "Sisaan",
col = "#27D3D3",
breaks=30,
prob = TRUE)
lines(density(res.cem), # density plot
lwd = 2, # thickness of line
col = "chocolate3")
Plot QQ-Norm
set.seed(1234)
res.cem1 <- as.numeric(res.cem)
qqnorm(res.cem1,datax=T, col="blue")
qqline(rnorm(length(res.cem1),mean(res.cem1),sd(res.cem1)),datax=T, col="red")
Uji Homoskedastisitas
Hipotesis:
H0 = Sisaan memiliki ragam homogen
H1 = Sisaan
tidak memiliki ragam homogen
(homos <- bptest(cem))##
## studentized Breusch-Pagan test
##
## data: cem
## BP = 1.0464, df = 4, p-value = 0.9027Uji Autokorelasi
Hipotesis:
H0 = Sisaan saling bebas
H1 = Sisaan tidak
saling bebas
(autokol <- pbgtest(cem))##
## Breusch-Godfrey/Wooldridge test for serial correlation in panel models
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 13.676, df = 3, p-value = 0.003381
## alternative hypothesis: serial correlation in idiosyncratic errorsdata.frame(Asumsi=c("Normalitas","Homoskesdastisitas","Non-Autokorelasi"),
PValue=c(normal$p.value,homos$p.value,autokol$p.value),
Keputusan= c("Terima H0","Terima H0","Tolak H0"),
Kesimpulan= c("Sisaan menyebar normal",
"Ragam sisaan homogen",
"Terdapat autokorelasi")) ## Asumsi PValue Keputusan Kesimpulan
## 1 Normalitas 0.32745172 Terima H0 Sisaan menyebar normal
## 2 Homoskesdastisitas 0.90267729 Terima H0 Ragam sisaan homogen
## 3 Non-Autokorelasi 0.00338087 Tolak H0 Terdapat autokorelasiBerdasarkan hasil uji diagnostik sisaan model CEM, didapatkan bahwa pengujian ketiga asumsi menghasilkan satu nilai p-value < 0.05, sehingga keputusannya adalah tolak H0 pada uji Autokorelasi. Dengan kata lain, asumsi normalitas terpenuhi dan memiliki ragam yang homogen serta sisaan tidak saling bebas.
Fixed Effect Model Within (FEM Within)
Fixed effect model adalah metode analisis data dalam statistik dan ekonometrika yang digunakan untuk memperkirakan hubungan antara variabel independen dan dependen dengan memperhitungkan efek tetap dari setiap individu atau unit pengamatan dalam sampel.
Dalam konteks penelitian seseorang, model efek tetap dapat digunakan untuk menganalisis data panel (data yang mengikuti subjek/subyek dari waktu ke waktu) atau data cross-section (data yang diambil pada satu titik waktu tertentu) untuk memeriksa hubungan antara variabel dependen dan independen dengan memperhitungkan perbedaan antar individu dalam karakteristik atau faktor lain yang mungkin memengaruhi hubungan tersebut.
Dalam model efek tetap, efek individu yang konstan dari setiap unit pengamatan dihitung dan dieliminasi, sehingga hanya variabel independen yang berubah dari waktu ke waktu atau antar individu yang dianggap sebagai faktor yang mempengaruhi variabel dependen. Metode ini sering digunakan dalam ekonomi, sosiologi, ilmu politik, dan bidang lain di mana data panel atau cross-section digunakan untuk menganalisis hubungan antara variabel.
fem.within <- plm(TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
model = "within",
index = c("Kab.Kota","Tahun"))
summary(fem.within)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## model = "within", index = c("Kab.Kota", "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.871451 -0.503443 0.041329 0.421842 1.398651
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS 0.446339 0.618656 0.7215 0.473983
## TPAK -0.212108 0.079037 -2.6837 0.009849 **
## PDRB -0.219655 0.035480 -6.1910 1.104e-07 ***
## PPM 0.429229 0.177325 2.4206 0.019170 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 104.68
## Residual Sum of Squares: 37.331
## R-Squared: 0.64337
## Adj. R-Squared: 0.42939
## F-statistic: 22.55 on 4 and 50 DF, p-value: 1.0918e-10Berdasarkan output di atas, model FEM Within hanya menghasilkan satu peubah yang tidak signifikan terhadap model. Selain itu, didapatkan juga koefisien determinasi dari model sebesar 64,34%. Hal tersebut menandakan bahwa keragaman yang dapat dijelaskan oleh model tergolong besar sehingga model FEM Within sudah baik digunakan untuk memodelkan kemiskinan di Indonesia.
res.fem.within <- residuals(fem.within)Uji Normalitas
Hipotesis:
H0 = Sisaan menyebar normal
H1 = Sisaan tidak
menyebar normal
(normal <- jarque.bera.test(res.fem.within))##
## Jarque Bera Test
##
## data: res.fem.within
## X-squared = 0.62825, df = 2, p-value = 0.7304Uji Homoskedastisitas
Hipotesis:
H0 = Sisaan memiliki ragam homogen
H1 = Sisaan
tidak memiliki ragam homogen
(homos <- bptest(fem.within))##
## studentized Breusch-Pagan test
##
## data: fem.within
## BP = 1.0464, df = 4, p-value = 0.9027Uji Autokorelasi
Hipotesis:
H0 = Sisaan saling bebas
H1 = Sisaan tidak
saling bebas
(autokol <- pbgtest(fem.within)) ##
## Breusch-Godfrey/Wooldridge test for serial correlation in panel models
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 13.534, df = 3, p-value = 0.003613
## alternative hypothesis: serial correlation in idiosyncratic errorsdata.frame(Asumsi=c("Normalitas","Homoskesdastisitas","Non-Autokorelasi"),
PValue=c(normal$p.value,homos$p.value,autokol$p.value),
Keputusan= c("Terima H0","Terima H0","Tolak H0"),
Kesimpulan= c("Sisaan menyebar normal",
"Ragam sisaan homogen",
"Terdapat autokorelasi")) ## Asumsi PValue Keputusan Kesimpulan
## 1 Normalitas 0.730427645 Terima H0 Sisaan menyebar normal
## 2 Homoskesdastisitas 0.902677288 Terima H0 Ragam sisaan homogen
## 3 Non-Autokorelasi 0.003613039 Tolak H0 Terdapat autokorelasiBerdasarkan hasil uji diagnostik sisaan model FEM Within, didapatkan bahwa ketiga asumsi menghasilkan satu nilai p-value < 0.05, sehingga keputusannya adalah tolak H0 pada uji Autokorelasi. Dengan kata lain, asumsi normalitas terpenuhi dan memiliki ragam yang homogen serta sisaan tidak saling bebas.
Fixed Effect Model dengan LSDV
Fixed effect model dengan Least Squares Dummy Variable (LSDV) adalah salah satu teknik estimasi dalam analisis data panel atau cross-section yang memperhitungkan efek tetap dari setiap individu atau unit pengamatan dalam sampel.
Dalam LSDV, variabel dummy dibuat untuk mewakili setiap individu atau unit pengamatan dalam sampel. Kemudian, variabel dummy ini dimasukkan ke dalam model regresi untuk memperkirakan efek tetap dari setiap individu atau unit pengamatan dalam sampel.
Metode ini dapat digunakan dalam analisis data panel atau cross-section untuk memeriksa hubungan antara variabel dependen dan independen dengan memperhitungkan efek tetap dari setiap individu atau unit pengamatan dalam sampel. LSDV juga berguna untuk menangani masalah heteroskedastisitas dan autokorelasi dalam data panel.
Dalam penelitian seseorang, metode ini dapat digunakan untuk memperkirakan efek tetap dari setiap individu atau unit pengamatan dalam sampel terhadap variabel dependen yang ingin diteliti. Metode ini dapat membantu untuk mengontrol perbedaan antar individu dalam karakteristik atau faktor lain yang mungkin memengaruhi hubungan tersebut, sehingga memungkinkan peneliti untuk lebih akurat menganalisis hubungan antara variabel independen dan dependen.
A. Pendugaan Dengan Memasukkan Intersep
fem.lsdv <- lm(TPT ~ RLS + TPAK + PDRB + PPM + factor(Tahun) + factor(`Kab/Kota`), data = data_adp)
summary(fem.lsdv)##
## Call:
## lm(formula = TPT ~ RLS + TPAK + PDRB + PPM + factor(Tahun) +
## factor(`Kab/Kota`), data = data_adp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.33122 -0.34488 -0.01369 0.31636 1.24842
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 32.60813 6.72012 4.852 1.33e-05 ***
## RLS -0.74490 0.56278 -1.324 0.191905
## TPAK -0.12807 0.06675 -1.919 0.060968 .
## PDRB -0.25720 0.10936 -2.352 0.022828 *
## PPM -1.69775 0.43388 -3.913 0.000287 ***
## factor(Tahun)2020 2.17069 0.89552 2.424 0.019170 *
## factor(Tahun)2021 3.83795 0.77657 4.942 9.81e-06 ***
## factor(`Kab/Kota`)Bandung Barat 8.46227 1.72886 4.895 1.15e-05 ***
## factor(`Kab/Kota`)Bekasi -0.74125 1.03374 -0.717 0.476816
## factor(`Kab/Kota`)Bogor 5.26072 0.71684 7.339 2.23e-09 ***
## factor(`Kab/Kota`)Ciamis -1.00718 0.95243 -1.057 0.295586
## factor(`Kab/Kota`)Cianjur 7.76423 1.74104 4.460 4.94e-05 ***
## factor(`Kab/Kota`)Cirebon 9.65905 2.02813 4.763 1.80e-05 ***
## factor(`Kab/Kota`)Garut 5.01310 1.54862 3.237 0.002191 **
## factor(`Kab/Kota`)Indramayu 9.09831 2.56785 3.543 0.000892 ***
## factor(`Kab/Kota`)Karawang 5.20398 0.98198 5.299 2.89e-06 ***
## factor(`Kab/Kota`)Kota Bandung -0.40050 1.49053 -0.269 0.789316
## factor(`Kab/Kota`)Kota Banjar -1.74516 0.67639 -2.580 0.012992 *
## factor(`Kab/Kota`)Kota Bekasi -0.01855 1.48356 -0.013 0.990076
## factor(`Kab/Kota`)Kota Bogor 4.34351 0.94232 4.609 3.01e-05 ***
## factor(`Kab/Kota`)Kota Cimahi 2.21934 1.19428 1.858 0.069264 .
## factor(`Kab/Kota`)Kota Cirebon 7.68268 1.46362 5.249 3.43e-06 ***
## factor(`Kab/Kota`)Kota Depok -4.70140 2.03143 -2.314 0.024975 *
## factor(`Kab/Kota`)Kota Sukabumi 4.23062 0.98373 4.301 8.31e-05 ***
## factor(`Kab/Kota`)Kota Tasikmalaya 10.43467 2.65004 3.938 0.000265 ***
## factor(`Kab/Kota`)Kuningan 12.27097 2.53087 4.849 1.35e-05 ***
## factor(`Kab/Kota`)Majalengka 5.56931 2.03299 2.739 0.008613 **
## factor(`Kab/Kota`)Pangandaran 1.17164 1.31233 0.893 0.376418
## factor(`Kab/Kota`)Purwakarta 4.48535 0.95755 4.684 2.34e-05 ***
## factor(`Kab/Kota`)Subang 4.94774 1.43668 3.444 0.001199 **
## factor(`Kab/Kota`)Sukabumi 0.88806 1.01332 0.876 0.385186
## factor(`Kab/Kota`)Sumedang 7.18394 1.50402 4.776 1.72e-05 ***
## factor(`Kab/Kota`)Tasikmalaya 4.51259 1.64930 2.736 0.008690 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7044 on 48 degrees of freedom
## Multiple R-squared: 0.946, Adjusted R-squared: 0.91
## F-statistic: 26.26 on 32 and 48 DF, p-value: < 2.2e-16B. Pendugaan Tanpa Memasukkan Intersep
fem.lsdv1 <- lm(TPT ~ RLS + TPAK + PDRB + PPM + factor(Tahun) + factor(`Kab/Kota`)-1, data = data_adp)
summary(fem.lsdv1)##
## Call:
## lm(formula = TPT ~ RLS + TPAK + PDRB + PPM + factor(Tahun) +
## factor(`Kab/Kota`) - 1, data = data_adp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.33122 -0.34488 -0.01369 0.31636 1.24842
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## RLS -0.74490 0.56278 -1.324 0.191905
## TPAK -0.12807 0.06675 -1.919 0.060968 .
## PDRB -0.25720 0.10936 -2.352 0.022828 *
## PPM -1.69775 0.43388 -3.913 0.000287 ***
## factor(Tahun)2019 32.60813 6.72012 4.852 1.33e-05 ***
## factor(Tahun)2020 34.77882 6.88387 5.052 6.75e-06 ***
## factor(Tahun)2021 36.44608 7.09657 5.136 5.07e-06 ***
## factor(`Kab/Kota`)Bandung Barat 8.46227 1.72886 4.895 1.15e-05 ***
## factor(`Kab/Kota`)Bekasi -0.74125 1.03374 -0.717 0.476816
## factor(`Kab/Kota`)Bogor 5.26072 0.71684 7.339 2.23e-09 ***
## factor(`Kab/Kota`)Ciamis -1.00718 0.95243 -1.057 0.295586
## factor(`Kab/Kota`)Cianjur 7.76423 1.74104 4.460 4.94e-05 ***
## factor(`Kab/Kota`)Cirebon 9.65905 2.02813 4.763 1.80e-05 ***
## factor(`Kab/Kota`)Garut 5.01310 1.54862 3.237 0.002191 **
## factor(`Kab/Kota`)Indramayu 9.09831 2.56785 3.543 0.000892 ***
## factor(`Kab/Kota`)Karawang 5.20398 0.98198 5.299 2.89e-06 ***
## factor(`Kab/Kota`)Kota Bandung -0.40050 1.49053 -0.269 0.789316
## factor(`Kab/Kota`)Kota Banjar -1.74516 0.67639 -2.580 0.012992 *
## factor(`Kab/Kota`)Kota Bekasi -0.01855 1.48356 -0.013 0.990076
## factor(`Kab/Kota`)Kota Bogor 4.34351 0.94232 4.609 3.01e-05 ***
## factor(`Kab/Kota`)Kota Cimahi 2.21934 1.19428 1.858 0.069264 .
## factor(`Kab/Kota`)Kota Cirebon 7.68268 1.46362 5.249 3.43e-06 ***
## factor(`Kab/Kota`)Kota Depok -4.70140 2.03143 -2.314 0.024975 *
## factor(`Kab/Kota`)Kota Sukabumi 4.23062 0.98373 4.301 8.31e-05 ***
## factor(`Kab/Kota`)Kota Tasikmalaya 10.43467 2.65004 3.938 0.000265 ***
## factor(`Kab/Kota`)Kuningan 12.27097 2.53087 4.849 1.35e-05 ***
## factor(`Kab/Kota`)Majalengka 5.56931 2.03299 2.739 0.008613 **
## factor(`Kab/Kota`)Pangandaran 1.17164 1.31233 0.893 0.376418
## factor(`Kab/Kota`)Purwakarta 4.48535 0.95755 4.684 2.34e-05 ***
## factor(`Kab/Kota`)Subang 4.94774 1.43668 3.444 0.001199 **
## factor(`Kab/Kota`)Sukabumi 0.88806 1.01332 0.876 0.385186
## factor(`Kab/Kota`)Sumedang 7.18394 1.50402 4.776 1.72e-05 ***
## factor(`Kab/Kota`)Tasikmalaya 4.51259 1.64930 2.736 0.008690 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7044 on 48 degrees of freedom
## Multiple R-squared: 0.9967, Adjusted R-squared: 0.9944
## F-statistic: 433.1 on 33 and 48 DF, p-value: < 2.2e-16Dengan membandingkan model dengan intersep dan tanpa intersep, dapat diketahui keduanya memiliki kesamaan dalam jumlah unit individu kabupaten/kota yang signifikan. Dari kabupaten/kota yang ada, terdapat beberapa kabupaten/kota seperti Bekasi, Ciamis, Kota Bandung, Kota Bekasi, Kota Cimahi, Pangandaran, dan Sukabumi yang tidak signifikan terhadap model.
Uji Chow
H0: Model Common Effect.
H1: Model Fixed Effect.
H0 ditolak
jika P-Value lebih kecil dari nilai alpha. Nilai
alpha yang digunakan sebesar 5%.
Memilih Antara CEM dan FEM
pooltest(cem,fem.within)##
## F statistic
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## F = 7.9419, df1 = 26, df2 = 50, p-value = 2.728e-10
## alternative hypothesis: unstabilityKeputusan : Karena p-value
(\2.728e-10) < alpha
(0.05), maka Tolak H0.
Kesimpulan : Dengan tingkat keyakinan
95%, kita yakin bahwa metode fixed effect
lebih baik digunakan daripada menggunakan metode
common effect.
Atau dapat juga dengan menggunakan
alternatif dibawah ini, apabila ingin membandingkan model OLS dengan
FEM.
OLS
ols <- lm(TPT ~ RLS + TPAK + PDRB + PPM, data=datapanel)
summary(ols)##
## Call:
## lm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.0842 -1.3777 0.0880 0.9973 3.4828
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 32.94244 4.63713 7.104 5.59e-10 ***
## RLS 0.21154 0.19954 1.060 0.292
## TPAK -0.39217 0.05303 -7.396 1.57e-10 ***
## PDRB -0.28539 0.05782 -4.936 4.60e-06 ***
## PPM 0.06124 0.08713 0.703 0.484
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.587 on 76 degrees of freedom
## Multiple R-squared: 0.5655, Adjusted R-squared: 0.5427
## F-statistic: 24.73 on 4 and 76 DF, p-value: 3.933e-13pFtest(fem.within, ols)##
## F test for individual effects
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## F = 7.9419, df1 = 26, df2 = 50, p-value = 2.728e-10
## alternative hypothesis: significant effectsP-Value yang kurang dari α (5%) berimplikasi Tolak H0 yang menandakan bahwa terdapat panel efek pada model. Berarti model FEM lebih baik digunakan dalam pemodelan dibandingkan hanya dengan OLS.
Pengaruh Spesifik Individu dan/atau Waktu
Berdasarkan hasil Uji Chouw diperoleh hasil bahwa model yang terpilih adalah model FEM. Sehingga lebih lanjut kita perlu mengetahui apakah terdapat pengaruh individu dan/atau pengaruh waktu.
Efek Individu dan Waktu
plmtest(fem.within,type = "bp", effect="twoways" )##
## Lagrange Multiplier Test - two-ways effects (Breusch-Pagan)
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 32.319, df = 2, p-value = 9.593e-08
## alternative hypothesis: significant effectsEfek Individu
plmtest(fem.within,type = "bp", effect="individual" )##
## Lagrange Multiplier Test - (Breusch-Pagan)
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 32.044, df = 1, p-value = 1.507e-08
## alternative hypothesis: significant effectsEfek Waktu
plmtest(fem.within,type = "bp", effect="time" )##
## Lagrange Multiplier Test - time effects (Breusch-Pagan)
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 0.27497, df = 1, p-value = 0.6
## alternative hypothesis: significant effectsBerdasarkan hasil di atas, diperoleh bahwa model memiliki pengaruh individu dan waktu sehingga model yang tepat digunakan adalah FEM Two ways
fem.two <- plm(TPT ~ RLS + TPAK + PDRB + PPM, datapanel, model = "within",
effect= "twoways", index = c("Kab.Kota","Tahun"))
summary(fem.two)## Twoways effects Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## effect = "twoways", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.331225 -0.344876 -0.013689 0.316364 1.248422
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS -0.744896 0.562779 -1.3236 0.1919050
## TPAK -0.128070 0.066745 -1.9188 0.0609681 .
## PDRB -0.257199 0.109358 -2.3519 0.0228283 *
## PPM -1.697752 0.433883 -3.9129 0.0002868 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 39.101
## Residual Sum of Squares: 23.813
## R-Squared: 0.39098
## Adj. R-Squared: -0.01503
## F-statistic: 7.70386 on 4 and 48 DF, p-value: 7.0383e-05Pengaruh Twoways
summary(fixef(fem.two, effect="twoways"))## Estimate
## Bandung-2019 32.60813
## Bandung-2020 34.77882
## Bandung-2021 36.44608
## Bandung Barat-2019 41.07040
## Bandung Barat-2020 43.24109
## Bandung Barat-2021 44.90835
## Bekasi-2019 31.86688
## Bekasi-2020 34.03758
## Bekasi-2021 35.70483
## Bogor-2019 37.86885
## Bogor-2020 40.03955
## Bogor-2021 41.70680
## Ciamis-2019 31.60095
## Ciamis-2020 33.77165
## Ciamis-2021 35.43890
## Cianjur-2019 40.37236
## Cianjur-2020 42.54305
## Cianjur-2021 44.21031
## Cirebon-2019 42.26718
## Cirebon-2020 44.43788
## Cirebon-2021 46.10513
## Garut-2019 37.62123
## Garut-2020 39.79192
## Garut-2021 41.45918
## Indramayu-2019 41.70644
## Indramayu-2020 43.87713
## Indramayu-2021 45.54439
## Karawang-2019 37.81211
## Karawang-2020 39.98280
## Karawang-2021 41.65005
## Kota Bandung-2019 32.20763
## Kota Bandung-2020 34.37833
## Kota Bandung-2021 36.04558
## Kota Banjar-2019 30.86297
## Kota Banjar-2020 33.03366
## Kota Banjar-2021 34.70091
## Kota Bekasi-2019 32.58958
## Kota Bekasi-2020 34.76027
## Kota Bekasi-2021 36.42753
## Kota Bogor-2019 36.95164
## Kota Bogor-2020 39.12233
## Kota Bogor-2021 40.78959
## Kota Cimahi-2019 34.82747
## Kota Cimahi-2020 36.99816
## Kota Cimahi-2021 38.66541
## Kota Cirebon-2019 40.29081
## Kota Cirebon-2020 42.46150
## Kota Cirebon-2021 44.12876
## Kota Depok-2019 27.90673
## Kota Depok-2020 30.07742
## Kota Depok-2021 31.74468
## Kota Sukabumi-2019 36.83875
## Kota Sukabumi-2020 39.00945
## Kota Sukabumi-2021 40.67670
## Kota Tasikmalaya-2019 43.04280
## Kota Tasikmalaya-2020 45.21349
## Kota Tasikmalaya-2021 46.88075
## Kuningan-2019 44.87910
## Kuningan-2020 47.04980
## Kuningan-2021 48.71705
## Majalengka-2019 38.17744
## Majalengka-2020 40.34813
## Majalengka-2021 42.01539
## Pangandaran-2019 33.77977
## Pangandaran-2020 35.95047
## Pangandaran-2021 37.61772
## Purwakarta-2019 37.09348
## Purwakarta-2020 39.26417
## Purwakarta-2021 40.93142
## Subang-2019 37.55587
## Subang-2020 39.72656
## Subang-2021 41.39382
## Sukabumi-2019 33.49619
## Sukabumi-2020 35.66688
## Sukabumi-2021 37.33413
## Sumedang-2019 39.79207
## Sumedang-2020 41.96276
## Sumedang-2021 43.63002
## Tasikmalaya-2019 37.12072
## Tasikmalaya-2020 39.29142
## Tasikmalaya-2021 40.95867res.fem.two <- residuals(fem.two)Normality
library(tseries)
jarque.bera.test(res.fem.two)##
## Jarque Bera Test
##
## data: res.fem.two
## X-squared = 0.15424, df = 2, p-value = 0.9258Histogram
hist(res.fem.two,
xlab = "Sisaan",
col = "#27D3D3",
breaks=30,
prob = TRUE)
lines(density(res.fem.two), # density plot
lwd = 2, # thickness of line
col = "chocolate3")
Plot QQ-Norm
set.seed(1234)
res.fem.two1 <- as.numeric(res.fem.two)
qqnorm(res.fem.two1,datax=T, col="blue")
qqline(rnorm(length(res.fem.two1),mean(res.fem.two1),sd(res.fem.two1)),datax=T, col="red")
Auto-Corelation
pbgtest(fem.two) #alternatif : Terdapat autokorelasi##
## Breusch-Godfrey/Wooldridge test for serial correlation in panel models
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 17.606, df = 3, p-value = 0.0005303
## alternative hypothesis: serial correlation in idiosyncratic errorsHomogen
bptest(fem.two) #alternatif : Ragam Tidak Homogen##
## studentized Breusch-Pagan test
##
## data: fem.two
## BP = 1.0464, df = 4, p-value = 0.9027Kesimpulan Pada Model FEM individual effect :
1. Kenormalan
sisaan Terpenuhi.
2. Kebebasan sisaan Belum Terpenuhi.
3.
Kehomogenan ragam Terpenuhi.
Uji Kebebasan pada Unit Cross Section
Biasanya pada data panel dengan series panjang, sering kali ditemukan adanya dependensi pada unit individu tetapi pada kasus data panel dengan series yang tidak terlalu panjang hal ini jarang terjadi. Pada bagian ini akan dicek apakah terdapat dependensi pada unit cross section di data panel yang digunakan.
Hipotesis untuk menguji ada tidaknya dependensi pada unit
individu yaitu :
H0 = tidak terdapat dependensi antar unit
individu.
H1 = terdapat dependensi antar unit individu.
Tolak H0 jika nilai P-Value yang
diperoleh < α (5%) yang berarti terdapat
dependensi pada unit individu.
pcdtest(fem.two, test = c("lm"), index=NULL,w =NULL )##
## Breusch-Pagan LM test for cross-sectional dependence in panels
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 510.09, df = 351, p-value = 5.636e-08
## alternative hypothesis: cross-sectional dependenceKarena nilai P-Value yang diperoleh < α (5%) maka tolak H0 yang berarti terdapat dependensi pada unit individu.
FEM INDIVIDUAL
Efek komponen sisaaan satu arah pada individu.
fem2 <- plm(TPT ~ RLS + TPAK + PDRB + PPM, datapanel, model = "within",effect= "individual",index = c("Kab.Kota","Tahun"))
summary(fem2)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## effect = "individual", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.871451 -0.503443 0.041329 0.421842 1.398651
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS 0.446339 0.618656 0.7215 0.473983
## TPAK -0.212108 0.079037 -2.6837 0.009849 **
## PDRB -0.219655 0.035480 -6.1910 1.104e-07 ***
## PPM 0.429229 0.177325 2.4206 0.019170 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 104.68
## Residual Sum of Squares: 37.331
## R-Squared: 0.64337
## Adj. R-Squared: 0.42939
## F-statistic: 22.55 on 4 and 50 DF, p-value: 1.0918e-10Pendugaan dengan model FEM efek komponen sisaan satu arah terhadap individual memberikan nilai R2-Adj sebesar 42,94%. Nilai dugaan koefisien capital masing-masing sebesar sebesar 0.446339 RLS, -0.212108 TPAK, -0.219655 PDRB, 0.429229 PPM. Nilai koefisien capital menunjukkan besarnya perubahan tingkat pengangguran terbuka pada rata-rata setiap kabupaten/kota, ketika (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) meningkat sebesar satu satuan unit. Nilai dari p-value < α (5%) menunjukkan bahwa (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) memiliki pengaruh yang signifikan terhadap tingkat pengangguran terbuka pada taraf 5%.
Pengaruh Individu
summary(fixef(fem2, effect="individual"))## Estimate Std. Error t-value Pr(>|t|)
## Bandung 14.6863 7.0000 2.0980 0.040976 *
## Bandung Barat 15.7769 6.4543 2.4444 0.018078 *
## Bekasi 17.9965 7.1973 2.5004 0.015727 *
## Bogor 18.8084 6.7071 2.8042 0.007163 **
## Ciamis 13.8859 6.8291 2.0333 0.047344 *
## Cianjur 17.3245 6.4491 2.6863 0.009780 **
## Cirebon 16.6153 6.3525 2.6156 0.011746 *
## Garut 14.2507 6.3051 2.2602 0.028196 *
## Indramayu 15.0302 6.3823 2.3550 0.022490 *
## Karawang 17.9015 6.5362 2.7388 0.008522 **
## Kota Bandung 18.1304 7.9644 2.2764 0.027136 *
## Kota Banjar 14.4303 6.9886 2.0648 0.044146 *
## Kota Bekasi 17.2682 8.0919 2.1340 0.037771 *
## Kota Bogor 17.5671 7.4513 2.3576 0.022349 *
## Kota Cimahi 18.5045 7.7076 2.4008 0.020120 *
## Kota Cirebon 15.6903 7.1960 2.1804 0.033960 *
## Kota Depok 16.7068 8.1024 2.0620 0.044429 *
## Kota Sukabumi 15.7418 6.9873 2.2529 0.028682 *
## Kota Tasikmalaya 12.3143 7.0826 1.7387 0.088246 .
## Kuningan 15.8751 6.3418 2.5032 0.015618 *
## Majalengka 12.4603 6.6061 1.8862 0.065086 .
## Pangandaran 13.5846 7.1766 1.8929 0.064165 .
## Purwakarta 16.7896 6.4850 2.5890 0.012574 *
## Subang 16.7337 6.4533 2.5931 0.012443 *
## Sukabumi 16.5914 6.3177 2.6262 0.011430 *
## Sumedang 15.5443 6.8883 2.2566 0.028433 *
## Tasikmalaya 13.6407 6.6103 2.0636 0.044269 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1FEM TIME
Pada fixed effect model dengan efek komponen sisaan satu
arah karena adanya pengaruh waktu, model panel yang dibangun
dengan memerhatikan spesifik waktu. Berikut merupakan contoh untuk
membangun model panel dengan efek komponen sisaan satu arah
pengaruh waktu. Pada fungsi plm dengan memberikan
keterangan effect = “time”
fem.time <- plm(TPT ~ RLS + TPAK + PDRB + PPM, datapanel, model = "within",effect= "time",index = c("Kab.Kota","Tahun"))
summary(fem.time)## Oneway (time) effect Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## effect = "time", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -3.67663 -1.09524 0.25236 0.96285 3.43890
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS 0.1284512 0.2068765 0.6209 0.536568
## TPAK -0.3881434 0.0520221 -7.4611 1.346e-10 ***
## PDRB -0.5121223 0.1706866 -3.0004 0.003673 **
## PPM 0.0069581 0.0914465 0.0761 0.939553
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 375.19
## Residual Sum of Squares: 176.07
## R-Squared: 0.53072
## Adj. R-Squared: 0.49267
## F-statistic: 20.9222 on 4 and 74 DF, p-value: 1.4374e-11Pendugaan dengan model FEM efek komponen sisaan satu arah terhadap waktu memberikan nilai R2-Adj sebesar 49,27%. Nilai dugaan koefisien capital masing-masing sebesar sebesar 0.1284512 RLS, -0.3881434 TPAK, -0.5121223 PDRB, 0.0069581 PPM. Nilai koefisien capital menunjukkan besarnya perubahan tingkat pengangguran terbuka pada rata-rata setiap kabupaten/kota, ketika (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) meningkat sebesar satu satuan unit. Nilai dari p-value < α (5%) menunjukkan bahwa (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) memiliki pengaruh yang signifikan terhadap jumlah kemiskinan pada taraf 5%. Sedangkan untuk melihat pengaruh spesifik waktu :
Pengaruh Waktu
summary(fixef(fem.time, effect="time"))## Estimate Std. Error t-value Pr(>|t|)
## 2019 34.8429 4.5797 7.6081 7.115e-11 ***
## 2020 33.3195 4.7398 7.0298 8.651e-10 ***
## 2021 35.2081 4.6829 7.5184 1.050e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1FEM Two Way
Pada model dengan efek komponen sisaan dua arah, model panel
yang dibangun dengan memperhatikan spesifik individu dan juga waktu.
Berikut merupakan contoh untuk membangun model panel dengan
efek komponen sisaan dua arah. Pada fungsi plm dengan
memberikan keterangan effect = “twoways”
fem.twoway <- plm(TPT ~ RLS + TPAK + PDRB + PPM, datapanel, model = "within",
effect= "twoways", index = c("Kab.Kota","Tahun"))
summary(fem.twoway)## Twoways effects Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## effect = "twoways", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.331225 -0.344876 -0.013689 0.316364 1.248422
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS -0.744896 0.562779 -1.3236 0.1919050
## TPAK -0.128070 0.066745 -1.9188 0.0609681 .
## PDRB -0.257199 0.109358 -2.3519 0.0228283 *
## PPM -1.697752 0.433883 -3.9129 0.0002868 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 39.101
## Residual Sum of Squares: 23.813
## R-Squared: 0.39098
## Adj. R-Squared: -0.01503
## F-statistic: 7.70386 on 4 and 48 DF, p-value: 7.0383e-05Berdasarkan output di atas, model FEM Two Way menghasilkan peubah RLS dan TPAK yang tidak signifikan terhadap model. Selain itu, didapatkan juga koefisien determinasi dari model sebesar -01,50%. Hal tersebut menandakan bahwa keragaman yang dapat dijelaskan oleh model tergolong sangat kecil sehingga model FEM Two Way belum baik digunakan untuk memodelkan kemiskinan di Indonesia.
Sedangkan untuk mengetahui nilai pengaruh spesifik individu dan waktu dilakukan dengan langkah berikut :
Pengaruh Waktu
summary(fixef(fem.twoway, effect="time"))## Estimate
## 2019 32.60813
## 2020 34.77882
## 2021 36.44608Pengaruh Individu
summary(fixef(fem.twoway, effect="individual"))## Estimate
## Bandung 32.60813
## Bandung Barat 41.07040
## Bekasi 31.86688
## Bogor 37.86885
## Ciamis 31.60095
## Cianjur 40.37236
## Cirebon 42.26718
## Garut 37.62123
## Indramayu 41.70644
## Karawang 37.81211
## Kota Bandung 32.20763
## Kota Banjar 30.86297
## Kota Bekasi 32.58958
## Kota Bogor 36.95164
## Kota Cimahi 34.82747
## Kota Cirebon 40.29081
## Kota Depok 27.90673
## Kota Sukabumi 36.83875
## Kota Tasikmalaya 43.04280
## Kuningan 44.87910
## Majalengka 38.17744
## Pangandaran 33.77977
## Purwakarta 37.09348
## Subang 37.55587
## Sukabumi 33.49619
## Sumedang 39.79207
## Tasikmalaya 37.12072Pengaruh Waktu dan Individu
summary(fixef(fem.twoway, effect="twoways"))## Estimate
## Bandung-2019 32.60813
## Bandung-2020 34.77882
## Bandung-2021 36.44608
## Bandung Barat-2019 41.07040
## Bandung Barat-2020 43.24109
## Bandung Barat-2021 44.90835
## Bekasi-2019 31.86688
## Bekasi-2020 34.03758
## Bekasi-2021 35.70483
## Bogor-2019 37.86885
## Bogor-2020 40.03955
## Bogor-2021 41.70680
## Ciamis-2019 31.60095
## Ciamis-2020 33.77165
## Ciamis-2021 35.43890
## Cianjur-2019 40.37236
## Cianjur-2020 42.54305
## Cianjur-2021 44.21031
## Cirebon-2019 42.26718
## Cirebon-2020 44.43788
## Cirebon-2021 46.10513
## Garut-2019 37.62123
## Garut-2020 39.79192
## Garut-2021 41.45918
## Indramayu-2019 41.70644
## Indramayu-2020 43.87713
## Indramayu-2021 45.54439
## Karawang-2019 37.81211
## Karawang-2020 39.98280
## Karawang-2021 41.65005
## Kota Bandung-2019 32.20763
## Kota Bandung-2020 34.37833
## Kota Bandung-2021 36.04558
## Kota Banjar-2019 30.86297
## Kota Banjar-2020 33.03366
## Kota Banjar-2021 34.70091
## Kota Bekasi-2019 32.58958
## Kota Bekasi-2020 34.76027
## Kota Bekasi-2021 36.42753
## Kota Bogor-2019 36.95164
## Kota Bogor-2020 39.12233
## Kota Bogor-2021 40.78959
## Kota Cimahi-2019 34.82747
## Kota Cimahi-2020 36.99816
## Kota Cimahi-2021 38.66541
## Kota Cirebon-2019 40.29081
## Kota Cirebon-2020 42.46150
## Kota Cirebon-2021 44.12876
## Kota Depok-2019 27.90673
## Kota Depok-2020 30.07742
## Kota Depok-2021 31.74468
## Kota Sukabumi-2019 36.83875
## Kota Sukabumi-2020 39.00945
## Kota Sukabumi-2021 40.67670
## Kota Tasikmalaya-2019 43.04280
## Kota Tasikmalaya-2020 45.21349
## Kota Tasikmalaya-2021 46.88075
## Kuningan-2019 44.87910
## Kuningan-2020 47.04980
## Kuningan-2021 48.71705
## Majalengka-2019 38.17744
## Majalengka-2020 40.34813
## Majalengka-2021 42.01539
## Pangandaran-2019 33.77977
## Pangandaran-2020 35.95047
## Pangandaran-2021 37.61772
## Purwakarta-2019 37.09348
## Purwakarta-2020 39.26417
## Purwakarta-2021 40.93142
## Subang-2019 37.55587
## Subang-2020 39.72656
## Subang-2021 41.39382
## Sukabumi-2019 33.49619
## Sukabumi-2020 35.66688
## Sukabumi-2021 37.33413
## Sumedang-2019 39.79207
## Sumedang-2020 41.96276
## Sumedang-2021 43.63002
## Tasikmalaya-2019 37.12072
## Tasikmalaya-2020 39.29142
## Tasikmalaya-2021 40.95867Pendugaan dengan Robustness
Pada data panel yang memiliki keberadaan heteroscedasticity,
dan autokorelasi maka untuk mengontrol hal tersebut, pendugaan model
dapat dilakukan dengan menggunakan robust covariance matrix
yaitu pada fungsi vcovHC. Identifikasi mengenai
heteroscedasticity dapat menggunakan uji Breucsh
Pagan. Terdapat beberapa kondisi diantaranya :
white1 : kondisi heteroskedasticity umum
tetapi tidak terdapat serial correlation.( Recommended
untuk random effects) white2 : kondisi
“white1” yang terbatas pada keragaman antar grup.(
Recommended untuk random effects)
arellano : Terdapat heteroskedasticity dan
serial correlation. ( Recommended untuk fixed
effects) Kondisi diatas dapat diaplikasikan dengan fungsi berikut :
HC0 -> heteroskedasticity consistent(
default ) HC1,HC2, HC3 –> Biasa
direkomendasikan untuk observasi yang sedikit/kecil.
HC3 memberikan bobot yang lebih sedikit apabila
terdapat amatan berpengaruh. HC4 -> Pada jumlah
observasi yang sedikit dengan adanya amatan berpengaruh.
Model Awal
coeftest(fem2)##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## RLS 0.446339 0.618656 0.7215 0.473983
## TPAK -0.212108 0.079037 -2.6837 0.009849 **
## PDRB -0.219655 0.035480 -6.1910 1.104e-07 ***
## PPM 0.429229 0.177325 2.4206 0.019170 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1coeftest(fem2, vcovHC)##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## RLS 0.446339 0.754396 0.5917 0.55675
## TPAK -0.212108 0.069928 -3.0332 0.00383 **
## PDRB -0.219655 0.028724 -7.6471 5.865e-10 ***
## PPM 0.429229 0.225743 1.9014 0.06302 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1coeftest(fem2, vcovHC(fem2, method = "arellano"))##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## RLS 0.446339 0.754396 0.5917 0.55675
## TPAK -0.212108 0.069928 -3.0332 0.00383 **
## PDRB -0.219655 0.028724 -7.6471 5.865e-10 ***
## PPM 0.429229 0.225743 1.9014 0.06302 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1coeftest(fem2, vcovHC(fem2, type = "HC3"))##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## RLS 0.446339 0.805184 0.5543 0.581823
## TPAK -0.212108 0.074530 -2.8459 0.006404 **
## PDRB -0.219655 0.030714 -7.1515 3.477e-09 ***
## PPM 0.429229 0.241527 1.7771 0.081627 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1t(sapply(c("HC0", "HC1", "HC2", "HC3", "HC4"), function(x) sqrt(diag(vcovHC(fem2, type = x)))))## RLS TPAK PDRB PPM
## HC0 0.7543956 0.06992755 0.02872386 0.2257432
## HC1 0.7737422 0.07172086 0.02946049 0.2315324
## HC2 0.7793120 0.07218116 0.02969840 0.2334742
## HC3 0.8051836 0.07452979 0.03071431 0.2415268
## HC4 0.7894930 0.07331606 0.03017538 0.2373523FEM FIRST DIFFERENCING
library(plm)
diff = function(x) {x-dplyr::lag(x)}
datapanel %<>%
group_by(Tahun) %>%
mutate(dTPT=diff(`TPT`),
dRLS=diff(`RLS`),
dTPAK=diff(`TPAK`),
dPDRB=diff(`PDRB`),
dPPM=diff(`PPM`)) %>%
ungroup()
fem_fd <- plm(dTPT ~ dRLS + dTPAK + dPDRB + dPPM , data = datapanel, effect = "individual",model = "within",index = c("Kab.Kota","Tahun"))
summary(fem_fd)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = dTPT ~ dRLS + dTPAK + dPDRB + dPPM, data = datapanel,
## effect = "individual", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 26, T = 3, N = 78
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.609423 -0.446059 -0.018828 0.456525 1.832886
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## dRLS -1.455234 0.541058 -2.6896 0.0098098 **
## dTPAK -0.110600 0.071551 -1.5457 0.1287355
## dPDRB -0.240271 0.098011 -2.4515 0.0179187 *
## dPPM -1.454074 0.407505 -3.5682 0.0008275 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 61.845
## Residual Sum of Squares: 38.99
## R-Squared: 0.36956
## Adj. R-Squared: -0.011339
## F-statistic: 7.03418 on 4 and 48 DF, p-value: 0.0001534Pendugaan dengan model FEM First Differencing memberikan nilai R2-Adj sebesar -01,13%. Nilai dugaan koefisien capital masing-masing sebesar -1.455234 dRLS, -0.110600 dTPAK, -0.240271 dPDRB, -1.454074 dPPM. Nilai koefisien capital (differencing) menunjukkan besarnya perubahan tingkat pengangguran terbuka pada rata-rata setiap kabupaten/kota, ketika (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) meningkat sebesar satu satuan unit. Nilai dari p-value < α (5%) menunjukkan bahwa (rata-rata lama sekolah, tingkat pengangguran terbuka, produk domestik regional bruto, dan persentase penduduk miskin) memiliki pengaruh yang signifikan terhadap jumlah kemiskinan pada taraf 5%. Koefisien capital (differencing) berbeda dibandingkan dengan pendekatan LSDV dan within. Ini karena koefisien dan standard error dari model yang dibedakan pertama hanya identik dengan hasil yang diperoleh sebelumnya ketika ada dua periode waktu. Untuk deret waktu yang lebih lama, koefisien dan standard error akan berbeda.
Melihat Fixed Effect dari FEM , nilai intersept yang berbeda beda pada masing-masing individu
fixef(fem_fd, type="level")## Bandung Barat Bekasi Bogor Ciamis
## 7.14591 -7.23137 4.75961 -6.76970
## Cianjur Cirebon Garut Indramayu
## 7.64061 1.96588 -4.14911 2.69911
## Karawang Kota Bandung Kota Banjar Kota Bekasi
## -1.96482 -2.59236 -3.42981 3.97475
## Kota Bogor Kota Cimahi Kota Cirebon Kota Depok
## 3.23530 -1.49766 3.82808 -9.84046
## Kota Sukabumi Kota Tasikmalaya Kuningan Majalengka
## 6.84330 4.59451 0.84376 -6.73660
## Pangandaran Purwakarta Subang Sukabumi
## -3.64115 3.87889 -0.50326 -3.32347
## Sumedang Tasikmalaya
## 6.26779 -3.38947FEM
fem <- plm(TPT ~ RLS + TPAK + PDRB + PPM, datapanel, model = "within",effect= "individual",index = c("Kab.Kota","Tahun"))
summary(fem)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
## effect = "individual", model = "within", index = c("Kab.Kota",
## "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.871451 -0.503443 0.041329 0.421842 1.398651
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## RLS 0.446339 0.618656 0.7215 0.473983
## TPAK -0.212108 0.079037 -2.6837 0.009849 **
## PDRB -0.219655 0.035480 -6.1910 1.104e-07 ***
## PPM 0.429229 0.177325 2.4206 0.019170 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 104.68
## Residual Sum of Squares: 37.331
## R-Squared: 0.64337
## Adj. R-Squared: 0.42939
## F-statistic: 22.55 on 4 and 50 DF, p-value: 1.0918e-10Pendugaan dengan model FEM memberikan nilai R^2-Adj sebesar 0.42939. Nilai dugaan koefisien RLS Sebesar 0.446339, TPAK sebesar -0.212108, PDRB sebesar -0.219655, dan PPM sebesar 0.429229. Nilai dari p-value < α (5%) menunjukkan bahwa x memiliki pengaruh yang signifikan terhadap Tingkat Pengangguran Terbuka pada taraf 5%.
Melihat Fixed Effect dari FEM , Nilai Intersep yang Berbeda-beda pada Masing-masing Individu
fixef(fem, type="level")## Bandung Bandung Barat Bekasi Bogor
## 14.686 15.777 17.997 18.808
## Ciamis Cianjur Cirebon Garut
## 13.886 17.325 16.615 14.251
## Indramayu Karawang Kota Bandung Kota Banjar
## 15.030 17.901 18.130 14.430
## Kota Bekasi Kota Bogor Kota Cimahi Kota Cirebon
## 17.268 17.567 18.505 15.690
## Kota Depok Kota Sukabumi Kota Tasikmalaya Kuningan
## 16.707 15.742 12.314 15.875
## Majalengka Pangandaran Purwakarta Subang
## 12.460 13.585 16.790 16.734
## Sukabumi Sumedang Tasikmalaya
## 16.591 15.544 13.641GLS: One Way Random Effect Model
Jika pada model FEM mengasumsikan bahwa nilai dari pengaruh spesifik individu bersifat fixed pada setiap unit, maka pada model random effect diasumsikan bahwa pengaruh tersebut bersifat acak. Dalam hal ini misalnya model one way individual pada data Jumlah Penduduk Miskin, maka seluruh Kab.Kota memiliki nilai rataan umum intersep yang sama, sedangkan perbedaan intersep antar masing-masing unit Kab.Kota direfleksikan pada error term-nya.
One Way Individual
model <- TPT ~ RLS + TPAK + PDRB + PPM
rem_gls_ind <- plm(model, data = datapanel,
index = c("Kab.Kota","Tahun"),
effect = "individual",
model = "random")
summary(rem_gls_ind)## Oneway (individual) effect Random Effect Model
## (Swamy-Arora's transformation)
##
## Call:
## plm(formula = model, data = datapanel, effect = "individual",
## model = "random", index = c("Kab.Kota", "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Effects:
## var std.dev share
## idiosyncratic 0.7466 0.8641 0.302
## individual 1.7226 1.3125 0.698
## theta: 0.6447
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.80819 -0.56784 -0.12780 0.57290 2.35318
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) 22.253938 4.988487 4.4611 8.156e-06 ***
## RLS 0.569395 0.239087 2.3815 0.01724 *
## TPAK -0.298693 0.059810 -4.9941 5.913e-07 ***
## PDRB -0.240203 0.034918 -6.8791 6.024e-12 ***
## PPM 0.218375 0.099103 2.2035 0.02756 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 147.1
## Residual Sum of Squares: 60.684
## R-Squared: 0.58747
## Adj. R-Squared: 0.56576
## Chisq: 108.23 on 4 DF, p-value: < 2.22e-16Selanjutnya output dirapikan dengan menggunakan library
stargazer, kableExtra, dan
broom.
library(stargazer)##
## Please cite as:## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.## R package version 5.2.3. https://CRAN.R-project.org/package=stargazerstargazer(rem_gls_ind, type='text')##
## ========================================
## Dependent variable:
## ---------------------------
## TPT
## ----------------------------------------
## RLS 0.569**
## (0.239)
##
## TPAK -0.299***
## (0.060)
##
## PDRB -0.240***
## (0.035)
##
## PPM 0.218**
## (0.099)
##
## Constant 22.254***
## (4.988)
##
## ----------------------------------------
## Observations 81
## R2 0.587
## Adjusted R2 0.566
## F Statistic 108.230***
## ========================================
## Note: *p<0.1; **p<0.05; ***p<0.01library(kableExtra)
library(broom)
kable(tidy(rem_gls_ind), digits=3,
caption="GLS_Random Effect Model_Individu")| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| (Intercept) | 22.254 | 4.988 | 4.461 | 0.000 |
| RLS | 0.569 | 0.239 | 2.382 | 0.017 |
| TPAK | -0.299 | 0.060 | -4.994 | 0.000 |
| PDRB | -0.240 | 0.035 | -6.879 | 0.000 |
| PPM | 0.218 | 0.099 | 2.204 | 0.028 |
Berdasarkan output di atas, model GLS REM Individu menghasilkan semua peubah signifikan terhadap model. Selain itu, didapatkan juga koefisien determinasi adjusted dari model sebesar 56,6%. Hal tersebut menandakan bahwa keragaman yang dapat dijelaskan oleh model cenderung besar sehingga model GLS REM Individu cukup baik digunakan untuk memodelkan tingkat penganguuran terbuka di Indonesia.
Pengecekan Pengaruh Spesifik Model
Selanjutnya, dilakukan pemeriksaan mengenai model random effect yang dibangun, apakah terdapat pengaruh efek individu, waktu, atau keduanya.
Efek Individu dan Waktu
plmtest(rem_gls_ind, type = "bp", effect = "twoways" )##
## Lagrange Multiplier Test - two-ways effects (Breusch-Pagan)
##
## data: model
## chisq = 32.319, df = 2, p-value = 9.593e-08
## alternative hypothesis: significant effectsEfek Individu
plmtest(rem_gls_ind,type = "bp", effect = "individu" )##
## Lagrange Multiplier Test - (Breusch-Pagan)
##
## data: model
## chisq = 32.044, df = 1, p-value = 1.507e-08
## alternative hypothesis: significant effectsEfek Waktu
plmtest(rem_gls_ind,type = "bp", effect = "time" )##
## Lagrange Multiplier Test - time effects (Breusch-Pagan)
##
## data: model
## chisq = 0.27497, df = 1, p-value = 0.6
## alternative hypothesis: significant effectsBerdasarkan pengujian pengaruh, didapatkan kesimpulan bahwa H0 ditolak pada model dengan pengaruh two-way, pengaruh individu, ataupun pengaruh waktu. Oleh karena itu, model random effect yang dibangun memiliki pengaruh two-way.
Walaupun memiliki pengaruh two-way, akan dilihat pengaruh baik individu maupun waktu untuk membuktikan two-way atau tidak. Adapun pengaruh individu yang terdapat pada model one way individual adalah sebagai berikut :
tidy_ranef_ind <- tidy(ranef(rem_gls_ind, effect="individual"))## Warning: 'tidy.numeric' is deprecated.
## See help("Deprecated")kable(tidy_ranef_ind, digits=3,
caption = "Pengaruh Acak Individu",
col.names = c("Kab.Kota", "Pengaruh Acak Individu"))| Kab.Kota | Pengaruh Acak Individu |
|---|---|
| Bandung | -1.474 |
| Bandung Barat | -0.027 |
| Bekasi | 1.042 |
| Bogor | 2.287 |
| Ciamis | -1.505 |
| Cianjur | 1.951 |
| Cirebon | 1.210 |
| Garut | -1.287 |
| Indramayu | 0.435 |
| Karawang | 1.720 |
| Kota Bandung | 0.922 |
| Kota Banjar | -1.572 |
| Kota Bekasi | 0.173 |
| Kota Bogor | 0.738 |
| Kota Cimahi | 1.261 |
| Kota Cirebon | -0.274 |
| Kota Depok | -0.741 |
| Kota Sukabumi | -0.911 |
| Kota Tasikmalaya | -2.375 |
| Kuningan | 0.625 |
| Majalengka | -2.068 |
| Pangandaran | -1.064 |
| Purwakarta | 0.553 |
| Subang | 1.187 |
| Sukabumi | 0.332 |
| Sumedang | 0.162 |
| Tasikmalaya | -1.298 |
Output diatas merupakan pengaruh acak dari setiap unit individu, nilai tersebut tersebut menunjukkan seberapa besar perbedaan nilai komponen error acak masing-masing unit indvidu terhadap nilai intersep umum.
Misalnya, nilai -1.474 pada kabupaten/kota 1 menunjukkan besarnya perbedaan komponen error acak kabupaten/kota 1 terhadap intersep umum. Artinya kabupaten/kota 1 memiliki intersep -1.474 yang lebih tinggi dari intersep umum. Begitu pula interpretasi untuk lainnya.
Diagnostik Model One Way REM
res.rem_ind <- residuals(rem_gls_ind)Normality (Jarque Bera Test)
Hipotesis:
H0 = Sisaan menyebar normal
H1 = Sisaan tidak
menyebar normal
library(tseries)
jarque.bera.test(res.rem_ind)##
## Jarque Bera Test
##
## data: res.rem_ind
## X-squared = 1.8211, df = 2, p-value = 0.4023Normality (Kolmogorov-Smirnov)
Hipotesis:
H0 = Sisaan menyebar normal
H1 = Sisaan tidak
menyebar normal
ks.test(res.rem_ind, "pnorm")##
## One-sample Kolmogorov-Smirnov test
##
## data: res.rem_ind
## D = 0.091169, p-value = 0.4832
## alternative hypothesis: two-sidedDalam menguji kenormalan, H0 ditolak baik menggunakan metode Jarque Bera Test maupun metode Kolmogorov-Smirnov sehingga mengimplikasikan bahwa asumsi kenormalan belum terpenuhi
Selain itu, kita akan melihat apakah pengujian asumsi sudah sesuai jika menggunakan metode secara eksploratif yaitu Histogram dan QQ-Plot.
Histogram
hist(res.rem_ind,
xlab = "Sisaan",
col = "#27D3D3",
breaks=30,
prob = TRUE)
lines(density(res.rem_ind), # density plot
lwd = 2, # thickness of line
col = "chocolate3")
Plot QQ-Norm
set.seed(1353)
res.rem_ind1 <- as.numeric(res.rem_ind)
qqnorm(res.rem_ind1,datax=T, col="blue")
qqline(rnorm(length(res.rem_ind1),mean(res.rem_ind1),sd(res.rem_ind1)),datax=T, col="red")
Auto-Correlation
Hipotesis:
H0 = Sisaan saling bebas
H1 = Sisaan tidak
saling bebas
adf.test(res.rem_ind, k=2) #alternatif : Terdapat autokorelasi## Warning in adf.test(res.rem_ind, k = 2): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: res.rem_ind
## Dickey-Fuller = -4.8605, Lag order = 2, p-value = 0.01
## alternative hypothesis: stationaryHomogen
Hipotesis:
H0 = Sisaan memiliki ragam homogen
H1 = Sisaan
memiliki ragam yang tidak homogen
bptest(rem_gls_ind)##
## studentized Breusch-Pagan test
##
## data: rem_gls_ind
## BP = 1.0464, df = 4, p-value = 0.9027Berdasarkan hasil uji diagnostik sisaan model REM Individual, didapatkan bahwa ketiga asumsi menghasilkan nilai p-value < 0.05, sehingga keputusannya adalah tolak H0. Dengan kata lain, asumsi normalitas belum terpenuhi, sisaan tidak saling bebas, dan ragam yang tidak homogen.
Pendugaan Model REM Individual
Pada package R, pendugaan model REM disediakan beberapa pilihan metode random yang dapat digunakan yaitu :
walhus: Wallace and Hussain (1969)swar: Swamy Arora (1972)amemiya1: Amemiya (1971)ht: Hausman and Taylor (1981)nerlove: Nerlove (1971)
rem.walhus <- update(rem_gls_ind, random.method = "walhus")
rem.amemiya <- update(rem_gls_ind, random.method = "amemiya")
rem.nerlove <- update(rem_gls_ind, random.method = "nerlove")
rem.models_ind <- list(swar = rem_gls_ind, walhus = rem.walhus,amemiya = rem.amemiya, nerlove = rem.nerlove)
sapply(rem.models_ind, coef)## swar walhus amemiya nerlove
## (Intercept) 22.2539377 23.7147695 19.6321495 17.8684766
## RLS 0.5693954 0.5232893 0.6454786 0.6857005
## TPAK -0.2986934 -0.3118769 -0.2743288 -0.2568723
## PDRB -0.2402032 -0.2441048 -0.2336119 -0.2292779
## PPM 0.2183753 0.1963715 0.2596414 0.2908538One Way Time
Pengujian pada model one way time dilakukan untuk membuktikan ada atau tidaknya pengaruh two-way sehingga juga mencoba melihat pengaruh waktu yang telah didapat, yaitu signifikan (tolak H0).
rem_gls_time <- plm(model, data = datapanel,
index = c("Kab.Kota","Tahun"),
effect = "time",
model = "random", random.method = "walhus")
summary(rem_gls_time)## Oneway (time) effect Random Effect Model
## (Wallace-Hussain's transformation)
##
## Call:
## plm(formula = model, data = datapanel, effect = "time", model = "random",
## random.method = "walhus", index = c("Kab.Kota", "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Effects:
## var std.dev share
## idiosyncratic 2.3260 1.5251 0.984
## time 0.0382 0.1955 0.016
## theta: 0.1677
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -3.17859 -1.30133 0.13558 0.96326 3.41608
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) 33.504701 4.630041 7.2364 4.608e-13 ***
## RLS 0.176768 0.200751 0.8805 0.3786
## TPAK -0.393177 0.052608 -7.4738 7.793e-14 ***
## PDRB -0.297669 0.066999 -4.4429 8.876e-06 ***
## PPM 0.042471 0.088164 0.4817 0.6300
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 420.62
## Residual Sum of Squares: 188.23
## R-Squared: 0.5525
## Adj. R-Squared: 0.52895
## Chisq: 93.8319 on 4 DF, p-value: < 2.22e-16kable(tidy(rem_gls_time), digits=3,
caption="GLS_Random Effect Model_Waktu")| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| (Intercept) | 33.505 | 4.630 | 7.236 | 0.000 |
| RLS | 0.177 | 0.201 | 0.881 | 0.379 |
| TPAK | -0.393 | 0.053 | -7.474 | 0.000 |
| PDRB | -0.298 | 0.067 | -4.443 | 0.000 |
| PPM | 0.042 | 0.088 | 0.482 | 0.630 |
tidy_ranef_time <- tidy(ranef(rem_gls_time, effect="time"))## Warning: 'tidy.numeric' is deprecated.
## See help("Deprecated")kable(tidy_ranef_time, caption = "Pengaruh Acak Waktu", col.names = c("Tahun", "Pengaruh Acak"))| Tahun | Pengaruh Acak |
|---|---|
| 2019 | -0.0696614 |
| 2020 | -0.0922757 |
| 2021 | 0.1619371 |
Berdasarkan hasil tersebut menunjukkan pengaruh efek waktu bernilai 0 pada seluruh periode amatan. Hal ini berarti efek waktu bersifat invariant atau tidak ada pengaruh spesifik waktu.
MLE (One Way REM)
Pendugaan model random effect model dengan menggunakan MLE dengan
menggunakan package lme4 dan lmerTest.
One Way Individual
library(lme4)## Loading required package: Matrix##
## Attaching package: 'Matrix'## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpacklibrary(lmerTest)##
## Attaching package: 'lmerTest'## The following object is masked from 'package:lme4':
##
## lmer## The following object is masked from 'package:stats':
##
## steprem_MLE_ind <- lmer(TPT ~ RLS + TPAK + PDRB + PPM+(1|Kab.Kota),
data=datapanel,
REML=FALSE)
summary(rem_MLE_ind)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: TPT ~ RLS + TPAK + PDRB + PPM + (1 | Kab.Kota)
## Data: datapanel
##
## AIC BIC logLik deviance df.resid
## 276.3 293.1 -131.2 262.3 74
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.91334 -0.55139 -0.08193 0.55548 2.19993
##
## Random effects:
## Groups Name Variance Std.Dev.
## Kab.Kota (Intercept) 1.8508 1.360
## Residual 0.7276 0.853
## Number of obs: 81, groups: Kab.Kota, 27
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 21.80091 4.84162 80.12791 4.503 2.25e-05 ***
## RLS 0.58328 0.23417 59.09307 2.491 0.0156 *
## TPAK -0.29456 0.05817 77.50555 -5.064 2.71e-06 ***
## PDRB -0.23903 0.03342 57.53700 -7.153 1.68e-09 ***
## PPM 0.22530 0.09655 73.90196 2.334 0.0223 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) RLS TPAK PDRB
## RLS -0.641
## TPAK -0.849 0.171
## PDRB -0.133 0.222 -0.044
## PPM -0.330 0.438 -0.038 0.325Kolom Groups menerangkan bahwa terdapat intersep unit individu yang bersifat random. Untuk dapat melihat nilai dari efek acak unit individu dengan menggunakan fungsi berikut :
coef(rem_MLE_ind)$Kab.Kota## (Intercept) RLS TPAK PDRB PPM
## Bandung 20.31787 0.5832841 -0.2945599 -0.2390336 0.2252982
## Bandung Barat 21.78039 0.5832841 -0.2945599 -0.2390336 0.2252982
## Bekasi 22.87311 0.5832841 -0.2945599 -0.2390336 0.2252982
## Bogor 24.12789 0.5832841 -0.2945599 -0.2390336 0.2252982
## Ciamis 20.27599 0.5832841 -0.2945599 -0.2390336 0.2252982
## Cianjur 23.77148 0.5832841 -0.2945599 -0.2390336 0.2252982
## Cirebon 23.02378 0.5832841 -0.2945599 -0.2390336 0.2252982
## Garut 20.51133 0.5832841 -0.2945599 -0.2390336 0.2252982
## Indramayu 22.22772 0.5832841 -0.2945599 -0.2390336 0.2252982
## Karawang 23.55386 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Bandung 22.73352 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Banjar 20.21730 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Bekasi 21.97025 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Bogor 22.54841 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Cimahi 23.08185 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Cirebon 21.50986 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Depok 21.06221 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Sukabumi 20.89537 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kota Tasikmalaya 19.35970 0.5832841 -0.2945599 -0.2390336 0.2252982
## Kuningan 22.42674 0.5832841 -0.2945599 -0.2390336 0.2252982
## Majalengka 19.68889 0.5832841 -0.2945599 -0.2390336 0.2252982
## Pangandaran 20.69107 0.5832841 -0.2945599 -0.2390336 0.2252982
## Purwakarta 22.37957 0.5832841 -0.2945599 -0.2390336 0.2252982
## Subang 23.00600 0.5832841 -0.2945599 -0.2390336 0.2252982
## Sukabumi 22.16655 0.5832841 -0.2945599 -0.2390336 0.2252982
## Sumedang 21.94720 0.5832841 -0.2945599 -0.2390336 0.2252982
## Tasikmalaya 20.47676 0.5832841 -0.2945599 -0.2390336 0.2252982ranef(rem_MLE_ind) #nilai yang disajikan merupakan perbedaan dengan base intersepnya## $Kab.Kota
## (Intercept)
## Bandung -1.4830430
## Bandung Barat -0.0205219
## Bekasi 1.0721947
## Bogor 2.3269751
## Ciamis -1.5249221
## Cianjur 1.9705652
## Cirebon 1.2228689
## Garut -1.2895875
## Indramayu 0.4268027
## Karawang 1.7529518
## Kota Bandung 0.9326063
## Kota Banjar -1.5836138
## Kota Bekasi 0.1693362
## Kota Bogor 0.7474952
## Kota Cimahi 1.2809377
## Kota Cirebon -0.2910560
## Kota Depok -0.7387077
## Kota Sukabumi -0.9055404
## Kota Tasikmalaya -2.4412135
## Kuningan 0.6258233
## Majalengka -2.1120246
## Pangandaran -1.1098416
## Purwakarta 0.5786551
## Subang 1.2050897
## Sukabumi 0.3656351
## Sumedang 0.1462843
## Tasikmalaya -1.3241493
##
## with conditional variances for "Kab.Kota"One Way Time
library(lme4)
library(lmerTest)
re_MLE_time = lmer(TPT ~ RLS + TPAK + PDRB + PPM+(1|Tahun),
data=datapanel,
REML=FALSE)
summary(re_MLE_time)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: TPT ~ RLS + TPAK + PDRB + PPM + (1 | Tahun)
## Data: datapanel
##
## AIC BIC logLik deviance df.resid
## 313.2 330.0 -149.6 299.2 74
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1848 -0.8614 0.1263 0.6478 2.1672
##
## Random effects:
## Groups Name Variance Std.Dev.
## Tahun (Intercept) 0.06242 0.2498
## Residual 2.30729 1.5190
## Number of obs: 81, groups: Tahun, 3
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 33.72204 4.48192 80.84411 7.524 6.52e-11 ***
## RLS 0.16344 0.19490 79.23953 0.839 0.4042
## TPAK -0.39349 0.05079 77.46707 -7.747 3.02e-11 ***
## PDRB -0.30451 0.06967 3.13534 -4.371 0.0202 *
## PPM 0.03518 0.08577 74.99248 0.410 0.6829
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) RLS TPAK PDRB
## RLS -0.723
## TPAK -0.867 0.306
## PDRB -0.042 0.050 -0.046
## PPM -0.484 0.681 0.078 0.103Kolom Groups menerangkan bahwa terdapat intersep unit individu yang bersifat random. Untuk dapat melihat nilai dari efek acak unit individu dengan menggunakan fungsi berikut :
ranef(re_MLE_time)## $Tahun
## (Intercept)
## 2019 -0.09039558
## 2020 -0.13753916
## 2021 0.22793474
##
## with conditional variances for "Tahun"Tambahan
Selain dengan package lmer, pemodelan data panel dengan
pendugaan MLE dapat dilakukan dengan menggunakan package
pglm dan nlme.
Pemodelan MLE Lainnya
library(pglm)## Loading required package: maxLik## Loading required package: miscTools##
## Please cite the 'maxLik' package as:
## Henningsen, Arne and Toomet, Ott (2011). maxLik: A package for maximum likelihood estimation in R. Computational Statistics 26(3), 443-458. DOI 10.1007/s00180-010-0217-1.
##
## If you have questions, suggestions, or comments regarding the 'maxLik' package, please use a forum or 'tracker' at maxLik's R-Forge site:
## https://r-forge.r-project.org/projects/maxlik/rem.mle <- pglm(TPT ~ RLS + TPAK + PDRB + PPM,
data = datapanel, model="random",
effect="individual",
family = gaussian)
summary(rem.mle)## --------------------------------------------
## Maximum Likelihood estimation
## Newton-Raphson maximisation, 5 iterations
## Return code 8: successive function values within relative tolerance limit (reltol)
## Log-Likelihood: -131.1508
## 7 free parameters
## Estimates:
## Estimate Std. error t value Pr(> t)
## (Intercept) 21.80088 5.19269 4.198 2.69e-05 ***
## RLS 0.58329 0.24102 2.420 0.0155 *
## TPAK -0.29456 0.06065 -4.857 1.19e-06 ***
## PDRB -0.23903 0.03376 -7.080 1.44e-12 ***
## PPM 0.22530 0.10075 2.236 0.0253 *
## sd.id 1.36043 0.22452 6.059 1.37e-09 ***
## sd.idios 0.85297 0.08247 10.343 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## --------------------------------------------Dengan menggunakan library nlme
library(nlme)##
## Attaching package: 'nlme'## The following object is masked from 'package:lme4':
##
## lmList## The following object is masked from 'package:dplyr':
##
## collapsereML <- lme(TPT ~ RLS + TPAK + PDRB + PPM, data = datapanel,
random = ~1 | Kab.Kota)
summary(reML)## Linear mixed-effects model fit by REML
## Data: datapanel
## AIC BIC logLik
## 289.9202 306.2353 -137.9601
##
## Random effects:
## Formula: ~1 | Kab.Kota
## (Intercept) Residual
## StdDev: 1.446746 0.8699816
##
## Fixed effects: TPT ~ RLS + TPAK + PDRB + PPM
## Value Std.Error DF t-value p-value
## (Intercept) 21.419854 5.006552 50 4.278364 0.0001
## RLS 0.594777 0.244155 50 2.436069 0.0185
## TPAK -0.291063 0.060251 50 -4.830830 0.0000
## PDRB -0.238062 0.034154 50 -6.970299 0.0000
## PPM 0.231171 0.100169 50 2.307817 0.0252
## Correlation:
## (Intr) RLS TPAK PDRB
## RLS -0.638
## TPAK -0.847 0.165
## PDRB -0.133 0.223 -0.046
## PPM -0.323 0.423 -0.042 0.329
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.88933023 -0.56452463 -0.07210676 0.54711085 2.13267733
##
## Number of Observations: 81
## Number of Groups: 27#random komponen
ranef(reML)## (Intercept)
## Bandung -1.48960166
## Bandung Barat -0.01498191
## Bekasi 1.09730261
## Bogor 2.35979445
## Ciamis -1.54097978
## Cianjur 1.98518845
## Cirebon 1.23310230
## Garut -1.29132857
## Indramayu 0.41874983
## Karawang 1.77994962
## Kota Bandung 0.94156177
## Kota Banjar -1.59238048
## Kota Bekasi 0.16698469
## Kota Bogor 0.75580756
## Kota Cimahi 1.29759939
## Kota Cirebon -0.30509114
## Kota Depok -0.73530326
## Kota Sukabumi -0.90004044
## Kota Tasikmalaya -2.49608900
## Kuningan 0.62597110
## Majalengka -2.14834132
## Pangandaran -1.14861934
## Purwakarta 0.59991033
## Subang 1.21951543
## Sukabumi 0.39397526
## Sumedang 0.13293880
## Tasikmalaya -1.34559471Perbandingan Model
Model dibandingkan salah satunya adalah dilihat dari AIC-nya. Agar dapat menjalankan AIC pada object PLM maka perlu dibuat fungsi logLik.plm terlebih dahulu.
logLik.plm <- function(object){
out <- -plm::nobs(object) * log(2 * var(object$residuals) * pi)/2 - deviance(object)/(2 * var(object$residuals))
attr(out,"df") <- nobs(object) - object$df.residual
attr(out,"nobs") <- plm::nobs(summary(object))
return(out)
}#AIC
AIC1 <- rbind(AIC(rem_gls_ind), AIC(reML))
#koefisien
koef1 <- rbind.data.frame(coef(rem_gls_ind), fixef(reML))
rownames(koef1) <- c("REM_GLS", "REM_MLE")
colnames(koef1) <- c("intersep", "INFRA", "PEKO", "INFL",
"PGGR")
koef1$AIC <- AIC1
koef1## intersep INFRA PEKO INFL PGGR AIC
## REM_GLS 22.25394 0.5693954 -0.2986934 -0.2402032 0.2183753 216.4844
## REM_MLE 21.41985 0.5947773 -0.2910629 -0.2380618 0.2311710 289.9202Berdasarkan hasil pendugaan dengan metode GLS dan MLE memberikan nilai dugaan koefisien yang tidak jauh berbeda. Walaupun begitu, pendugaan dengan metode GLS cenderung lebih baik karena memiliki nilai AIC yang lebih kecil yaitu 216.484 dibandingkan dengan metode MLE yang memiliki AIC sebesar 289.920.
Pemilihan FEM dan REM (Model Individu)
Pemilihan FEM dan REM dilakukan dengan menggunakan Uji Hausman. Uji
Hausman sendiri memiliki default yang berisi
phtest(fixed effect model, random effect model) sehingga
:
Hipotesis :
H0 : Random Effect Model
(REM)
H1 : Fixed Effect Model (FEM)
FEM dengan REM GLS
kable(tidy(phtest(fem, rem_gls_ind)), caption="Uji Hausman FEM_ind VS REM_ind GLS")| statistic | p.value | parameter | method | alternative |
|---|---|---|---|---|
| 7.762078 | 0.1006928 | 4 | Hausman Test | one model is inconsistent |
FEM dengan REM MLE package pglm
kable(tidy(phtest(fem, rem.mle)),caption="Uji Hausman FEM_ind VS REM_ind MLE")| statistic | p.value | parameter | method | alternative |
|---|---|---|---|---|
| 10.62814 | 0.0310769 | 4 | Hausman Test | one model is inconsistent |
Berdasarkan pengujian Hausman pada kedua perbandingan, yaitu FEM vs REM GLS dan FEM vs REM MLE memberikan hasil tolak H0, sehingga model FEM lebih tepat digunakan daripada menggunakan model REM (baik menggunakan metode GLS maupun MLE).
Uji Hausman Manual
# fem dengan rem_MLE_ind
# matrix covariance dari masing-masing model
V1 <- fem$vcov
V2 <- rem_MLE_ind@vcov_beta[2,2]
# Hitung perbedaan slope koefisien masing-masing model
bdiff <- coef(fem) - coef(rem_gls_ind)[2]
# Hitung perbedaan dari matrix covariance
Vdiff <- V1 - V2
# Hitung Hausman statistic
H <- t(bdiff) %*% solve(Vdiff) %*% bdiff
# Hitung nilai pvalue
pval <- 1 - pchisq(H, length(bdiff))
# Print hasil Hausman test
cat("Hausman test:<br>")## Hausman test:<br>cat("H = ", H, "<br>")## H = 15.67926 <br>cat("p-value = ", pval, "<br>")## p-value = 0.003481182 <br>Karena nilai p-value yang diperoleh lebih kecil dari alpha (5%), maka Terima H0 yang menandakan bahwa model yang dipilih dari perhitungan manual adalah model FEM. Namun, perlu dicatat lebih lanjut bahwa perhitungan manual cenderung memiliki banyak kesalahan dibandingkan dengan pengujian secara langsung sehingga hasil kurang konsisten.
Selain itu, jika terindikasi terdapat gejala heteroskedastisitas, pengujian Hausman dapat dilakukan dengan menggunakan uji Hausman robust dengan menyertakan matriks robust covariance untuk mendapatkan hasil pengujian yang konsisten.
kable(tidy(phtest(fem, rem_gls_ind, method = "aux", vcov = vcovHC)),caption="Uji Hausman Robust FEM_ind VS REM_ind GLS")| statistic | p.value | parameter | method | alternative |
|---|---|---|---|---|
| 7.762078 | 0.1006928 | 4 | Hausman Test | one model is inconsistent |
Berdasarkan pengujian robust Hausman, keputusan yang diambil adalah Terima H0 yang berimplikasi model yang lebih baik adalah model REM.
Kemudian, tahapan selanjutnya adalah perlu memeriksa antara model pengaruh acak dengan common effect model dengan uji pengganda lagrange.
cem <- plm(model, datapanel, model="pooling")
c.tw <- plmtest(cem, effect = "twoways", type = c("bp"))
c.ind <- plmtest(cem, effect = "individual", type = c("bp"))
c.time <-plmtest(cem, effect = "time", type = c("bp"))
tidy_tests <- bind_rows(tidy(c.tw), tidy(c.ind), tidy(c.time))
colnames(tidy_tests) <- c("Chisq", "p-value", "df","Jenis-Uji", "H1")
kable(tidy_tests, digits=2,caption = "Hasil Uji LM Breusch-Pagan")| Chisq | p-value | df | Jenis-Uji | H1 |
|---|---|---|---|---|
| 32.32 | 0.0 | 2 | Lagrange Multiplier Test - two-ways effects (Breusch-Pagan) | significant effects |
| 32.04 | 0.0 | 1 | Lagrange Multiplier Test - (Breusch-Pagan) | significant effects |
| 0.27 | 0.6 | 1 | Lagrange Multiplier Test - time effects (Breusch-Pagan) | significant effects |
Berdasarkan pengecekan dengan uji pengganda lagrange, dapat ditunjukkan bahwa terdapat efek panel (efek two ways, dan efek individu signifikan pada taraf 5%) sehingga pemodelan dengan pooling saja tidak lah cukup.
Dengan demikian model yang terpilih adalah model Two Ways REM.
GLS : Two Ways Random Effect Model
rem_gls_two<- plm(model, data = datapanel,
index = c("Kab.Kota", "Tahun"),
effect = "twoways", model = "random", random.method = "walhus")
summary(rem_gls_two)## Twoways effects Random Effect Model
## (Wallace-Hussain's transformation)
##
## Call:
## plm(formula = model, data = datapanel, effect = "twoways", model = "random",
## random.method = "walhus", index = c("Kab.Kota", "Tahun"))
##
## Balanced Panel: n = 27, T = 3, N = 81
##
## Effects:
## var std.dev share
## idiosyncratic 0.71722 0.84689 0.305
## individual 1.54036 1.24111 0.654
## time 0.09779 0.31271 0.042
## theta: 0.6335 (id) 0.5378 (time) 0.4711 (total)
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -1.95584 -0.51412 -0.11184 0.65524 2.19514
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) 27.278421 5.087325 5.3620 8.229e-08 ***
## RLS 0.222026 0.260290 0.8530 0.3937
## TPAK -0.302118 0.057075 -5.2933 1.201e-07 ***
## PDRB -0.269474 0.062831 -4.2888 1.796e-05 ***
## PPM 0.021265 0.119113 0.1785 0.8583
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 98.264
## Residual Sum of Squares: 56.763
## R-Squared: 0.42235
## Adj. R-Squared: 0.39194
## Chisq: 55.5666 on 4 DF, p-value: 2.4717e-11Nilai intersep pada model pengaruh acak sebesar 27.278421, nilai ini menggambarkan rataan umum intersep model. Koefisien RLS, TPAK, PDRB, PPM masing-masing secara berurutan sebesar 0.222026, -0.302118, -0.269474, 0.021265 dengan nilai R^2-Adj sebesar 39,19%. Nilai dari p-value < α (5%) menunjukkan bahwa x memiliki pengaruh yang signifikan terhadap Kemiskinan pada model REM Two Ways dengan taraf nyata 5%.
Diagnostik Model Two Ways REM
Uji Normalitas
Hipotesis :
H0 = Sisaan menyebar normal
H1 = Sisaan tidak
menyebar normal
res.rem.two <- residuals(rem_gls_two)
#Normality
library(tseries)
jarque.bera.test(res.rem.two)##
## Jarque Bera Test
##
## data: res.rem.two
## X-squared = 0.86827, df = 2, p-value = 0.6478#Kolmogorov-smirnov
ks.test(res.rem.two, "pnorm")##
## One-sample Kolmogorov-Smirnov test
##
## data: res.rem.two
## D = 0.081034, p-value = 0.6324
## alternative hypothesis: two-sidedUji normalitas baik dengan menggunakan metode Jarque Bera Test maupun metode Kolmogorov-Smirnov memberikan hasil Terima H0 yang bermakna bahwa asumsi kenormalan terpenuhi. Berikutnya kita cek sebaran secara eksploratif untuk menguatkan hasil uji yang didapat dengan Histogram dan QQ-Plot.
#Histogram
hist(res.rem.two,
xlab = "sisaan",
col = "#27D3D3",
breaks=30,
prob = TRUE)
lines(density(res.rem.two), # density plot
lwd = 2, # thickness of line
col = "chocolate3")
#Plot QQ-Norm
set.seed(1353)
res.rem1 <- as.numeric(res.rem.two)
qqnorm(res.rem1, datax=T, col="blue")
qqline(rnorm(length(res.rem1),mean(res.rem1),sd(res.rem1)), datax=T, col="red")
Uji Autokorelasi
Hipotesis:
H0 = Sisaan saling bebas
H1 = Sisaan tidak
saling bebas
#Autocorrelation
adf.test(res.rem.two, k=2) #alternatif : Terdapat autokorelasi## Warning in adf.test(res.rem.two, k = 2): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: res.rem.two
## Dickey-Fuller = -5.2479, Lag order = 2, p-value = 0.01
## alternative hypothesis: stationaryUji Homoskedastisitas
Hipotesis:
H0 = Sisaan memiliki ragam homogen
H1 = Sisaan
memiliki ragam yang tidak homogen
#Homogen
bptest(rem_gls_two)##
## studentized Breusch-Pagan test
##
## data: rem_gls_two
## BP = 1.0464, df = 4, p-value = 0.9027Berdasarkan hasil uji diagnostik sisaan model REM Two Ways, didapatkan bahwa ketiga pengujian asumsi menghasilkan suatu asumsi bernilai p-value < 0.05, sehingga keputusannya adalah Tolak H0 pada uji Autokorelasi. Dengan kata lain, asumsi normalitas terpenuhi, sisaan tidak saling bebas, dan ragam yang homogen.
Pendugaan Model REM
Pada package R, pendugaan model rem disediakan beberapa pilihan metode random yang dapat digunakan yaitu :
’walhus’ : Wallace and Hussain (1969)
’swar’ : Swamy Arora (1972)
’amemiya1’ : Amemiya (1971)
’ht’ : Hausman and Taylor (1981)
’nerlove’ : Nerlove (1971)rem.walhus <- update(rem_gls_two, random.method = "walhus")
rem.amemiya <- update(rem_gls_two, random.method = "amemiya")
rem.nerlove <- update(rem_gls_two, random.method = "nerlove")library(texreg)## Version: 1.38.6
## Date: 2022-04-06
## Author: Philip Leifeld (University of Essex)
##
## Consider submitting praise using the praise or praise_interactive functions.
## Please cite the JSS article in your publications -- see citation("texreg").##
## Attaching package: 'texreg'## The following object is masked from 'package:tidyr':
##
## extractscreenreg(list("Swamy-Arora" = rem_gls_two,"Wallace-Hussain" = rem.walhus, "Amemiya" = rem.amemiya, "Nerlove"= rem.nerlove ), digits = 5)##
## ======================================================================
## Swamy-Arora Wallace-Hussain Amemiya Nerlove
## ----------------------------------------------------------------------
## (Intercept) 27.27842 *** 27.27842 *** 31.93686 *** 33.71042 ***
## (5.08733) (5.08733) (5.83552) (5.95884)
## RLS 0.22203 0.22203 -0.56142 -0.66839
## (0.26029) (0.26029) (0.42124) (0.43447)
## TPAK -0.30212 *** -0.30212 *** -0.17410 ** -0.16045 **
## (0.05708) (0.05708) (0.05681) (0.05610)
## PDRB -0.26947 *** -0.26947 *** -0.27046 ** -0.26845 **
## (0.06283) (0.06283) (0.09167) (0.09190)
## PPM 0.02127 0.02127 -0.70482 ** -0.91211 ***
## (0.11911) (0.11911) (0.24752) (0.27198)
## ----------------------------------------------------------------------
## s_idios 0.84689 0.84689 0.67672 0.54221
## s_id 1.24111 1.24111 4.11724 4.21452
## s_time 0.31271 0.31271 1.56592 1.92447
## R^2 0.42235 0.42235 0.29142 0.30910
## Adj. R^2 0.39194 0.39194 0.25413 0.27274
## Num. obs. 81 81 81 81
## ======================================================================
## *** p < 0.001; ** p < 0.01; * p < 0.05Pada tabel diatas dapat diamati bahwa pendugaan REM Two Ways Effect dengan metode random swar, walhus, anemiya , maupun nerlove memberikan dugaan parameter tidak jauh berbeda.
Efek Two Ways dapat diekstrak dengan fungsi
ranef() sama halnya pada efek One Way.
#Efek Individu pada REM Two Ways
tidy_ranef_ind <- tidy(ranef(rem_gls_two, effect="individual"))## Warning: 'tidy.numeric' is deprecated.
## See help("Deprecated")kable(tidy_ranef_ind, digits=3,caption = "Pengaruh Acak Individu", col.names = c("Kab.Kota", "Pengaruh Acak Individu"))| Kab.Kota | Pengaruh Acak Individu |
|---|---|
| Bandung | -1.627 |
| Bandung Barat | 0.213 |
| Bekasi | 0.562 |
| Bogor | 2.065 |
| Ciamis | -1.850 |
| Cianjur | 1.843 |
| Cirebon | 1.334 |
| Garut | -1.293 |
| Indramayu | 0.502 |
| Karawang | 1.468 |
| Kota Bandung | 0.812 |
| Kota Banjar | -1.882 |
| Kota Bekasi | 0.221 |
| Kota Bogor | 0.928 |
| Kota Cimahi | 1.277 |
| Kota Cirebon | 0.251 |
| Kota Depok | -1.029 |
| Kota Sukabumi | -0.731 |
| Kota Tasikmalaya | -1.411 |
| Kuningan | 1.100 |
| Majalengka | -1.806 |
| Pangandaran | -1.141 |
| Purwakarta | 0.357 |
| Subang | 0.894 |
| Sukabumi | -0.244 |
| Sumedang | 0.439 |
| Tasikmalaya | -1.253 |
#Efek Waktu
tidy_ranef_time <- tidy(ranef(rem_gls_two, effect="time"))## Warning: 'tidy.numeric' is deprecated.
## See help("Deprecated")kable(tidy_ranef_time, digits=3,caption = "Pengaruh Acak Waktu", col.names = c("Tahun", "Pengaruh Acak Waktu"))| Tahun | Pengaruh Acak Waktu |
|---|---|
| 2019 | -0.261 |
| 2020 | -0.143 |
| 2021 | 0.404 |
Pemilihan FEM dan REM (Model Two Ways)
Hipotesis :
H0 : Random Effect Model
(REM)
H1 : Fixed Effect Model (FEM)
phtest(fem.twoway,rem_gls_two)##
## Hausman Test
##
## data: TPT ~ RLS + TPAK + PDRB + PPM
## chisq = 32.002, df = 4, p-value = 1.911e-06
## alternative hypothesis: one model is inconsistentBerdasarkan uji Hausman, nilai p-value < alpha (5%) sehingga Tolak H0. Jadi, FEM Two Ways merupakan model yang tepat.
Final Model
final_model <- fem.twoway
kable(tidy(final_model),caption="Model Pengaruh Acak Efek Individu")| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| RLS | -0.7448960 | 0.5627794 | -1.323602 | 0.1919050 |
| TPAK | -0.1280698 | 0.0667451 | -1.918790 | 0.0609681 |
| PDRB | -0.2571995 | 0.1093581 | -2.351902 | 0.0228283 |
| PPM | -1.6977524 | 0.4338827 | -3.912929 | 0.0002868 |
tidy_fixef_twoways <- tidy(fixef(fem.twoway, effect="twoways"))## Warning: 'tidy.numeric' is deprecated.
## See help("Deprecated")kable(tidy_fixef_twoways,digits=3, caption = "Pengaruh Acak Waktu dan Individu", col.names = c("Kab.Kota", "Pengaruh Acak Waktu dan Individu"))| Kab.Kota | Pengaruh Acak Waktu dan Individu |
|---|---|
| Bandung-2019 | 32.60813 |
| Bandung-2020 | 34.77882 |
| Bandung-2021 | 36.44608 |
| Bandung Barat-2019 | 41.07040 |
| Bandung Barat-2020 | 43.24109 |
| Bandung Barat-2021 | 44.90835 |
| Bekasi-2019 | 31.86688 |
| Bekasi-2020 | 34.03758 |
| Bekasi-2021 | 35.70483 |
| Bogor-2019 | 37.86885 |
| Bogor-2020 | 40.03955 |
| Bogor-2021 | 41.70680 |
| Ciamis-2019 | 31.60095 |
| Ciamis-2020 | 33.77165 |
| Ciamis-2021 | 35.43890 |
| Cianjur-2019 | 40.37236 |
| Cianjur-2020 | 42.54305 |
| Cianjur-2021 | 44.21031 |
| Cirebon-2019 | 42.26718 |
| Cirebon-2020 | 44.43788 |
| Cirebon-2021 | 46.10513 |
| Garut-2019 | 37.62123 |
| Garut-2020 | 39.79192 |
| Garut-2021 | 41.45918 |
| Indramayu-2019 | 41.70644 |
| Indramayu-2020 | 43.87713 |
| Indramayu-2021 | 45.54439 |
| Karawang-2019 | 37.81211 |
| Karawang-2020 | 39.98280 |
| Karawang-2021 | 41.65005 |
| Kota Bandung-2019 | 32.20763 |
| Kota Bandung-2020 | 34.37833 |
| Kota Bandung-2021 | 36.04558 |
| Kota Banjar-2019 | 30.86297 |
| Kota Banjar-2020 | 33.03366 |
| Kota Banjar-2021 | 34.70091 |
| Kota Bekasi-2019 | 32.58958 |
| Kota Bekasi-2020 | 34.76027 |
| Kota Bekasi-2021 | 36.42753 |
| Kota Bogor-2019 | 36.95164 |
| Kota Bogor-2020 | 39.12233 |
| Kota Bogor-2021 | 40.78959 |
| Kota Cimahi-2019 | 34.82747 |
| Kota Cimahi-2020 | 36.99816 |
| Kota Cimahi-2021 | 38.66541 |
| Kota Cirebon-2019 | 40.29081 |
| Kota Cirebon-2020 | 42.46150 |
| Kota Cirebon-2021 | 44.12876 |
| Kota Depok-2019 | 27.90673 |
| Kota Depok-2020 | 30.07742 |
| Kota Depok-2021 | 31.74468 |
| Kota Sukabumi-2019 | 36.83875 |
| Kota Sukabumi-2020 | 39.00945 |
| Kota Sukabumi-2021 | 40.67670 |
| Kota Tasikmalaya-2019 | 43.04280 |
| Kota Tasikmalaya-2020 | 45.21349 |
| Kota Tasikmalaya-2021 | 46.88075 |
| Kuningan-2019 | 44.87910 |
| Kuningan-2020 | 47.04980 |
| Kuningan-2021 | 48.71705 |
| Majalengka-2019 | 38.17744 |
| Majalengka-2020 | 40.34813 |
| Majalengka-2021 | 42.01539 |
| Pangandaran-2019 | 33.77977 |
| Pangandaran-2020 | 35.95047 |
| Pangandaran-2021 | 37.61772 |
| Purwakarta-2019 | 37.09348 |
| Purwakarta-2020 | 39.26417 |
| Purwakarta-2021 | 40.93142 |
| Subang-2019 | 37.55587 |
| Subang-2020 | 39.72656 |
| Subang-2021 | 41.39382 |
| Sukabumi-2019 | 33.49619 |
| Sukabumi-2020 | 35.66688 |
| Sukabumi-2021 | 37.33413 |
| Sumedang-2019 | 39.79207 |
| Sumedang-2020 | 41.96276 |
| Sumedang-2021 | 43.63002 |
| Tasikmalaya-2019 | 37.12072 |
| Tasikmalaya-2020 | 39.29142 |
| Tasikmalaya-2021 | 40.95867 |
Model Akhir
\[ \begin{aligned} I\hat nv_{it} &= -0.7448960RLS-0.1280698TPAK-0.2571995PDRB-1.6977524PPM\\ \end{aligned} \]
\(RLS\)
Setiap penurunan satu satuan RLS(Rata-Rata Lama Sekolah di Atas 15 Tahun) maka akan menyebabkan Tingkat Pengangguran Terbuka turun sebesar -0.7448960 dengan mengganggap peubah lain konstan (apabila memasukkan peubah lain).
$TPAK
Setiap penurunan satu satuan TPAK(Tingkat Partisipasi Angkatan Kerja) maka akan menyebabkan Tingkat Pengangguran Terbuka turun sebesar -0.1280698 dengan mengganggap peubah lain konstan (apabila memasukkan peubah lain).
\(PDRB\)
Setiap penurunan satu satuan PDRB(Produk Domestik Regional Bruto) maka akan menyebabkan Tingkat Pengangguran Terbuka turun sebesar -0.2571995 dengan mengganggap peubah lain konstan (apabila memasukkan peubah lain).
\(PPM\)
Setiap penurunan satu satuan PPM(Persentase Penduduk Miskin) maka akan menyebabkan Tingkat Pengangguran Terbuka turun sebesar -1.6977524 dengan mengganggap peubah lain konstan (apabila memasukkan peubah lain).