Analisis Regresi Data Panel dalam Pemodelan Angka Buta Huruf di Indonesia Tahun 2016-2020

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).

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.

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.

Library

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.abhifull <- read.csv("C:/SEM6/ADP/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("C:/SEM6/ADP/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

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

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")

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`.

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

Pemodelan Regresi

Hipotesis yang diuji adalah sebagai berikut. H0 : Tidak terdeteksi gejala multikolinearitas (nilai VIF < 10) H1 : 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.

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 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.

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. H0 : Model gabungan H1 : 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. H0 : Model pengaruh acak dua arah H1 : 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. H0 : Tidak ada pengaruh individu H1 : 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. H0 : Tidak ada pengaruh waktu H1 : 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. H0 : Tidak ada pengaruh individu dan waktu H1 : 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. H0 : Sisaan menyebar normal H1 : 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. H0 : Sisaan saling bebas H1 : 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

# 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. H0 : Sisaan memiliki ragam homogen H1 : 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. H0 : Sisaan menyebar normal H1 : 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. H0 : Sisaan saling bebas H1 : 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. H0 : Sisaan memiliki ragam homogen H1 : 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. H0 : Sisaan menyebar normal H1 : 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. H0 : Sisaan saling bebas H1 : 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. H0 : Sisaan memiliki ragam homogen H1 : 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

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) #ln

## # 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) #boxcox

## # 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

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.

log(ABH) = 3.936 - 0.039(Angka Kesakitan_X7) + (pengaruh spesifik unit waktu ke-t)

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

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.

ADPPPP

Dhea Puspita Adinda

2023-06-08

Tinjauan Pustaka

Model Gabungan

Model Pengaruh Tetap

Model Pengaruh Acak

Library

Input Data

Eksplorasi Data

Pemodelan Regresi

Pemeriksaan Multikolinieritas

Penentuan Model Estimasi

Model Gabungan

Model Pengaruh Tetap

Model Pengaruh Tetap: Pengaruh Dua Arah

Model Pengaruh Tetap: Pengaruh Individu

Model Pengaruh Tetap: Pengaruh Waktu

Model Pengaruh Acak

Model Pengaruh Acak: Pengaruh Individu

Pemilihan Metode

Model Pengaruh Acak: Pengaruh Waktu

Model Pengaruh Acak: Pengaruh Dua Arah

Uji Chow

Uji Hausman

Uji Langrange Multiplier

Pengaruh Individu

Pengaruh Waktu

Pengaruh Individu dan Waktu

Uji Asumsi Sisaan

Sebelum Transformasi

1. Normalitas

Autokorelasi

Homogenitas

Penanganan Kondisi Tidak Standar

Transformasi Logaritma Natural

Transformasi Box-Cox

Uji Asumsi Sisaan Setelah Transformasi ln

1. Normalitas

2. Autokorelasi

2. Homogenitas

Uji Asumsi Sisaan Setelah Transformasi Box-Cox

1. Normalitas

2. Autokorelasi

2. Homogenitas

Ukuran Kebaikan Model

Interpretasi Model

Simpulan

Daftar Pustaka