Project ADP

Afris Setiya Intan Amanda

12/5/2023

Analisis Regresi Data Panel dalam Pemodelan Angka Buta Huruf di Indonesia Tahun 2016-2020 Anggota Kelompok 11 STA1353 :
1. Afris Setiya Intan Amanda    (G1401201018)
2. Duwi Sulistianingsih      (G1401201033)
3. Muhammad Habib Ghozi    (G1401201065)
4. Naura Tirza Ardhani      (G1401201073)
5. Dhea Puspita Adinda    (G1401201090)

Outline

Pendahuluan

Tinjauan Pustaka

Metodologi

Hasil dan Pembahasan

Simpulan

Daftar Pustaka

Pendahuluan

Terlampir pada PPT.

Tinjauan Pustaka

Model Gabungan

Model gabungan atau common effect model (CEM) atau pooled least square (PLS) merupakan pendekatan model data panel yang sederhana karena mengkombinasikan data time series dan cross section, tanpa memperhatikan pengaruh spesifik waktu maupun individu. Koefisien regresi (intersep ataupun kemiringan) diasumsikan konstan antar individu dan waktu. Metode ini bisa menggunakan pendekatan ordinary least square (OLS) atau metode kuadrat terkecil (MKT) untuk mengestimasi model data panel (Nandita et al. 2019). Persamaan model gabungan dinyatakan sebagai berikut.

\[y_{it}=\alpha+x^{'}_{it} Q+\epsilon_{it}\]

Model Pengaruh Tetap

Model pengaruh tetap atau fixed effect model (FEM) merupakan model yang mengasumsikan antara unit individu atau waktu memiliki perilaku yang berbeda, terlihat dari nilai intersep yang berbeda untuk setiap unit individu atau waktu, tetapi konstan pada nilai koefisien kemiringan dan koefisien regresi antara unit individu maupun waktu (Gujarati dan Porter 2009). Pendugaan parameter model pengaruh tetap dapat menggunakan Metode Kuadrat Terkecil Peubah Boneka (least square dummy variable) dan MKT. Persamaan model pengaruh tetap dinyatakan sebagai berikut.

\[y_{it}=\alpha^{*}+x^{'}_{it} Q+\epsilon_{it}\]

Model Pengaruh Acak

Menurut Baltagi (2011), model pengaruh acak atau random effect model digunakan ketika individu amatan mengikuti kaidah pengacakan dari sejumlah populasi yang besar sehingga pengaruh pada setiap individu bersifat acak. Pendugaan parameter pada model pengaruh acak yaitu dengan metode kuadrat terkecil terampat (generalize least square). Model ini memiliki asumsi bahwa tidak ada korelasi antara pengaruh spesifik individu dan pengaruh spesifik waktu dengan peubah bebas sehingga komponen sisaan dari kedua pengaruh spesifik dimasukkan ke dalam model. Persamaan model random effect model dapat dinyatakan sebagai berikut

\[y_{it}=\alpha^{*}+x^{'}_{it} Q+W_{it}\] \[W_{it}=\mu_{i}+\epsilon_{it}\] \[atau\] \[W_{it}=\delta_{t}+\epsilon_{it}\]

Uji Chow

Menurut Ghozali dan Ratmono (2013), uji chow digunakan untuk memilih pendekatan yang lebih baik antara model gabungan dengan model pengaruh tetap.Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Model gabungan
$H_1$ : Model pengaruh tetap
Statistik uji Chow adalah sebagai berikut:
\[F_{c}=\frac{(RSS_1-RSS_2)/(n-1)}{RSS_2/(nT-n-p)}\]
dengan $RSS_1$ adalah jumlah kuadrat galat hasil pendugaan untuk model gabungan, sedangkan $RSS_2$ adalah jumlah kuadrat galat untuk hasil pendugaan model pengaruh tetap, $n$ adalah banyaknya unit individu, $T$ adalah banyaknya unit waktu, dan $p$ adalah banyaknya peubah penjelas. Keputusan tolak $H_0$ (model pengaruh tetap terpilih) apabila$F_0 > _{n-1,n(𝑇−1)−𝑝}$.

Uji Hausman

Uji spesifikasi Hausman membandingkan model pengaruh tetap dan model pengaruh acak. Jika hipotesis nol yang menyatakan tidak ada korelasi antara pengaruh individu dengan regresor tidak ditolak, model pengaruh random disarankan daripada pengaruh tetap (Susanti 2013). Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Model pengaruh acak
$H_1$ : Model pengaruh tetap
Statistik uji Hausman adalah sebagai berikut:
\[W=(b-𝜷)^{'}[var(b-𝜷)]-1(b-𝜷)\]
dengan $b$ merupakan vektor koefisien kemiringan model pengaruh tetap dan $β$ adalah vektor koefisien kemiringan model pengaruh acak. Keputusan tolak C apabila $W > X^{2}_{(p-1;α)^{2}}$ atau p-value < α dengan $p$ adalah banyaknya peubah penjelas. Jika tolak $H_0$ , berarti model pengaruh tetap terpilih, jika tak tolak $H_0$ , berarti model pengaruh acak terpilih (Greene 2003).

Multikolinieritas

Uji multikolinearitas dilakukan untuk mengetahui korelasi di antara peubah penjelas. Pendeteksian pada multikolinearitas dapat dilakukan dengan memeriksa nilai variance inflation factor (VIF). Suatu model regresi dikatakan bebas dari multikolinearitas jika mempunyai nilai VIF disekitar angka 1. Jika nilai VIF lebih dari 10, maka terjadi pengaruh serius yang menyebabkan peubah yang mengalami hal tersebut tidak dapat digunakan dalam analisis (Akinwande et al. 2015).

Metodologi

Data

Data yang digunakan berasal dari Data Publikasi Badan Pusat Statistik pada situs web bps.go.id. Data tersebut merupakan data sekunder mengenai Angka Buta Huruf di Indonesia tahun 2016-2020. Data tersebut terdiri dari 34 individu dan 7 peubah bebas. Adapun peubah-peubah yang digunakan dalam penelitian ini, dapat dilihat pada Tabel 1.

Prosedur

Analisis data dilakukan dengan bantuan software R Studio dan Google Data Studio. Prosedur analisis data pada penelitian ini dapat dilihat pada diagram alir yang tertera pada Gambar di bawah ini.

Hasil dan Pembahasan

Library

Library yang digunakan dalam data yang dianalisis ini adalah sebagai berikut.

library(dplyr)

## ## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats': ## ## filter, lag

## The following objects are masked from 'package:base': ## ## intersect, setdiff, setequal, union

library(ggplot2) library(gplots)

## ## Attaching package: 'gplots'

## The following object is masked from 'package:stats': ## ## lowess

library(ggcorrplot) library(DataExplorer) library(lme4)

## Loading required package: Matrix

library(readr) library(readxl) library(car)

## Loading required package: carData

## ## Attaching package: 'car'

## The following object is masked from 'package:dplyr': ## ## recode

library(lmerTest)

## ## Attaching package: 'lmerTest'

## The following object is masked from 'package:lme4': ## ## lmer

## The following object is masked from 'package:stats': ## ## step

library(plm)

## ## Attaching package: 'plm'

## The following objects are masked from 'package:dplyr': ## ## between, lag, lead

library(kableExtra)

## ## Attaching package: 'kableExtra'

## The following object is masked from 'package:dplyr': ## ## group_rows

library(broom) library(performance) library(tseries)

## Registered S3 method overwritten by 'quantmod': ## method from ## as.zoo.data.frame zoo

library(lmtest)

## Loading required package: zoo

## ## Attaching package: 'zoo'

## The following objects are masked from 'package:base': ## ## as.Date, as.Date.numeric

Input Data

Data yang digunakan pada penelitian ini adalah sebagai berikut.

data.abhifull <- read.csv("D:/MY COLLEGE/SEMESTER 6/ADP/DATA/ABHIFULL.csv", sep=";") kableExtra::kable(head(data.abhifull), caption = 'Subset Data ABH di Indonesia 2016-2020')

Subset Data ABH di Indonesia 2016-2020

Provinsi Tahun X1 X2 X3 X4 X5 X6 X7 Y

ACEH 2016 22835.29 848.44 35.23 52.97 98.16 7.57 25.78 2.26

ACEH 2017 23362.90 872.61 44.83 54.20 98.54 6.57 24.85 1.99

ACEH 2018 24013.79 839.49 56.89 59.05 99.10 6.34 30.18 5.04

ACEH 2019 24842.30 819.44 65.16 57.75 99.12 6.17 29.33 0.79

ACEH 2020 25018.28 814.91 71.97 59.60 99.03 6.59 27.12 5.68

SUMATERA UTARA 2016 32885.09 1455.95 40.44 54.28 96.57 5.84 22.88 2.06

data.abhitinggi <- read.csv("D:/MY COLLEGE/SEMESTER 6/ADP/DATA/ABHITINGGI.csv", sep=";") kableExtra::kable(head(data.abhitinggi), caption = 'Subset Data ABH Tinggi di Indonesia 2016-2020')

Subset Data ABH Tinggi di Indonesia 2016-2020

Provinsi Tahun X1 X2 X3 X4 X5 X6 X7 Y

DKI JAKARTA 2016 149831.90 384.30 76.96 75.78 97.01 6.12 30.45 1.19

DKI JAKARTA 2017 157636.60 389.69 85.70 76.99 97.64 7.14 27.05 1.92

DKI JAKARTA 2018 165768.99 373.12 89.04 76.16 98.03 6.65 28.83 2.38

DKI JAKARTA 2019 174812.51 365.55 93.33 78.42 98.12 6.54 29.28 1.04

DKI JAKARTA 2020 170099.68 480.86 93.24 77.57 98.05 10.95 33.80 7.21

JAWA BARAT 2016 26923.51 4224.33 49.43 60.99 97.82 8.89 28.32 1.15

str(data.abhifull)

## 'data.frame': 170 obs. of 10 variables: ## $ Provinsi: chr "ACEH" "ACEH" "ACEH" "ACEH" ... ## $ Tahun : int 2016 2017 2018 2019 2020 2016 2017 2018 2019 2020 ... ## $ X1 : num 22835 23363 24014 24842 25018 ... ## $ X2 : num 848 873 839 819 815 ... ## $ X3 : num 35.2 44.8 56.9 65.2 72 ... ## $ X4 : num 53 54.2 59 57.8 59.6 ... ## $ X5 : num 98.2 98.5 99.1 99.1 99 ... ## $ X6 : num 7.57 6.57 6.34 6.17 6.59 5.84 5.6 5.55 5.39 6.91 ... ## $ X7 : num 25.8 24.9 30.2 29.3 27.1 ... ## $ Y : num 2.26 1.99 5.04 0.79 5.68 2.06 3.22 4.91 0.78 5.54 ...

str(data.abhitinggi)

## 'data.frame': 40 obs. of 10 variables: ## $ Provinsi: chr "DKI JAKARTA" "DKI JAKARTA" "DKI JAKARTA" "DKI JAKARTA" ... ## $ Tahun : int 2016 2017 2018 2019 2020 2016 2017 2018 2019 2020 ... ## $ X1 : num 149832 157637 165769 174813 170100 ... ## $ X2 : num 384 390 373 366 481 ... ## $ X3 : num 77 85.7 89 93.3 93.2 ... ## $ X4 : num 75.8 77 76.2 78.4 77.6 ... ## $ X5 : num 97 97.6 98 98.1 98 ... ## $ X6 : num 6.12 7.14 6.65 6.54 10.95 ... ## $ X7 : num 30.4 27.1 28.8 29.3 33.8 ... ## $ Y : num 1.19 1.92 2.38 1.04 7.21 1.15 1.16 2.12 0.97 7.15 ...

Eksplorasi Data

Line Chart

ggplot(data.abhitinggi, aes(Tahun, Y, colour=Provinsi)) + geom_line()

Boxplot

par(mfrow=c(3,3)) boxplot(data.abhifull$Y, xlab = "ABH") boxplot(data.abhifull$X1, xlab = "PDRB") boxplot(data.abhifull$X2, xlab = "Miskin") boxplot(data.abhifull$X3, xlab = "Akses Internet") boxplot(data.abhifull$X4, xlab = "Telepon") boxplot(data.abhifull$X5, xlab = "APM SD") boxplot(data.abhifull$X4, xlab = "Penganguran") boxplot(data.abhifull$X4, xlab = "Sakit")

Histogram

plot_histogram(data = data.abhifull[, c(-1,-2)], nrow = 3, ncol = 3, geom_histogram_args = list (fill="#D27685"))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Matriks Korelasi (Heatmap)

korelasi <- cor(data.abhifull[, c(-1,-2)]) ggcorrplot(korelasi, type="lower", lab = TRUE)

Geospasial

Visualisasi geospasial dengan bantuan software Google Data Studio dan grafik iteratif dapat diakses melalui tautan sebagai berikut: Geospasial ABH Indonesia

Grafik Geospasial

Plot Means & Coplot

plotmeans(Y ~ Provinsi, main="Keragaman Antar Individu",data.abhifull)

## Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped ## Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped

## Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped ## Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped

plotmeans(Y ~ Tahun, main="Keragaman Antar Waktu", data.abhifull)

## Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped

## Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, : zero- ## length arrow is of indeterminate angle and so skipped

coplot(Y ~ Tahun|Provinsi, type="b", data=data.abhifull)

Gambar di atas menunjukkan bahwa nilai keragaman pada data terlihat terdapat keragaman yang terjadi baik antara individu dan/atau antara waktu.

Pemodelan Regresi

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Tidak terdeteksi gejala multikolinearitas (nilai VIF < 10)
$H_1$ : Terdeteksi gejala multikolineeritas (nilai VIF > 10)

ols <- lm(Y~X1+X2+X3+X4+X5+X6+X7,data=data.abhifull) summary(ols)

## ## Call: ## lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.9136 -2.6333 -0.9656 1.0235 21.9516 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 2.245e+01 1.080e+01 2.079 0.039205 * ## X1 -1.925e-05 1.718e-05 -1.120 0.264246 ## X2 -3.975e-04 3.455e-04 -1.150 0.251695 ## X3 1.511e-01 4.104e-02 3.682 0.000314 *** ## X4 -2.200e-01 9.563e-02 -2.301 0.022692 * ## X5 -8.303e-02 1.399e-01 -0.593 0.553682 ## X6 2.757e-02 2.277e-01 0.121 0.903813 ## X7 -1.929e-01 8.191e-02 -2.355 0.019709 * ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 4.439 on 162 degrees of freedom ## Multiple R-squared: 0.1325, Adjusted R-squared: 0.09499 ## F-statistic: 3.534 on 7 and 162 DF, p-value: 0.001468

Pemeriksaan Multikolinieritas

multicollinearity(ols)

## # Check for Multicollinearity ## ## Low Correlation ## ## Term VIF VIF 95% CI Increased SE Tolerance Tolerance 95% CI ## X1 2.46 [2.00, 3.13] 1.57 0.41 [0.32, 0.50] ## X2 1.26 [1.11, 1.59] 1.12 0.80 [0.63, 0.90] ## X3 3.55 [2.82, 4.57] 1.88 0.28 [0.22, 0.36] ## X6 1.53 [1.31, 1.92] 1.24 0.65 [0.52, 0.76] ## X7 1.78 [1.49, 2.24] 1.33 0.56 [0.45, 0.67] ## X4 5.76 [4.49, 7.51] 2.40 0.17 [0.13, 0.22] ## ## Moderate Correlation ## ## Term VIF VIF 95% CI Increased SE Tolerance Tolerance 95% CI ## X5 2.09 [1.73, 2.65] 1.45 0.48 [0.38, 0.58]

Berdasarkan uji multikolinieritas, dapat disimpulkan bahwa semua model regresi klasik dengan pendekatan OLS tidak terdeteksi adanya gejala multikolinieritas karena nilai VIF < 10 dan nilia tolerance > 0.1 maka terima H0.

Penentuan Model Estimasi

Model Gabungan

cem <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data=data.abhifull, model = "pooling") summary(cem)

## Pooling Model ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## model = "pooling") ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -5.91360 -2.63326 -0.96562 1.02354 21.95161 ## ## Coefficients: ## Estimate Std. Error t-value Pr(>|t|) ## (Intercept) 2.2450e+01 1.0799e+01 2.0789 0.0392046 * ## X1 -1.9249e-05 1.7182e-05 -1.1203 0.2642459 ## X2 -3.9749e-04 3.4554e-04 -1.1504 0.2516948 ## X3 1.5114e-01 4.1043e-02 3.6824 0.0003144 *** ## X4 -2.2000e-01 9.5630e-02 -2.3006 0.0226922 * ## X5 -8.3026e-02 1.3989e-01 -0.5935 0.5536822 ## X6 2.7566e-02 2.2775e-01 0.1210 0.9038134 ## X7 -1.9291e-01 8.1909e-02 -2.3552 0.0197094 * ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 3680 ## Residual Sum of Squares: 3192.5 ## R-Squared: 0.13247 ## Adj. R-Squared: 0.094988 ## F-statistic: 3.53398 on 7 and 162 DF, p-value: 0.0014685

model_performance(cem, metrics = "all")

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## ----------------------------------------------------------------- ## 999.007 | 1000.132 | 1027.229 | 0.132 | 0.095 | 4.333 | 4.439

Berdasarkan persamaan model gabungan di atas, peubah X3, X4, dan X7 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.13247 dan Adj. R-Squared sebesar 0.094988. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 9.5% keragaman peubah respon.

Model Pengaruh Tetap

Model Pengaruh Tetap: Pengaruh Dua Arah

fem.two <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, model = "within", effect= "twoways", index = c("Provinsi","Tahun")) summary(fem.two)

## Twoways effects Within Model ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "twoways", model = "within", index = c("Provinsi", ## "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -7.7208 -1.3892 -0.1254 1.2470 14.9641 ## ## Coefficients: ## Estimate Std. Error t-value Pr(>|t|) ## X1 -0.00045164 0.00013463 -3.3546 0.001053 ** ## X2 -0.00104572 0.00365101 -0.2864 0.775031 ## X3 -0.45078919 0.14802214 -3.0454 0.002834 ** ## X4 0.72769482 0.34871526 2.0868 0.038939 * ## X5 1.38282549 1.65972481 0.8332 0.406341 ## X6 -0.52209977 0.55404018 -0.9424 0.347831 ## X7 -0.04805540 0.17949361 -0.2677 0.789350 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2100.4 ## Residual Sum of Squares: 1721 ## R-Squared: 0.18066 ## Adj. R-Squared: -0.10775 ## F-statistic: 3.93734 on 7 and 125 DF, p-value: 0.00064602

model_performance(fem.two, metrics = "all")

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 891.962 | 892.856 | 917.048 | 0.181 | -0.108 | 3.182 | 3.249

Berdasarkan persamaan model pengaruh tetap dua arah di atas, peubah X1, X3, dan X4 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.18066 dan Adj. R-Squared sebesar 0.10775. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 10.8% keragaman peubah respon.

Model Pengaruh Tetap: Pengaruh Individu

fem.inv <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, model = "within", effect= "individual", index = c("Provinsi","Tahun")) summary(fem.inv)

## Oneway (individual) effect Within Model ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "individual", model = "within", index = c("Provinsi", ## "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -7.52935 -2.49602 -0.28879 1.46908 15.84982 ## ## Coefficients: ## Estimate Std. Error t-value Pr(>|t|) ## X1 -0.00050243 0.00013775 -3.6475 0.0003831 *** ## X2 0.00046444 0.00390009 0.1191 0.9053947 ## X3 -0.11001831 0.07825728 -1.4059 0.1621715 ## X4 0.85027027 0.33582426 2.5319 0.0125450 * ## X5 4.19894195 1.62881845 2.5779 0.0110633 * ## X6 0.50461915 0.48691307 1.0364 0.3019712 ## X7 -0.14612133 0.18526322 -0.7887 0.4317211 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2869.6 ## Residual Sum of Squares: 2118.3 ## R-Squared: 0.26182 ## Adj. R-Squared: 0.032931 ## F-statistic: 6.53642 on 7 and 129 DF, p-value: 1.3177e-06

model_performance(fem.inv, metrics = "all")

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 927.276 | 928.171 | 952.363 | 0.262 | 0.033 | 3.530 | 3.605

Berdasarkan persamaan model pengaruh tetap individu di atas, peubah X1, X4, dan X5 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.262 dan Adj. R-Squared sebesar 0.033. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 3.3% keragaman peubah respon.

Model Pengaruh Tetap: Pengaruh Waktu

fem.wkt <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, model = "within", effect= "time", index = c("Provinsi","Tahun")) summary(fem.wkt)

## Oneway (time) effect Within Model ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "time", model = "within", index = c("Provinsi", ## "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -9.94438 -2.17911 -0.37428 1.40159 21.71972 ## ## Coefficients: ## Estimate Std. Error t-value Pr(>|t|) ## X1 -7.8964e-06 1.5885e-05 -0.4971 0.61980 ## X2 -2.6346e-04 3.1714e-04 -0.8307 0.40737 ## X3 -5.1285e-02 8.5111e-02 -0.6026 0.54767 ## X4 -3.2275e-02 1.1878e-01 -0.2717 0.78620 ## X5 -7.3608e-03 1.2889e-01 -0.0571 0.95453 ## X6 8.2209e-02 2.1352e-01 0.3850 0.70074 ## X7 -1.4700e-01 7.8869e-02 -1.8638 0.06421 . ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2910.7 ## Residual Sum of Squares: 2570.9 ## R-Squared: 0.11676 ## Adj. R-Squared: 0.055263 ## F-statistic: 2.98369 on 7 and 158 DF, p-value: 0.005741

model_performance(fem.wkt, metrics = "all")

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 960.195 | 961.090 | 985.282 | 0.117 | 0.055 | 3.889 | 3.971

Berdasarkan persamaan model pengaruh tetap waktu di atas, tidak ada peubah penjelas berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.11676 dan Adj. R-Squared sebesar 0.055263. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 5.5% keragaman peubah respon.

Perbandingan Model FEM

Berdasarkan perbandingan antara model pengaruh acak dua arah di atas, model pengaruh tetap dua arah terpilih berdasarkan ukuran nilai AIC, AICc, R-Squared, dan Adj. R-Squared untuk dilakukan tahapan analisis selanjutnya.

Model Pengaruh Acak

model <- Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7

Model Pengaruh Acak: Pengaruh Individu

rem.gls.inv <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, index = c("Provinsi","Tahun"), effect = "individual", model = "random",random.method = "nerlove") summary(rem.gls.inv)

## Oneway (individual) effect Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "individual", model = "random", random.method = "nerlove", ## index = c("Provinsi", "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 12.46 3.53 0.024 ## individual 497.00 22.29 0.976 ## theta: 0.9294 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -6.38227 -2.52275 -0.53031 1.46063 15.61696 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) -1.1042e+02 8.9113e+01 -1.2391 0.2153055 ## X1 -3.2037e-04 9.2787e-05 -3.4528 0.0005549 *** ## X2 -1.3314e-04 2.5938e-03 -0.0513 0.9590612 ## X3 -5.3715e-02 6.1631e-02 -0.8716 0.3834470 ## X4 9.6985e-01 2.8088e-01 3.4529 0.0005545 *** ## X5 7.7642e-01 9.6728e-01 0.8027 0.4221587 ## X6 4.7550e-01 4.3703e-01 1.0880 0.2765814 ## X7 -1.9240e-01 1.6484e-01 -1.1672 0.2431308 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2873.7 ## Residual Sum of Squares: 2285.9 ## R-Squared: 0.20455 ## Adj. R-Squared: 0.17018 ## Chisq: 41.6594 on 7 DF, p-value: 6.0472e-07

model_performance(rem.gls.inv)

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 942.218 | 943.343 | 970.440 | 0.205 | 0.170 | 3.667 | 3.756

Berdasarkan persamaan model pengaruh acak individu di atas, peubah X1 dan X4 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.20455 dan Adj. R-Squared sebesar 0.17018. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 17% keragaman peubah respon.

Pemilihan Metode

Pada package R, pendugaan model pengaruh acak disediakan beberapa pilihan metode yang dapat digunakan adalah sebagai berikut.
1.walhus : Wallace and Hussain (1969)
2. swar : Swamy Arora (1972)
3. amemiya : Amemiya (1971)
4. ht : Hausman and Taylor (1981)
5. nerlove : Nerlove (1971)

rem.swar <- update(rem.gls.inv, random.method = "swar") rem.walhus <- update(rem.gls.inv, random.method = "walhus") rem.amemiya <- update(rem.gls.inv, random.method = "amemiya") rem.nerlove <- update(rem.gls.inv, random.method = "nerlove") rem.ht <- update(rem.gls.inv, random.method = "ht")

compare_performance(rem.swar, rem.walhus, rem.amemiya, rem.nerlove, rem.ht, metrics="all")

## Some of the nested models seem to be identical

## # Comparison of Model Performance Indices ## ## Name | Model | AIC (weights) | AICc (weights) | BIC (weights) | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------------------------------------------------- ## rem.swar | plm | 999.0 (<.001) | 1000.1 (<.001) | 1027.2 (<.001) | 0.132 | 0.095 | 4.333 | 4.439 ## rem.walhus | plm | 999.0 (<.001) | 1000.1 (<.001) | 1027.2 (<.001) | 0.132 | 0.095 | 4.333 | 4.439 ## rem.amemiya | plm | 943.7 (0.241) | 944.9 (0.241) | 972.0 (0.241) | 0.198 | 0.163 | 3.683 | 3.773 ## rem.nerlove | plm | 942.2 (0.518) | 943.3 (0.518) | 970.4 (0.518) | 0.205 | 0.170 | 3.667 | 3.756 ## rem.ht | plm | 943.7 (0.241) | 944.9 (0.241) | 972.0 (0.241) | 0.198 | 0.163 | 3.683 | 3.773

Berdasarkan perbandingan metode di atas, metode nerlove terpilih berdasarkan ukuran nili AIC, AICc, R-Squared, dan Adj. R-Squared.

Model Pengaruh Acak: Pengaruh Waktu

rem.gls.wkt <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, index = c("Provinsi","Tahun"), effect = "time", model = "random",random.method = "nerlove") summary(rem.gls.wkt)

## Oneway (time) effect Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "time", model = "random", random.method = "nerlove", ## index = c("Provinsi", "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 15.123 3.889 0.628 ## time 8.972 2.995 0.372 ## theta: 0.7827 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -9.09184 -2.14776 -0.60897 1.02578 22.13454 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) 1.5785e+01 1.0185e+01 1.5499 0.12117 ## X1 -9.8086e-06 1.5786e-05 -0.6214 0.53436 ## X2 -2.9076e-04 3.1564e-04 -0.9212 0.35695 ## X3 -4.2675e-03 7.6715e-02 -0.0556 0.95564 ## X4 -8.0270e-02 1.1217e-01 -0.7156 0.47424 ## X5 -1.9263e-02 1.2827e-01 -0.1502 0.88062 ## X6 6.8841e-02 2.1244e-01 0.3241 0.74590 ## X7 -1.6106e-01 7.7777e-02 -2.0707 0.03838 * ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2947.1 ## Residual Sum of Squares: 2622.1 ## R-Squared: 0.11025 ## Adj. R-Squared: 0.071809 ## Chisq: 20.0745 on 7 DF, p-value: 0.005411

model_performance(rem.gls.wkt)

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 965.551 | 966.676 | 993.773 | 0.110 | 0.072 | 3.927 | 4.023

Berdasarkan persamaan model pengaruh acak waktu di atas, peubah X7 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.11025 dan Adj. R-Squared sebesar 0.071809. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 7.2% keragaman peubah respon.

Model Pengaruh Acak: Pengaruh Dua Arah

rem.gls.two <- plm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, index = c("Provinsi","Tahun"), effect = "twoways", model = "random",random.method = "nerlove") summary(rem.gls.two)

## Twoways effects Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data.abhifull, ## effect = "twoways", model = "random", random.method = "nerlove", ## index = c("Provinsi", "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 10.123 3.182 0.039 ## individual 216.079 14.700 0.830 ## time 34.090 5.839 0.131 ## theta: 0.9037 (id) 0.9069 (time) 0.8777 (total) ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -8.48893 -1.58003 -0.20336 1.06741 15.42534 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) 4.4154e+00 6.7125e+01 0.0658 0.9475539 ## X1 -2.0752e-04 7.4234e-05 -2.7955 0.0051824 ** ## X2 -4.9949e-04 1.9766e-03 -0.2527 0.8005011 ## X3 -3.6343e-01 1.2034e-01 -3.0200 0.0025279 ** ## X4 9.2360e-01 2.6257e-01 3.5175 0.0004356 *** ## X5 -2.1340e-01 7.2793e-01 -0.2932 0.7694028 ## X6 -3.7671e-01 4.7674e-01 -0.7902 0.4294241 ## X7 -9.3815e-02 1.5592e-01 -0.6017 0.5473907 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 2114.6 ## Residual Sum of Squares: 1854.4 ## R-Squared: 0.12307 ## Adj. R-Squared: 0.08518 ## Chisq: 22.7359 on 7 DF, p-value: 0.0018948

rem.walhus1 <- update(rem.gls.two, random.method = "walhus") rem.amemiya1 <- update(rem.gls.two, random.method = "amemiya") rem.nerlove1 <- update(rem.gls.two, random.method = "nerlove")

compare_performance(rem.walhus1, rem.amemiya1, rem.nerlove1, metrics="all")

## Some of the nested models seem to be identical

## # Comparison of Model Performance Indices ## ## Name | Model | AIC (weights) | AICc (weights) | BIC (weights) | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------------------------------------------------- ## rem.walhus1 | plm | 971.6 (<.001) | 972.7 (<.001) | 999.8 (<.001) | 0.105 | 0.066 | 3.998 | 4.095 ## rem.amemiya1 | plm | 909.2 (0.216) | 910.4 (0.216) | 937.5 (0.216) | 0.112 | 0.074 | 3.328 | 3.409 ## rem.nerlove1 | plm | 906.7 (0.784) | 907.8 (0.784) | 934.9 (0.784) | 0.123 | 0.085 | 3.303 | 3.383

model_performance(rem.gls.two)

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 906.653 | 907.778 | 934.875 | 0.123 | 0.085 | 3.303 | 3.383

Berdasarkan persamaan model pengaruh acak dua arah di atas, peubah X1, X3, dan X4 berpengaruh signifikan terhadap Peubah Y pada taraf nyata 5% dengan R-squared sebesar 0.12307 dan Adj. R-Squared sebesar 0.08518. Artinya, peubah penjelas model mampu menjelaskan hanya sebesar 8.5% keragaman peubah respon.

Perbandingan Model FEM

Berdasarkan perbandingan antara model pengaruh acak dua arah di atas, model pengaruh acak dua arah terpilih berdasarkan ukuran nilai AIC, AICc, R-Squared, dan Adj. R-Squared untuk dilakukan tahapan analisis selanjutnya.

Uji Chow

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Model gabungan
$H_1$ : Model pengaruh tetap dua arah
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

pooltest(cem,fem.two)

## ## F statistic ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## F = 2.8887, df1 = 37, df2 = 125, p-value = 6.115e-06 ## alternative hypothesis: unstability

Dengan tingkat keyakinan 95%, kita yakin bahwa model pengaruh tetap dua arah merupakan model yang lebih baik untuk digunakan dibandingkan oleh model gabungan.

Uji Hausman

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Model pengaruh acak dua arah
$H_1$ : Model pengaruh tetap dua arah
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

phtest(fem.two, rem.gls.two)

## ## Hausman Test ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## chisq = 6.7816, df = 7, p-value = 0.452 ## alternative hypothesis: one model is inconsistent

Dengan tingkat keyakinan 95%, kita yakin bahwa model pengaruh acak dua arah merupakan model yang lebih baik untuk digunakan dibandingkan oleh model pengaruh tetap dua arah.

Uji Langrange Multiplier

Pengaruh Individu

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Tidak ada pengaruh individu
$H_1$ : Ada pengaruh individu
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

plmtest(rem.gls.two, type = "bp", effect="individual")

## ## Lagrange Multiplier Test - (Breusch-Pagan) ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## chisq = 1.6738, df = 1, p-value = 0.1958 ## alternative hypothesis: significant effects

Pengaruh Waktu

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Tidak ada pengaruh waktu
$H_1$ : Ada pengaruh waktu
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

plmtest(rem.gls.two, type = "bp", effect="time")

## ## Lagrange Multiplier Test - time effects (Breusch-Pagan) ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## chisq = 52.45, df = 1, p-value = 4.413e-13 ## alternative hypothesis: significant effects

Pengaruh Individu dan Waktu

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Tidak ada pengaruh individu dan waktu
$H_1$ : Ada pengaruh individu dan waktu
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

plmtest(rem.gls.two, type = "bp", effect="twoways" )

## ## Lagrange Multiplier Test - two-ways effects (Breusch-Pagan) ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## chisq = 54.124, df = 2, p-value = 1.767e-12 ## alternative hypothesis: significant effects

Berdasarkan hasil uji pengaruh di atas, model memiliki pengaruh individu dan waktu, maka model yang paling tepat digunakan adalah model pengaruh acak: pengaruh waktu pada taraf nyata 5%.

Uji Asumsi Sisaan

Sebelum Transformasi

1. Normalitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan menyebar normal
$H_1$ : Sisaan tidak menyebar normal
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

res.rem.gls.wkt <- residuals(rem.gls.wkt) jarque.bera.test(res.rem.gls.wkt) # Tolak H0

## ## Jarque Bera Test ## ## data: res.rem.gls.wkt ## X-squared = 777.23, df = 2, p-value < 2.2e-16

#Histogram hist(res.rem.gls.wkt, xlab = "sisaan", col = "light blue", breaks=30, prob = TRUE) lines(density(res.rem.gls.wkt), # density plot lwd = 2, # thickness of line col = "blue")

#Plotqqnorm set.seed(1353) res.rem.gls.wkt1 <- as.numeric(res.rem.gls.wkt) qqnorm(res.rem.gls.wkt1,datax=T, col="blue") qqline(rnorm(length(res.rem.gls.wkt1),mean(res.rem.gls.wkt1),sd(res.rem.gls.wkt1)),datax=T, col="red")

Berdasarkan uji formal dan secara eksploratif di atas, terlihat bahwa sebaran data pada model pengaruh acak: pengaruh waktu sisaan tidak mengikuti sebaran normal.

Autokorelasi

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan saling bebas
$H_1$ : Sisaan tidak saling bebas
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

# Durbin Watson pdwtest(rem.gls.wkt) # Tak tolak H0

## ## Durbin-Watson test for serial correlation in panel models ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## DW = 1.8851, p-value = 0.118 ## alternative hypothesis: serial correlation in idiosyncratic errors

# Durbin Watson pdwtest(rem.gls.wkt) # Tak tolak H0

## ## Durbin-Watson test for serial correlation in panel models ## ## data: Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## DW = 1.8851, p-value = 0.118 ## alternative hypothesis: serial correlation in idiosyncratic errors

# Augmented DF adf.test(res.rem.gls.wkt, k=2) # Tolak H0

## Warning in adf.test(res.rem.gls.wkt, k = 2): p-value smaller than printed p- ## value

## ## Augmented Dickey-Fuller Test ## ## data: res.rem.gls.wkt ## Dickey-Fuller = -7.134, Lag order = 2, p-value = 0.01 ## alternative hypothesis: stationary

Berdasarkan uji formal di atas, terlihat bahwa sebaran data pada model pengaruh acak: pengaruh waktu sisaan tidak saling bebas berdasarkan Uji ADF.

Homogenitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan memiliki ragam homogen
$H_1$ : Sisaan tidak memiliki ragam homogen
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

bptest(rem.gls.wkt)

## ## studentized Breusch-Pagan test ## ## data: rem.gls.wkt ## BP = 19.236, df = 7, p-value = 0.00748

Berdasarkan uji formal di atas, terlihat bahwa sebaran data pada model pengaruh acak: pengaruh waktu sisaan tidak memiliki ragam homogen.

Penanganan Kondisi Tidak Standar

Transformasi Logaritma Natural

rem.gls.wkt2 <- plm(log(Y)~log(X1)+log(X2)+log(X3)+log(X5)+log(X6)+log(X7), data.abhifull, index = c("Provinsi","Tahun"), effect = "time", model = "random",random.method = "nerlove") summary(rem.gls.wkt2)

## Oneway (time) effect Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = log(Y) ~ log(X1) + log(X2) + log(X3) + log(X5) + ## log(X6) + log(X7), data = data.abhifull, effect = "time", ## model = "random", random.method = "nerlove", index = c("Provinsi", ## "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 0.7248 0.8513 0.679 ## time 0.3433 0.5859 0.321 ## theta: 0.7582 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -2.649621 -0.468945 0.038415 0.561874 2.143141 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) 22.163781 11.863190 1.8683 0.061723 . ## log(X1) -0.392400 0.190947 -2.0550 0.039877 * ## log(X2) -0.036411 0.074403 -0.4894 0.624573 ## log(X3) 0.650257 0.602251 1.0797 0.280271 ## log(X5) -3.449608 2.727809 -1.2646 0.206012 ## log(X6) 0.212650 0.237561 0.8951 0.370713 ## log(X7) -1.269112 0.473924 -2.6779 0.007409 ** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 140.99 ## Residual Sum of Squares: 125.89 ## R-Squared: 0.1071 ## Adj. R-Squared: 0.074234 ## Chisq: 19.5515 on 6 DF, p-value: 0.0033269

Transformasi Box-Cox

Ketika nilai lamda mendekati nol maka transformasi Box-cox mirip dengan transformasi logaritma natural yang dapat diterapkan pada peubah angka buta huruf (Y) (Montgomery et al. 2015).

y.trans <- caret::BoxCoxTrans(data.abhifull$Y) data2 <- cbind(data.abhifull, y.new = predict(y.trans, data.abhifull$Y)) print(y.trans) # Y^lamda = Y^0 = log(Y)

## Box-Cox Transformation ## ## 170 data points used to estimate Lambda ## ## Input data summary: ## Min. 1st Qu. Median Mean 3rd Qu. Max. ## 0.130 1.173 2.075 3.976 6.223 28.980 ## ## Largest/Smallest: 223 ## Sample Skewness: 2.78 ## ## Estimated Lambda: 0 ## With fudge factor, Lambda = 0 will be used for transformations

rem.gls.wkt3 <- plm(y.new ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data2, index = c("Provinsi","Tahun"), effect = "time", model = "random",random.method = "nerlove") summary(rem.gls.wkt3)

## Oneway (time) effect Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = y.new ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data2, ## effect = "time", model = "random", random.method = "nerlove", ## index = c("Provinsi", "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 0.7325 0.8559 0.676 ## time 0.3506 0.5921 0.324 ## theta: 0.7594 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -2.520405 -0.465669 -0.012494 0.485558 2.123710 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) 3.9357e+00 2.2362e+00 1.7600 0.0784 . ## X1 -3.1387e-06 3.4749e-06 -0.9032 0.3664 ## X2 -6.3002e-05 6.9501e-05 -0.9065 0.3647 ## X3 1.6069e-02 1.6562e-02 0.9703 0.3319 ## X4 -2.8438e-02 2.4456e-02 -1.1628 0.2449 ## X5 -1.3834e-02 2.8242e-02 -0.4898 0.6242 ## X6 5.6388e-02 4.6770e-02 1.2056 0.2280 ## X7 -3.9584e-02 1.7097e-02 -2.3152 0.0206 * ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 140.96 ## Residual Sum of Squares: 127.23 ## R-Squared: 0.097451 ## Adj. R-Squared: 0.058452 ## Chisq: 17.4916 on 7 DF, p-value: 0.014487

rem.gls.wkt4 <- plm(log(Y) ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data2, index = c("Provinsi","Tahun"), effect = "time", model = "random",random.method = "nerlove") summary(rem.gls.wkt4)

## Oneway (time) effect Random Effect Model ## (Nerlove's transformation) ## ## Call: ## plm(formula = log(Y) ~ X1 + X2 + X3 + X4 + X5 + X6 + X7, data = data2, ## effect = "time", model = "random", random.method = "nerlove", ## index = c("Provinsi", "Tahun")) ## ## Balanced Panel: n = 34, T = 5, N = 170 ## ## Effects: ## var std.dev share ## idiosyncratic 0.7325 0.8559 0.676 ## time 0.3506 0.5921 0.324 ## theta: 0.7594 ## ## Residuals: ## Min. 1st Qu. Median 3rd Qu. Max. ## -2.520405 -0.465669 -0.012494 0.485558 2.123710 ## ## Coefficients: ## Estimate Std. Error z-value Pr(>|z|) ## (Intercept) 3.9357e+00 2.2362e+00 1.7600 0.0784 . ## X1 -3.1387e-06 3.4749e-06 -0.9032 0.3664 ## X2 -6.3002e-05 6.9501e-05 -0.9065 0.3647 ## X3 1.6069e-02 1.6562e-02 0.9703 0.3319 ## X4 -2.8438e-02 2.4456e-02 -1.1628 0.2449 ## X5 -1.3834e-02 2.8242e-02 -0.4898 0.6242 ## X6 5.6388e-02 4.6770e-02 1.2056 0.2280 ## X7 -3.9584e-02 1.7097e-02 -2.3152 0.0206 * ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Total Sum of Squares: 140.96 ## Residual Sum of Squares: 127.23 ## R-Squared: 0.097451 ## Adj. R-Squared: 0.058452 ## Chisq: 17.4916 on 7 DF, p-value: 0.014487

Uji Asumsi Sisaan Setelah Transformasi ln

1. Normalitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan menyebar normal
$H_1$ : Sisaan tidak menyebar normal
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

res.rem.gls.wkt2 <- residuals(rem.gls.wkt2) jarque.bera.test(res.rem.gls.wkt2) #Tak tolak H0

## ## Jarque Bera Test ## ## data: res.rem.gls.wkt2 ## X-squared = 5.6596, df = 2, p-value = 0.05902

#Histogram hist(res.rem.gls.wkt2, xlab = "sisaan", col = "light blue", breaks=30, prob = TRUE) lines(density(res.rem.gls.wkt2), # density plot lwd = 2, # thickness of line col = "blue")

#Plotqqnorm set.seed(1353) res.rem.gls.wkt22 <- as.numeric(res.rem.gls.wkt2) qqnorm(res.rem.gls.wkt22,datax=T, col="blue") qqline(rnorm(length(res.rem.gls.wkt22),mean(res.rem.gls.wkt22),sd(res.rem.gls.wkt22)),datax=T, col="red")

Berdasarkan histogram dan Q-Q plot di atas, terlihat bahwa sebaran data pada model pengaruh acak waktu sisaan mengikuti sebaran normal.

2. Autokorelasi

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan saling bebas
$H_1$ : Sisaan tidak saling bebas
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

# Durbin Watson pdwtest(rem.gls.wkt2) # Tak tolak H0

## ## Durbin-Watson test for serial correlation in panel models ## ## data: log(Y) ~ log(X1) + log(X2) + log(X3) + log(X5) + log(X6) + log(X7) ## DW = 1.9776, p-value = 0.3035 ## alternative hypothesis: serial correlation in idiosyncratic errors

adf.test(res.rem.gls.wkt22, k=2) #Tolak H0

## Warning in adf.test(res.rem.gls.wkt22, k = 2): p-value smaller than printed p- ## value

## ## Augmented Dickey-Fuller Test ## ## data: res.rem.gls.wkt22 ## Dickey-Fuller = -6.2527, Lag order = 2, p-value = 0.01 ## alternative hypothesis: stationary

2. Homogenitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan memiliki ragam homogen
$H_1$ : Sisaan tidak memiliki ragam homogen
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

bptest(rem.gls.wkt2) # Tolak H0

## ## studentized Breusch-Pagan test ## ## data: rem.gls.wkt2 ## BP = 13.62, df = 6, p-value = 0.03418

Uji Asumsi Sisaan Setelah Transformasi Box-Cox

1. Normalitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan menyebar normal
$H_1$ : Sisaan tidak menyebar normal
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

res.rem.gls.wkt3 <- residuals(rem.gls.wkt3) jarque.bera.test(res.rem.gls.wkt3) # Tak Tolak H0

## ## Jarque Bera Test ## ## data: res.rem.gls.wkt3 ## X-squared = 5.1278, df = 2, p-value = 0.077

#Histogram hist(res.rem.gls.wkt3, xlab = "sisaan", col = "light blue", breaks=30, prob = TRUE) lines(density(res.rem.gls.wkt3), # density plot lwd = 2, # thickness of line col = "blue")

#Plotqqnorm set.seed(1353) res.rem.gls.wkt33 <- as.numeric(res.rem.gls.wkt3) qqnorm(res.rem.gls.wkt33,datax=T, col="blue") qqline(rnorm(length(res.rem.gls.wkt33),mean(res.rem.gls.wkt33),sd(res.rem.gls.wkt33)),datax=T, col="red")

Berdasarkan histogram dan Q-Q plot di atas, terlihat bahwa sebaran data pada model pengaruh acak waktu sisaan mengikuti sebaran normal.

2. Autokorelasi

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan saling bebas
$H_1$ : Sisaan tidak saling bebas
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

# Durbin Watson pdwtest(rem.gls.wkt3) # Tak tolak H0

## ## Durbin-Watson test for serial correlation in panel models ## ## data: y.new ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 ## DW = 1.9812, p-value = 0.2906 ## alternative hypothesis: serial correlation in idiosyncratic errors

adf.test(res.rem.gls.wkt33, k=2) # Tolak H0

## Warning in adf.test(res.rem.gls.wkt33, k = 2): p-value smaller than printed p- ## value

## ## Augmented Dickey-Fuller Test ## ## data: res.rem.gls.wkt33 ## Dickey-Fuller = -6.4739, Lag order = 2, p-value = 0.01 ## alternative hypothesis: stationary

2. Homogenitas

Hipotesis yang diuji adalah sebagai berikut.
$H_0$ : Sisaan memiliki ragam homogen
$H_1$ : Sisaan tidak memiliki ragam homogen
Tolak H0 jika P-value < α. Nilai α yang digunakan sebesar 5%.

bptest(rem.gls.wkt3) #Tak Tolak H0

## ## studentized Breusch-Pagan test ## ## data: rem.gls.wkt3 ## BP = 13.837, df = 7, p-value = 0.05415

Hasil Uji Formal Asumsi Sisaan

Ukuran Kebaikan Model

compare_performance(rem.gls.wkt, rem.gls.wkt2, rem.gls.wkt3, rem.gls.wkt4, metrics="all")

## Warning: When comparing models, please note that probably not all models were fit ## from same data.

## # Comparison of Model Performance Indices ## ## Name | Model | AIC (weights) | AICc (weights) | BIC (weights) | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------------------------------------------------- ## rem.gls.wkt | plm | 965.6 (<.001) | 966.7 (<.001) | 993.8 (<.001) | 0.110 | 0.072 | 3.927 | 4.023 ## rem.gls.wkt2 | plm | 739.1 (<.001) | 740.0 (<.001) | 764.2 (<.001) | 0.107 | 0.074 | 0.861 | 0.879 ## rem.gls.wkt3 | plm | 451.2 (>.999) | 452.3 (>.999) | 479.4 (>.999) | 0.097 | 0.058 | 0.865 | 0.886 ## rem.gls.wkt4 | plm | 742.9 (<.001) | 744.0 (<.001) | 771.1 (<.001) | 0.097 | 0.058 | 0.865 | 0.886

performance(rem.gls.wkt2)

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 739.095 | 739.990 | 764.182 | 0.107 | 0.074 | 0.861 | 0.879

performance(rem.gls.wkt3)

## # Indices of model performance ## ## AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma ## --------------------------------------------------------------- ## 451.169 | 452.294 | 479.391 | 0.097 | 0.058 | 0.865 | 0.886

Hasil Ukuran Kebaikan Kandidat Model

Interpretasi Model

Model Pengaruh Acak Waktu Setelah Transformasi Box-Cox

Analisis yang telah dilakukan menunjukkan bahwa Model Pengaruh Acak Waktu Setelah Transformasi Box-Cox dengan nilai lamda = 0 merupakan model terbaik dibandingkan model-model data panel lainnya. Persamaan model tersebut dinyatakan sebagai berikut.

Simpulan

Berdasarkan analisis data panel yang telah dilakukan menunjukkan bahwa Model Pengaruh Acak Waktu Setelah Transformasi Box-Cox dengan nilai lamda=0 merupakan model terbaik. Peubah penjelas yang signifikan, yaitu angka kesakitan (X7). Selain itu, model regesi ini merupakan model layak karena memenuhi asumsi normalitas, non autokorelasi, dan homoskedastisitas dengan hasil akurasi model tergolong cukup baik.

Daftar Pustaka

Akinwande MO, Dikko HG, Samson A. 2015. Variance inflation factor: as a condition for the inclusion of suppressor variable(s) in regression analysis. Open Journal of Statistics 5: 754-767.

Baltagi BH. 2011. Econometrics. Ed ke-5. New York(US) : Springer.

Ghozali I dan Ratmono D. 2013. Analisis Multivariat dan Ekonometrika, Teori, Konsep, dan Aplikasi dengan Eviews 8. Semarang: Badan Penerbit Universitas Diponegoro.

Greene WH. 2003. Econometric Analysis. 5th ed. New Jersey (NJ): Prentice Hall International.

Gujarati DN dan Porter DC. 2009. Basic Econometric 5th Edition. New York(US): The McGraw-Hil Companies.

Jaya GNM, Sunengsih N. Tahun terbit. Kajian Analisis Regresi dengan data Panel. Seminar Nasional Penelitian, Pendidikan dan Penerapan MIPA; 2009 Mei 16; Yogyakarta, Indonesia. Yogyakarta: Fakultas MIPA, Universitas Negeri Yogyakarta,51-58.

Nandita DA, Alamsyah LB, Jati EP, dan Widodo E. 2019. Regresi data panel untuk mengetahui faktor-faktor yang mempengaruhi PDRB di Provinsi DIY tahun 2011-2015. Indonesian Journal of Applied Statistics. 2(1): 42-52.

Montgomery DC, Jennings CL, Kulahci M. 2015. Introduction to Time Series Analysis and Forecasting. New Jersey(US):John Wiley & Sons.

Susanti S. 2013. Pengaruh produk domestik regional bruto, pengangguran dan indeks pembangunan manusia terhadap kemiskinan di Jawa Barat dengan menggunakan analisis data panel. Jurnal Matematika Integratif, ISSN, pp.1412-6184.

Provinsi	Tahun	X1	X2	X3	X4	X5	X6	X7	Y
ACEH	2016	22835.29	848.44	35.23	52.97	98.16	7.57	25.78	2.26
ACEH	2017	23362.90	872.61	44.83	54.20	98.54	6.57	24.85	1.99
ACEH	2018	24013.79	839.49	56.89	59.05	99.10	6.34	30.18	5.04
ACEH	2019	24842.30	819.44	65.16	57.75	99.12	6.17	29.33	0.79
ACEH	2020	25018.28	814.91	71.97	59.60	99.03	6.59	27.12	5.68
SUMATERA UTARA	2016	32885.09	1455.95	40.44	54.28	96.57	5.84	22.88	2.06

Provinsi	Tahun	X1	X2	X3	X4	X5	X6	X7	Y
DKI JAKARTA	2016	149831.90	384.30	76.96	75.78	97.01	6.12	30.45	1.19
DKI JAKARTA	2017	157636.60	389.69	85.70	76.99	97.64	7.14	27.05	1.92
DKI JAKARTA	2018	165768.99	373.12	89.04	76.16	98.03	6.65	28.83	2.38
DKI JAKARTA	2019	174812.51	365.55	93.33	78.42	98.12	6.54	29.28	1.04
DKI JAKARTA	2020	170099.68	480.86	93.24	77.57	98.05	10.95	33.80	7.21
JAWA BARAT	2016	26923.51	4224.33	49.43	60.99	97.82	8.89	28.32	1.15