Analisis Regresi Data Panel dalam Pemodelan Persentase Tingkat Pengangguran di Pulau Jawa Tahun 2015-2018

Pendahuluan

Latar Belakang

Pulau Jawa adalah salah satu pulau terpadat di Indonesia dengan tingkat urbanisasi yang tinggi dan aktivitas ekonomi yang maju. Namun, di balik kemajuan tersebut, masih ada permasalahan serius yang perlu diatasi, yaitu tingkat pengangguran terbuka. Tingkat pengangguran terbuka mengacu pada proporsi individu yang siap bekerja, aktif mencari pekerjaan, tetapi tidak dapat menemukan pekerjaan yang sesuai.

Penelitian ini bertujuan untuk menganalisis tingkat pengangguran terbuka di pulau Jawa, dengan menggunakan metode analisis deskriptif dan regresi data panel. Data yang akan digunakan adalah data tingkat pengangguran terbuka dari Badan Pusat Statistik (BPS) dan data faktor-faktor yang mempengaruhinya seperti besaran serapan investasi suatu wilayah, pendidikan, keterampilan, dan informasi demografi di suatu wilayah seperti proporsi jenis kelamin dan satatus perkawinan.

Hasil penelitian ini diharapkan dapat memberikan pemahaman yang lebih baik tentang tingkat pengangguran terbuka di pulau Jawa, faktor-faktor yang mempengaruhinya, dan implikasinya terhadap kebijakan ketenagakerjaan. Penelitian ini diharapkan dapat menjadi dasar bagi pemerintah dan pemangku kepentingan lainnya dalam mengembangkan kebijakan yang efektif untuk mengurangi tingkat pengangguran terbuka dan menciptakan kesempatan kerja yang lebih baik bagi masyarakat.


Tentang Data

Pada penelitian kali ini akan menggunakan beberapa dataset berikut:

  1. Tingkat Penganguran Terbuka: Tingkat Pengangguran Terbuka Menurut Provinsi (Persen), 2015-2023

  2. Persentase Satatus Perkawinan di Suatu Wilayah: Persentase Penduduk Berumur 10 Tahun ke Atas menurut Provinsi, Jenis Kelamin, dan Status Perkawinan, 2009-2018

  3. Pendidikan: Tingkat Penyelesaian Pendidikan Menurut Jenjang Pendidikan dan Provinsi 2015-2022

  4. Persentase Jenis Kelamin: Persentase Penduduk menurut Provinsi dan Jenis Kelamin (Persen), 2015-2018

  5. Keterampilan TIK: Proporsi Remaja Dan Dewasa Usia 15-59 Tahun Dengan Keterampilan Teknologi Informasi Dan Komputer (TIK) Menurut Provinsi (Persen), 2015-2022

  6. Tingkat Investasi di Suatu Wilayah: Realisasi Investasi Penanaman Modal Dalam Negeri Menurut Provinsi (Investasi) (Milyar Rupiah), 2015-2022

Dikarenakan terdapat perbedaan rentang waktu antar data, maka pada penelitian kali ini akan diambil rentang waktu yang dimiliki oleh keseluruhan data yakni tahun 2015 s.d 2018.

Dengan kondisi data yang diperoleh dari web resmi BPS memiliki tampilan multiindeks dan dataset yang terpisah berdasarkan jenjang waktu, maka sebelum digunakan pada penelitian ini terlebih dahulu telah dilakukan penggabungan dataset, penyesuaian indeks data dan pemilihan kolom tertentu yang akan digunakan menggunakan MS.Excel.


Teori Dasar

Data Panel

Salah satu aspek yang harus diperhatikan dalam melakukan pemodelan statistika adalah jenis data yang digunakan. Salah satu jenis data berdasarkan waktu pengumpulannya adalah data panel. Data panel merupakan gabungan data cross section dan time series (runtun/ deret waktu). Dengan kata lain, data panel merupakan data dari beberapa individu sama yang diamati dalam kurun waktu tertentu. Terkadang data panel disebut juga data longitudinal.

Jika kita memiliki T periode waktu (t = 1, 2, …, T) dan N jumlah individu (i = 1, 2, …, N), maka dengan data panel kita akan memiliki total unit observasi sebanyak NT.

Jika jumlah unit waktu sama untuk setiap individu, maka data disebut balanced panel. Jika sebaliknya, yakni jumlah unit waktu berbeda untuk setiap individu, maka disebut unbalanced panel

Ilustrasi Balanced & Unbalanced Data Panel


Pemodelan

Terdapat tiga pendekatan yang digunakan dalam model panel yaitu Common/ Polled Effects, Fixed Effects dan Random Effects.

1. Model Common/ Pooled Effects (CEM)

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. Untuk membuat model CEM di R dapat menggunakan fungsi plm() dari package plm dengan parameter sebagai berikut:

  • formula = Target ~ Prediktor
  • data = berupa dataframe
  • index = c(“kolom_individu”,“kolom_waktu”)
  • model = “pooling”

2. Fixed Effects

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. Untuk membuat model CEM di R dapat menggunakan fungsi plm() dari package plm dengan parameter sebagai berikut:

  • formula = Target ~ Prediktor
  • data = berupa dataframe
  • index = c(“kolom_individu”,“kolom_waktu”)
  • model = “within”

3. Random Effects

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 (generalized 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. Untuk membuat model CEM di R dapat menggunakan fungsi plm() dari package plm dengan parameter sebagai berikut:

  • formula = Target ~ Prediktor
  • data = berupa dataframe
  • index = c(“kolom_individu”,“kolom_waktu”)
  • model = “random”

Pemilihan Model Terbaik

Untuk memilih model mana yang terbaik dari ketiga model diatas, perlu dilakukan pengujian hipotesis berikut ini:

Uji Chow

Menurut Ghozali dan Ratmono (2013), uji chow digunakan untuk memilih pendekatan yang lebih baik antara model gabungan dengan model pengaruh tetap. Untuk melakukan uji Chow dapat menggunakan fungsi pooltest(model_cem, model_fem), dengan hipotesis yang diujikan adalah sebagai berikut.

  • H0 : Model gabungan
  • H1 : Model pengaruh tetap

Keputusan tolak H0 (model pengaruh tetap terpilih) apabila nilai p-value < α.


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). Untuk melakukan uji Chow di R dapat menggunakan fungsi phtest(model_rem, model_fem), dengan Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Model pengaruh acak
  • H1 : Model pengaruh tetap

Keputusan tolak H0 (model pengaruh tetap terpilih) apabila nilai p-value < α.


Pengujian Lanjutan

Apabila model gabungan atau model pengaruh acak yang terpilih sebagai model terbaik, perlu dilakukan uji lanjutan untuk memeriksa apakah model dipengaruhi oleh individu/waktu/individu & waktu. Untuk memeriksa pengaruh tersebut digunakan Uji Lagrange Multiplier.

Di R untuk melakukan Uji Lagrange Multiplier dapat dilakukan dengan menggunakan fungsi plmtest() dari package plm dengan parameter:

  • x = model terpilih
  • type = "bp" , menggunakan Breusch Pagan test
  • effect =
    • “time” : untuk menguji pengaruh waktu
    • “individual” : untuk menguji pengaruh individu
    • “twoways” : untuk menguji pengaruh individu & waktu

berikut hipotesis pengujian untuk pemeriksaan setiap pengaruh:

Pemeriksaan Pengaruh Individu & Waktu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh individu & waktu
  • H1 : Ada pengaruh individu & waktu

Pemeriksaan Pengaruh Individu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh individu
  • H1 : Ada pengaruh individu

Pemeriksaan Pengaruh Waktu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh waktu
  • H1 : Ada pengaruh waktu

dengan ketentuan untuk ketiga hipotesis adalah H0 ditolak jika P-value < α. Dimana nilai α yang umum digunakan adalah sebesar 5%.


Pengujian Asumsi

dikarenakan analisis data panel menggunakan konsep regresi dan time series maka ada beberapa Asumsi yang perlu dipenuhi sebagai berikut:

Asumsi model linear regression:

  1. No Multicollinearity (VIF)
  2. Normality of Residuals
  3. Homogenitas of Residuals

Asumsi Time Series:

  1. No Autocorelation

Note: Pengujian Asumsi hanya perlu dilakukan untuk model yang akan digunakan

1. Pemeriksaan Multikolinieritas

Multicollinearity adalah kondisi adanya korelasi antar prediktor yang kuat. Hal ini tidak diinginkan karena menandakan prediktor redundan pada model, yang seharusnya dapat dipilih salah satu saja dari variable yang hubungannya amat kuat tersebut. Harapannya tidak terjadi multicollinearity. Pemeriksaan Multikolinieritas bisa dilakukan dengan menggunakan fungsi fungsi vif(), dengan ketentuan

  • nilai VIF > 10: terjadi multicollinearity pada model
  • nilai VIF < 10: tidak terjadi multicollinearity pada model

Note: Pemeriksaan asumsi multikolinieritas bisa dilakukan diawal sebelum dilakukan pemodelan data panel, dengan cara terlebih dahulu dilakukan pembuatan model regresi dengan fungsi lm() dan akan diujikan menggunakan fungsi vif().

2. Pengujian Asumsi Normalitas Residual

Model diharapkan menghasilkan error yang berdistribusi normal. Dengan begitu, error lebih banyak berkumpul di sekitar angka nol. Pengujian asumsi normalitas di R dapat dilakukan dengan menggunakan fungsi shapiro.test(). dengan hipotesis yang diujikan adalah sebagai berikut.

  • H0 : Sisaan berdistribusi normal
  • H1 : Sisaan tidak berdistribusi normal

dengan ketentuan H0 ditolak jika P-value < α. Dimana nilai α yang umum digunakan adalah sebesar 5%.

3. Pengujian Asumsi Homogenitas Residual

Diharapkan error yang dihasilkan oleh model menyebar secara acak atau dapat dikatakan variasi konstan. Pengujian asumsi Homogenitas di R dapat dilakukan dengan menggunakan fungsi bptest(), dengan hipotesis yang diujikan adalah sebagai berikut.

  • H0 : Sisaan memiliki ragam homogen
  • H1 : Sisaan tidak memiliki ragam homogen

dengan ketentuan H0 ditolak jika P-value < α. Dimana nilai α yang umum digunakan adalah sebesar 5%.

4. Pengujian Asumsi Autokorelasi

Untuk mengecek ada/tidaknya bisa menggunakan uji Ljung-box dengan menggunakan fungsi Box.test(residual model, type = "Ljung-Box), dengan hipotesis yang diujikan adalah sebagai berikut.

  • H0 : tidak terjadi autokorelasi pada sisaan
  • H1 : terjadi autokorelasi pada sisaan

dengan ketentuan H0 ditolak jika P-value < α. Dimana nilai α yang umum digunakan adalah sebesar 5%.


Alur Kerja Data Panel

Analisis data panel di R dapat dilakukan dengan mengikuti alur kerja berikut ini:
Alur Kerja Analisis Data Panel di R

Analisis Data

Persiapan Data

Load packages

Berikut adalah beberapa Packages yang akan digunakan untuk analisis data kali ini.

#Packages untuk pengolahan dataframe
library(dplyr)
library(tidyr)
library(lubridate)

#Packages untuk membuat visualisasi
library(ggcorrplot)
library(ggplot2)
library(plotly)

#Packages untuk melakukan analisis
library(plm)
library(lfe)
library(lmtest)
library(car)
library(tseries)
library(MLmetrics)

Data Preparation

Dikarenakan setiap variabel yang akan digunakan dalam file csv yang terpisah, maka pada tahapan Read data kali ini akan dilakukan beberapa proses untuk setiap variabel yang akan digunakan, yakni:

  1. Read data setiap variabel
  2. Pengubahan nama kolom
  3. Mengambil data provinsi di pulau Jawa
  4. Melakukan pengubahan struktur data (baris menjadi kolom) untuk memperoleh informasi keterangan waktu
  5. Mengambil data yang diperlukan
  6. Mengulangi langkah 1 s.d 5 untuk setiap variabel
  7. Menggabungkan setiap variabel menjadi sebuah dataframe
# ---- Data Tingkat Pengangguran ----
pengangguran <- read.csv("data/penganguran.csv", sep = ";") %>%
  # mengubah nama kolom
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  
  #  mengambil data provinsi yang berada di pulau Jawa
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  
  # melakukan transformasi untuk mengambil informasi tahun per wilayah
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "Y"
  ) %>%
  select(Provinsi, Tahun, Y)

# ---- Data Persentase Perkawinan ----
perkawinan <- read.csv("data/perkawinan.csv", sep = ";") %>%
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "X1"
  ) %>%
  select(X1)

# ---- Data Persentase Lulus SMA ----
pendidikan <- read.csv("data/persenlulusSMA.csv", sep = ";") %>%
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "X2"
  ) %>%
  select(X2)

# ---- Data Persentase Pria ----
gender <- read.csv("data/Persentase_pria.csv", sep = ";") %>%
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "X3"
  ) %>%
  select(X3)

# ---- Data Persentase Kemampuan TIK ----
TIK <- read.csv("data/TIK.csv", sep = ";") %>%
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "X4"
  ) %>%
  select(X4)

# ---- Data Tingkat Investasi Wilayah ----
investasi <- read.csv("data/investasi.csv", sep = ";") %>%
  mutate(
    "2015" = X2015,
    "2016" = X2016,
    "2017" = X2017,
    "2018" = X2018
  ) %>%
  filter(
    Provinsi %in% "BANTEN" |
      Provinsi %in% "DKI JAKARTA" |
      Provinsi %in% "JAWA BARAT" |
      Provinsi %in% "JAWA TENGAH" |
      Provinsi %in% "DI YOGYAKARTA" |
      Provinsi %in% "JAWA TIMUR"
  ) %>%
  pivot_longer(
    cols = c("2015", "2016", "2017", "2018"),
    names_to = "Tahun",
    values_to = "X5"
  ) %>%
  select(X5)

menggabungkan setiap variabel menjadi sebuah dataframe dan disimpan kedalam objek bernama data_tpt.

data_tpt <- cbind(pengangguran,perkawinan,
                  pendidikan,gender,TIK, 
                  investasi) %>% 
  mutate(Tahun = as.factor(Tahun),
         Provinsi = as.factor(Provinsi))

Pemeriksaan Struktur Data

glimpse(data_tpt)
#> Rows: 24
#> Columns: 8
#> $ Provinsi <fct> DKI JAKARTA, DKI JAKARTA, DKI JAKARTA, DKI JAKARTA, JAWA BARA…
#> $ Tahun    <fct> 2015, 2016, 2017, 2018, 2015, 2016, 2017, 2018, 2015, 2016, 2…
#> $ Y        <dbl> 7.23, 6.65, 6.54, 6.65, 8.72, 8.23, 8.04, 8.23, 4.99, 4.47, 4…
#> $ X1       <dbl> 35.80, 35.20, 35.53, 34.62, 32.23, 32.75, 31.93, 31.66, 28.57…
#> $ X2       <dbl> 74.10, 74.74, 78.25, 83.48, 48.53, 55.03, 48.32, 61.04, 43.86…
#> $ X3       <dbl> 50.24, 50.25, 50.21, 50.15, 50.73, 50.70, 50.68, 50.66, 49.60…
#> $ X4       <dbl> 53.25, 58.40, 71.39, 77.14, 29.98, 34.84, 46.09, 55.91, 24.54…
#> $ X5       <dbl> 15512.7, 12216.9, 47262.3, 49097.4, 26272.9, 30360.2, 38390.6…

keterangan setiap kolom:

  • Provinsi: Nama provinsi di Jawa Barat
  • Tahun: Tahun 2015 s.d 2018
  • Y: Tingkat pengangguran wilayah (persen)
  • X1: Persentase Penduduk Berumur 10 Tahun ke Atas Berstatus Kawin menurut Provinsi (persen)
  • X2: Tingkat Penyelesaian Pendidikan Menurut Jenjang Pendidikan SMA dan Provinsi (Persen)
  • X3: Persentase Penduduk Berjenis Kelamin Laki-Laki menurut Provinsi (Persen)
  • X4: Proporsi Remaja Dan Dewasa Usia 15-59 Tahun Dengan Keterampilan Teknologi Informasi Dan Komputer (TIK) Menurut Provinsi (Persen)
  • X5: Realisasi Investasi Penanaman Modal Dalam Negeri Menurut Provinsi (Milyar Rupiah)

Pemeriksaan Keseimbangan Data

Untuk memeriksa apakah data kita sudah seimbang dapat dilakukan dengan fungsi is.pbalanced() yang telah tersedia pada packages plm.

is.pbalanced(data_tpt)
#> [1] TRUE

dari hasil pemeriksaan keseimbangan data dapat diketahui bahwa data yang akan kita sudah seimbang (Balanced).


Eksplorasi Data

Ringkasan Data

data_tpt %>% summary()
#>           Provinsi  Tahun         Y               X1              X2       
#>  BANTEN       :4   2015:6   Min.   :3.180   Min.   :27.32   Min.   :43.86  
#>  DI YOGYAKARTA:4   2016:6   1st Qu.:4.030   1st Qu.:29.15   1st Qu.:52.72  
#>  DKI JAKARTA  :4   2017:6   Median :5.765   Median :31.40   Median :60.37  
#>  JAWA BARAT   :4   2018:6   Mean   :5.973   Mean   :31.45   Mean   :63.25  
#>  JAWA TENGAH  :4            3rd Qu.:8.140   3rd Qu.:33.74   3rd Qu.:75.62  
#>  JAWA TIMUR   :4            Max.   :9.550   Max.   :35.80   Max.   :85.53  
#>        X3              X4              X5         
#>  Min.   :49.34   Min.   :24.54   Min.   :  294.6  
#>  1st Qu.:49.43   1st Qu.:34.00   1st Qu.:12374.0  
#>  Median :49.88   Median :45.79   Median :21968.2  
#>  Mean   :50.05   Mean   :45.84   Mean   :23877.7  
#>  3rd Qu.:50.69   3rd Qu.:56.27   3rd Qu.:36215.0  
#>  Max.   :51.02   Max.   :77.14   Max.   :49097.4

Berdasarkan ringkasan diatas dapat kita ketahui beberapa hal berikut:

  • [Y] : Rata-rata tingkat pengangguran terbuka di pulau Jawa adalah 5.973% dengan nilai terrendah di 3.180% (DI Yogyakarta) dan nilai tertinggi di 9.55% (Banten)
  • [X1] : Rata-rata persentase status perkawinan di pulau Jawa adalah 31.45% artinya penduduk berusia 10 tahun ke atas lebih banyak yang belum menikah dibandingkan dengan yang sudah menikah.
  • [X2] : Rata-rata persentase penyelesaian pendidikan tingkat SMA di pulau Jawa adalah 63.25%, artinya 63.35% penduduk di pulau Jawa sudah mengenyam pendidikan hingga tingkat SMA.
  • [X3] : Rata-rata persentase penduduk dengan jenis kelamin Laki-Laki di pulau Jawa adalah 50.05%, artinya jumlah penduduk Pria dan Wanita di pulau Jawa seimbang.
  • [X4] : Rata-rata persentase penduduk dengan kemampuan TIK di pulau jawa adalah sebesar 45.84%, artinya penduduk usia 15-59 tahun di pulau Jawa lebih banyak yang tidak memiliki kemampuan TIK.
  • [X5] : Terdapat wilayah yang memiliki tingkat realisasi investasi sangat kecil yaitu sebesar 294.6 Milyar (DI Yogyakarta) dibandingkan dengan rata-rata realisasi investasi sebesar 23877.7 Milyar
#wilayah dengan nilai tingkat pengangguran terbuka terrendah
data_tpt[data_tpt$Y==min(data_tpt$Y),]
Provinsi Tahun Y X1 X2 X3 X4 X5
15 DI YOGYAKARTA 2017 3.18 30.51 85.53 49.59 57.37 294.6
#wilayah dengan nilai tingkat pengangguran terbuka tertinggi
data_tpt[data_tpt$Y==max(data_tpt$Y),]
Provinsi Tahun Y X1 X2 X3 X4 X5
21 BANTEN 2015 9.55 33.37 52.95 51.02 31.48 10709.9
#wilayah dengan tingkat realisasi investasi terrendah
data_tpt[data_tpt$X5==min(data_tpt$X5),]
Provinsi Tahun Y X1 X2 X3 X4 X5
15 DI YOGYAKARTA 2017 3.18 30.51 85.53 49.59 57.37 294.6

Matriks Korelasi (Heatmap)

untuk mengetahui seberapa besar tingkat hubungan antar variabel prediktor terhadap variabel target dapat kita gunakan fungsi ggcorrplot.

data_tpt %>% select(c(-1,-2)) %>% cor() %>% 
ggcorrplot(type = "lower",lab = TRUE)

Berdasarkan hasil plot heatmap diatas, dapat diketahui bahwa:

  • Terdapat dua variabel yang memiliki korelasi tinggi terhadap persentase tingkat pengangguran terbuka yakni Variabel Persentase Penduduk Berjenis Kelamin Laki-Laki menurut Provinsi dan Persentase Penduduk Berumur 10 Tahun ke Atas berstatus kawin menurut Provinsi.
  • Sedangkan ketiga variabel lainnya memiliki korelasi yang relatif kecil terhadap persentase tingkat pengangguran terbuka.
  • terdapat indikasi terjadi multikolinieritas antara variabel X2 dengan X4 dan variabel X1 dengan X3.

Pemodelan

Cross-Validation

Tahapan cross validation baik akan selalu dilakukan sebelum pembuatan model, data akan dibagi menjadi data train dan data test. Dikarenakan data panel memiliki informasi keterangan waktu maka pembagian data tidak boleh diambil secara acak melainkan dibagi dengan cara dipisah secara berurutan.

Data Train akan menggunakan data yang paling awal Data Test akan menggunakan data yang paling akhir

untuk melakukannya kita bisa menggunakan bantuan fungsi filter()

#membuat data train
data_train <- data_tpt %>% filter( Tahun != 2018)

#membuat data test
data_test <- data_tpt %>% filter( Tahun %in% 2018)

Pemeriksaan Asumsi Multikolinieritas

Dikarenakan pada hasil pemeriksaan korelasi pada tahapan EDA sebelumnya menunjukkan adanya indikasi multikolinieritas antar variabel prediktor, maka akan dilakukan pemeriksaan asumsi multikolinieritas terlebih dahulu dengan cara pembuatan model regresi dengan fungsi lm() dan dilanjutkan pengujian menggunakan fungsi vif().

nilai VIF > 10: terjadi multicollinearity pada model nilai VIF < 10: tidak terjadi multicollinearity pada model

lm(Y~X1+X2+X3+X4+X5, data_train) %>% vif()
#>       X1       X2       X3       X4       X5 
#> 8.840779 6.035948 6.322624 5.607408 1.709643

Berdasarkan hasil pengujian VIF diatas didapati keseluruhan variabel prediktor memiliki nilai VIF < 10, artinya tidak terjadi multikolinieritas pada model. Sehingga seluruh variabel prediktor dapat digunakan untuk pembuatan model.

Penentuan Model Estimasi

Pembuatan Model

Untuk setiap pembuatan model akan digunakan fungsi plm() dari package plm dengan parameter sebagai berikut:

  • formula = Target ~ Prediktor
  • data = berupa dataframe
  • index = c(“kolom_individu”,“kolom_waktu”)
  • model =
    • “pooling” : untuk model CEM
    • “within” : untuk model FEM
    • “random” : untuk model REM

dengan variabel target

  • variabel target : [Y]
  • Variabel prediktor : [X1 s.d X5]

Model Gabungan (CEM)

membuat model CEM dan disimpan kedalam objek cem

cem <- plm(Y~X1+X2+X3+X4+X5,
           data = data_train,
           index = c("Provinsi","Tahun"),
           model = "pooling")

Model Pengaruh Tetap (FEM)

membuat model FEM dengan memberikan parameter tambahan effect = "twoways" untuk memasukan pengaruh individu dan waktu kemudian disimpan kedalam objek fem.two

fem.two <- plm(Y~X1+X2+X3+X4+X5,
           data = data_train,
           index = c("Provinsi","Tahun"),
           model = "within")

Uji Chow

Uji chow dilakukan untuk memilih model terbaik antara model gabungan (cem) dengan model fixed effec (fem). untuk melakukan uji Chow dapat menggunakan fungsi pooltest(model_cem, model_fem)

Hipotesis yang diuji adalah sebagai berikut:

  • H0 : Model gabungan
  • H1 : Model pengaruh tetap

H0 ditolak jika P-value < α. Nilai α yang digunakan sebesar 5%.

pooltest(cem,fem.two)
#> 
#>  F statistic
#> 
#> data:  Y ~ X1 + X2 + X3 + X4 + X5
#> F = 1.8609, df1 = 5, df2 = 7, p-value = 0.2198
#> alternative hypothesis: unstability

Berdasarkan hasil pengujian diatas kita peroleh nilai p-value > α (5%), artinya dengan tingkat keyakinan 95% kita yakin bahwa model pengaruh gabungan merupakan model yang lebih baik untuk digunakan dibandingkan oleh model tetap. Sehingga kita tidak perlu membuat model random. dan dapat dilanjutkan kepada tahapan selanjutnya untuk melihat apakah terdapat pengaruh dari Individu / Waktu / Individu & Waktu terhadap model.


Pemeriksaan Pengaruh Individu & Waktu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh individu & waktu
  • H1 : Ada pengaruh individu & waktu
plmtest(cem, effect = "twoways", type = "bp")
#> 
#>  Lagrange Multiplier Test - two-ways effects (Breusch-Pagan)
#> 
#> data:  Y ~ X1 + X2 + X3 + X4 + X5
#> chisq = 12.006, df = 2, p-value = 0.002471
#> alternative hypothesis: significant effects

Pemeriksaan Pengaruh Individu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh individu
  • H1 : Ada pengaruh individu
plmtest(cem, effect = "individual", type = "bp")
#> 
#>  Lagrange Multiplier Test - (Breusch-Pagan)
#> 
#> data:  Y ~ X1 + X2 + X3 + X4 + X5
#> chisq = 1.423, df = 1, p-value = 0.2329
#> alternative hypothesis: significant effects

Pemeriksaan Pengaruh Waktu

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Tidak ada pengaruh waktu
  • H1 : Ada pengaruh waktu
plmtest(cem, effect = "time", type = "bp")
#> 
#>  Lagrange Multiplier Test - time effects (Breusch-Pagan)
#> 
#> data:  Y ~ X1 + X2 + X3 + X4 + X5
#> chisq = 10.583, df = 1, p-value = 0.001141
#> alternative hypothesis: significant effects

Berdasarkan hasil pengujian Breusch Pagan diperoleh kesimpulan bahwa pada model gabungan terdapat pengaruh dua arah. Namun, setelah diuji pengaruh individu dan waktu, hanya terdapat pengaruh waktu. Sehingga model yang terbentuk adalah model gabungan dengan pengaruh satu arah, yaitu pengaruh waktu.

#model gabungan dengan pengaruh waktu
cem.time <- plm(Y~X1+X2+X3+X4+X5,
           data = data_train,
           index = c("Provinsi","Tahun"),
           model = "pooling",
           effect = "time")

Pengujian Asumsi

Normalitas

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Sisaan menyebar normal
  • H1 : Sisaan tidak menyebar normal

H0 ditolak jika P-value < α. Nilai α yang digunakan sebesar 5%.

cem.time$residuals %>% shapiro.test()
#> 
#>  Shapiro-Wilk normality test
#> 
#> data:  .
#> W = 0.96132, p-value = 0.6273

Berdasarkan hasil pengujian normalitas diperoleh nilai p-value > 0.05, artinya sisaan menyebar secara normal.

Homogenitas

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : Sisaan memiliki ragam homogen
  • H1 : Sisaan tidak memiliki ragam homogen

H0 ditolak jika P-value < α. Nilai α yang digunakan sebesar 5%.

cem.time %>% bptest()
#> 
#>  studentized Breusch-Pagan test
#> 
#> data:  .
#> BP = 6.3898, df = 5, p-value = 0.2701

Berdasarkan hasil pengujian homogenitas diperoleh nilai p-value > 0.05, artinya sisaan memiliki ragam yang homogen.

Autokorelasi

Hipotesis yang diuji adalah sebagai berikut.

  • H0 : tidak terjadi autokorelasi pada sisaan
  • H1 : terjadi autokorelasi pada sisaan

H0 ditolak jika P-value < α. Nilai α yang digunakan sebesar 5%.

cem.time$residuals %>% Box.test(type = "Ljung-Box")
#> 
#>  Box-Ljung test
#> 
#> data:  .
#> X-squared = 1.5905, df = 1, p-value = 0.2073

Berdasarkan hasil pengujian autokorelasi diperoleh nilai p-value > 0.05, artinya tidak terjadi autokorelasi antar sisaan.


Prediksi

Untuk melakukan prediksi akan kita gunakan fungsi predict() dengan parameter:

  • object = nama model yang kita gunakan
  • newdata = data baru yang akan kita prediksi
predict(cem.time,newdata = data_test)
#>        1        2        3        4        5        6 
#> 5.819528 7.041081 4.003161 2.469587 2.984266 7.935911

Untuk menguji apakah model yang kita miliki sudah baik dalam memprediksi data baru maka kita akan evaluasi dengan menggunakan nilai error, salah satu metric eror yang biasa digunakan adalah MAPE. Kita dapat melakukannya menggunakan fungsi MAPE() dengan parameter:

  • y_pred = nilai hasil prediksi
  • y_true = nilai target asli
MAPE(y_pred = predict(cem.time,newdata = data_test),
    y_true = data_test$Y)
#> [1] 0.1567975

Diperoleh nilai MAPE sebesar 0.1567975 artinya model kita dalam memprediksi data baru mengalami kesalahan prediksi sebesar 15.68% atau dapat dikatakan model kita sudah cukup baik dalam memprediksi nilai baru.


Interpretasi Model

summary(cem.time)
#> Pooling Model
#> 
#> Call:
#> plm(formula = Y ~ X1 + X2 + X3 + X4 + X5, data = data_train, 
#>     effect = "time", model = "pooling", index = c("Provinsi", 
#>         "Tahun"))
#> 
#> Balanced Panel: n = 6, T = 3, N = 18
#> 
#> Residuals:
#>      Min.   1st Qu.    Median   3rd Qu.      Max. 
#> -0.646827 -0.128901  0.033248  0.239362  0.474528 
#> 
#> Coefficients:
#>                    Estimate      Std. Error t-value    Pr(>|t|)    
#> (Intercept) -110.1820834249   14.9320138557 -7.3789 0.000008514 ***
#> X1             0.3443627134    0.0916460360  3.7575    0.002733 ** 
#> X2            -0.0129411224    0.0155167856 -0.8340    0.420569    
#> X3             2.1502739738    0.3358762129  6.4020 0.000033925 ***
#> X4            -0.0458480802    0.0158858977 -2.8861    0.013675 *  
#> X5             0.0000175280    0.0000072536  2.4165    0.032527 *  
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Total Sum of Squares:    77.762
#> Residual Sum of Squares: 1.5284
#> R-Squared:      0.98034
#> Adj. R-Squared: 0.97215
#> F-statistic: 119.704 on 5 and 12 DF, p-value: 0.00000000082424

Berdasarkan output diatas dapat kita peroleh informasi sebagai berikut:

  • Bahwa variabel yang mempengaruhi persentase tingkat pengangguran terbuka adalah Rata-rata persentase status perkawinan, Rata-rata persentase penduduk dengan jenis kelamin Laki-Laki dan tingkat realisasi investasi. Hal ini ditunjukkan dengan nilai p-value untuk variabel X1,X3, X4 dan X5 yang kurang dari 0.05.
  • Estimasi koefisien X1 bernilai positif artinya semakin tinggi persentase status perkawinan di suatu wilayah akan meningkatkan persentase tingkat pengagguran. Hal ini linier dengan kondisi dilapangan dimana banyak penyedia lapangan pekerjaan memberikan syarat kepada pelamar untuk mendapatkan pekerjaan adalah belum menikah/belum kawin.
  • Estimasi koefisien X2 bernilai negatif artinya semakin tinggi persentase penyelesaian pendidikan tingkat SMA di suatu wilayah akan mengurangi persentase tingkat pengangguran. Hal ini linier dengan kondisi di lapangan dimana banyak penyedia lapangan pekerjaan memberikan syarat kepada pelamar untuk mendapatkan pekerjaan adalah memiliki pendidikan minimal tingkat SMA .
  • Estimasi koefisien X3 bernilai positif artinya semakin tinggi persentase penduduk laki-laki di suatu wilayah akan meningkatkan persentase tingkat pengaguran.
  • Estimasi koefisien X4 bernilai negatif artinya semakin tinggi persentase penduduk dengan kemampuan TIK di suatu wilayah akan mengurangi persentase tingkat pengangguran. Hal ini linier dengan kondisi di lapangan dimana banyak penyedia lapangan pekerjaan memberikan syarat kepada pelamar untuk mendapatkan pekerjaan adalah mampu dalam mengoperasikan komputer.
  • Estimasi koefisien X5 bernilai positif artinya semakin tinggi tingkat realisasi investasi di suatu wilayah akan meningkatkan persentase tingkat pengangguran.
  • Diperoleh nilai adj. R-squared sebesar 0.97215, artinya model dapat menjelaskan persentase tingkat pengangguran dengan baik sebesar 96.67%

Kesimpulan & Saran

Dari serangkaian proses analisis yang telah dilakukan, dapat kita peroleh kesimpulan sebagai berikut:

  1. Berdasarkan tahapan EDA kita peroleh informasi bahwa rata-rata tingkat pengangguran terbuka di pulau Jawa adalah 5.973% dengan nilai terrendah sebesar 3.180% untuk wilayah DI Yogyakarta dan nilai tertinggi sebesar 9.55% untuk wilayah Banten.
  2. Diperoleh model terbaik untuk memprediksi tingkat pengangguran adalah model Common Effect Model (CEM) dengan pengaruh waktu.
  3. Dari kelima variabel prediktor yang digunakan hanya variabel persentase penyelesaian tingkat pendidikan SMA yang tidak signifikan berpengaruh terhadap persentase tingkat pengangguran suatu wilayah.
  4. Melihat dari hasil pengamatan pada tahapan EDA dan uji signifikansi prediktor. Diperoleh variabel penting yakni persentase penduduk dengan kemampuan TIK di suatu wilayah positif yang berpengaruh dalam menurunkan persentase tingkat pengangguran di suatu wilayah.
  5. Diperoleh nilai adj. R-squared sebesar 0.97215, artinya model dapat menjelaskan persentase tingkat pengangguran dengan baik sebesar 96.67%.
  6. Diperoleh nilai MAPE sebesar 0.1567975 artinya model kita dalam memprediksi data baru mengalami kesalahan prediksi sebesar 15.68%. Sehingga dapat dikatakan model yang kita miliki sudah cukup baik dalam memprediksi nilai baru di masa mendatang.

Dari beberapa kesimpulan diatas terdapat beberapa saran sebagai berikut:

  1. Pada penelitan kali ini masih didapati kendala dimana terbatasnya data yang digunakan, yakni 6 wilayah dan 4 periode waktu. Hal ini dikhawatirkan akan membuat model kurang baik dalam mempelajari data. Sehingga diharapkan pada penelitan selanjutnya bisa mencoba dilakukan analisis di wilayah lainnya seperti seluruh provinsi di Indonesia dan atau menambah rentang waktu dari data.
  2. Memperhatikan kesimpulan poin 4, diharapkan pemerintah bisa lebih mengedepankan pembekalan kemampuan pengoperasian komputer bagi penduduk di suatu wilayah dikarenakan hal ini dapat mendorong penurunan persentase tingkat pengagguran di suatu wilayah.