Library

library(readxl)
## Warning: package 'readxl' was built under R version 4.4.2
library(GGally)
## Warning: package 'GGally' was built under R version 4.4.3
## Loading required package: ggplot2
## Warning: package 'ggplot2' was built under R version 4.4.2
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
library(corrplot)
## Warning: package 'corrplot' was built under R version 4.4.3
## corrplot 0.95 loaded
library(nortest)
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.4.3
## 
## 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(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.4.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.4.2
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(car) 
## Warning: package 'car' was built under R version 4.4.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.4.2
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
library(olsrr)
## Warning: package 'olsrr' was built under R version 4.4.3
## 
## Attaching package: 'olsrr'
## The following object is masked from 'package:datasets':
## 
##     rivers
library(MASS)
## Warning: package 'MASS' was built under R version 4.4.3
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:olsrr':
## 
##     cement
## The following object is masked from 'package:dplyr':
## 
##     select
library(psych)
## Warning: package 'psych' was built under R version 4.4.3
## 
## Attaching package: 'psych'
## The following object is masked from 'package:car':
## 
##     logit
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
library(car)

Import Data

dt <- read_excel("C:/Users/LENOVO/Downloads/Data tugas praktikum anreg.xlsx")
head(dt)
## # A tibble: 6 × 6
##   IPM     TPT  TPAK   PPM  SAML    RKMRT
##   <chr> <dbl> <dbl> <dbl> <dbl>    <dbl>
## 1 72.8   6.17  63.5 14.8   89.7 2773116.
## 2 72.71  6.16  69.5  8.33  92.1 2653705.
## 3 73.26  6.28  69.3  6.04  85.2 2995544.
## 4 73.52  4.37  63.9  6.84  90.1 2953900.
## 5 72.14  4.59  67.8  7.7   79.2 2604722.
## 6 70.9   4.63  69.3 12.0   86.4 2488706.

Penjelasan Variabel

Y _> Index Pembangunan Manusia,

X1 _> Tingkat Pengangguran Terbuka

X2 _> Tingkat Partisipasi Angkatan Kerja (%)

X3 _> Persentase Penduduk Miskin

X4 _> Persentase Sumber Air Minum Layak

X5 _> Rata-rata konsumsi Makanan Rumah Tangga

dt_numeric <- data.frame(lapply(dt, function(x) as.numeric(as.character(x))))

Eksplorasi Data

Peubah Penjelas vs Peubah Respon

TPT

plot(dt_numeric$IPM, dt_numeric$TPT)

TPAK

plot(dt_numeric$IPM, dt_numeric$TPAK)

PPM

plot(dt_numeric$IPM, dt_numeric$PPM)

SAML

plot(dt_numeric$IPM, dt_numeric$SAML)

RKMRT

plot(dt_numeric$IPM, dt_numeric$RKMRT)

Hubungan Antar Peubah

cor_matrix <- cor(dt_numeric, use = "complete.obs")

corrplot(cor_matrix, method = "color", type = "lower",
         col = colorRampPalette(c("red", "white", "blue"))(200),
         addCoef.col = "black", tl.col = "black", tl.srt = 35)

library(GGally)
library(psych)
pairs.panels(dt_numeric)

ggpairs(dt_numeric,
        upper = list(continuous = wrap('cor', size = 3)),
        title = "Matriks Scatterplot Data")

Berdasarkan hasil eksplorasi korelasi antara peubah penjelas dan peubah terikat dapat disimpulkan bahwa

IPM dan TPT memiliki hubungan positif sedang, sehingga ketika TPT naik maka IPM pun akan naik. Akan tetapi hasil korelasi ini berbanding terbalik dengan konsep IPM.

IPM dan TPAT memiliki hbungan negatif lemah, sehingga ketika TPAT naik maka IPM turun, Akan tetapi hasil korelasi ini berbanding terbalik dengan konsep IPM.

IPM dan PPM memiliki hubungan negatif kuat, sehingga ketika PPM naik maka IPM turun

IPM dan SAML memiliki hubungan postif kuat, sehingga ketika SAML naik maka IPM naik

IPM dan RKMRT hampir tidak memiliki hubungan karena mendekati nol

Permodelan Awal

model1 = lm(formula = IPM ~ ., data = dt_numeric)
summary(model1)
## 
## Call:
## lm(formula = IPM ~ ., data = dt_numeric)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.6440 -1.2412 -0.2265  0.3090  7.6984 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  5.450e+01  1.224e+01   4.454 0.000123 ***
## TPT          3.488e-01  3.535e-01   0.987 0.332293    
## TPAK         2.424e-02  1.497e-01   0.162 0.872528    
## PPM         -3.804e-01  1.005e-01  -3.787 0.000742 ***
## SAML         1.948e-01  6.111e-02   3.188 0.003508 ** 
## RKMRT        3.489e-07  8.224e-07   0.424 0.674586    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.533 on 28 degrees of freedom
## Multiple R-squared:  0.6423, Adjusted R-squared:  0.5784 
## F-statistic: 10.06 on 5 and 28 DF,  p-value: 1.368e-05

Pada permodelan awal terbentuk model regresi linear berganda :

IPM = 54.5 + 0.3448 TPT + 0.0242TPAK - 0.3804PPM + 0.1948 + 0.0000003RKMRT

Multikolinieritas

library(car)
vif(model1)
##      TPT     TPAK      PPM     SAML    RKMRT 
## 1.646228 1.528020 1.450515 1.185360 1.150418

Dari hasi VIF tidak terdapat peubah yang mengalami multikolinieritas (VIF > 10).

Seleksi Peubah

bs <- ols_step_best_subset(model1)
bs
##        Best Subsets Regression        
## --------------------------------------
## Model Index    Predictors
## --------------------------------------
##      1         PPM                     
##      2         PPM SAML                
##      3         TPT PPM SAML            
##      4         TPT PPM SAML RKMRT      
##      5         TPT TPAK PPM SAML RKMRT 
## --------------------------------------
## 
##                                                    Subsets Regression Summary                                                    
## ---------------------------------------------------------------------------------------------------------------------------------
##                        Adj.        Pred                                                                                           
## Model    R-Square    R-Square    R-Square     C(p)        AIC        SBIC        SBC         MSEP       FPE       HSP       APC  
## ---------------------------------------------------------------------------------------------------------------------------------
##   1        0.4537      0.4367      0.3878    12.7596    173.4736    76.0824    178.0527    291.4631    9.0757    0.2765    0.6145 
##   2        0.6211      0.5966      0.5458     1.6604    163.0370    67.3842    169.1424    208.9137    6.6790    0.2046    0.4523 
##   3        0.6398      0.6037      0.5457     2.1985    163.3185    68.3919    170.9503    205.4654    6.7388    0.2079    0.4563 
##   4        0.6420      0.5926      0.5239     4.0262    165.1101    70.6629    174.2683    211.5035    7.1109    0.2214    0.4815 
##   5        0.6423      0.5784      0.5008     6.0000    167.0783    73.0701    177.7628    219.1318    7.5466    0.2376    0.5110 
## ---------------------------------------------------------------------------------------------------------------------------------
## AIC: Akaike Information Criteria 
##  SBIC: Sawa's Bayesian Information Criteria 
##  SBC: Schwarz Bayesian Criteria 
##  MSEP: Estimated error of prediction, assuming multivariate normality 
##  FPE: Final Prediction Error 
##  HSP: Hocking's Sp 
##  APC: Amemiya Prediction Criteria

Berdasarkan best subset dapat disimpulkan bahwa hanya PPM dan SAML (model 2) yang dapat digunakan.

Backward

null<-lm(IPM ~ 1, data=dt_numeric) 
full<-lm(IPM ~ ., data=dt_numeric)
 
step(full, scope=list(lower=null, upper=full),data=dt_numeric, direction='backward', trace=0)
## 
## Call:
## lm(formula = IPM ~ PPM + SAML, data = dt_numeric)
## 
## Coefficients:
## (Intercept)          PPM         SAML  
##     57.4838      -0.4035       0.2127

Berdasarkan metode Backward dapat disimpulkan bahwa hanya PPM dan SAML yang dapat digunakan dalam model.

Foward

step(null, scope=list(lower=null, upper=full),data=data, direction='forward', trace=0)
## 
## Call:
## lm(formula = IPM ~ PPM + SAML, data = dt_numeric)
## 
## Coefficients:
## (Intercept)          PPM         SAML  
##     57.4838      -0.4035       0.2127

Berdasarkan metode Foward dapat disimpulkan bahwa hanya PPM dan SAML yang dapat digunakan dalam model.

Stepwise (Campuran)

step(null, scope=list(lower=null, upper=full),data=data, direction='both', trace=0)
## 
## Call:
## lm(formula = IPM ~ PPM + SAML, data = dt_numeric)
## 
## Coefficients:
## (Intercept)          PPM         SAML  
##     57.4838      -0.4035       0.2127

Berdasarkan metode stepwise dapat disimpulkan bahwa hanya PPM dan SAML yang dapat digunakan dalam model.

Kesimpulan

IPM = 57.4838 - 0.4035 PPM + 0.2127 SAML

model2 <- lm(IPM ~ PPM + SAML, data=dt_numeric)
vif(model2)
##      PPM     SAML 
## 1.096144 1.096144

Uji Asumsi

Uji Formal

Nilai Harapan Sisaaan = 0

t.test(model2$residuals,mu = 0,conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  model2$residuals
## t = -1.348e-16, df = 33, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.8377975  0.8377975
## sample estimates:
##     mean of x 
## -5.551115e-17
  • H₀ : Mean residual = 0 (rata-rata sisa = 0)
  • H₁ : Mean residual ≠ 0

Berdasarkan hasil uji One Sample t-test terhadap residual, diperoleh p-value sebesar 1 (lebih besar dari 0.05). Dengan demikian, tidak terdapat cukup bukti untuk menolak hipotesis nol. Oleh karena itu, dapat disimpulkan bahwa rata-rata sisa adalah nol, sehingga asumsi nilai harapan residual = 0 terpenuhi.

Ragam Sisaan Homogen

bptest(model2)
## 
##  studentized Breusch-Pagan test
## 
## data:  model2
## BP = 0.65066, df = 2, p-value = 0.7223
  • H₀: Varians residual homogen (homoskedastisitas ada)
  • H₁ : Varians residual tidak homogen (terjadi heteroskedastisitas)

Berdasarkan uji Breusch-Pagan, diperoleh p-value sebesar 0.7223 (lebih besar dari 0.05). Dengan demikian, tidak terdapat cukup bukti untuk menolak hipotesis nol. Oleh karena itu, dapat disimpulkan bahwa model memenuhi asumsi homoskedastisitas (ragam residual homogen).

Sisaan Saling Bebas

dwtest(model2)
## 
##  Durbin-Watson test
## 
## data:  model2
## DW = 1.8346, p-value = 0.2494
## alternative hypothesis: true autocorrelation is greater than 0
  • H₀: Tidak ada autokorelasi pada residual
  • H₁: Ada autokorelasi positif

Berdasarkan uji Durbin–Watson, diperoleh statistik DW = 1.8346 dengan p-value sebesar 0.2494 (lebih besar dari 0.05). Dengan demikian, tidak terdapat cukup bukti untuk menolak hipotesis nol. Oleh karena itu, dapat disimpulkan bahwa model memenuhi asumsi residual saling bebas (tidak ada autokorelasi pada sisa).

Normalitas Sisaan

shapiro.test(model2$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  model2$residuals
## W = 0.87899, p-value = 0.001337

Berdasarkan uji Shapiro-Wilk, diperoleh statistik W = 0.87899 dengan p-value sebesar 0.001337 (lebih kecil dari 0.05). Dengan demikian, terdapat cukup bukti untuk menolak hipotesis nol. Oleh karena itu, dapat disimpulkan bahwa residual model tidak berdistribusi normal.

Box - Cox

bc_model <- boxcox(IPM ~ PPM + SAML, data = dt,
         lambda = seq(-2, 2, by = 0.1))

interpretasi

(optimal_lambda <- bc_model$x[which.max(bc_model$y)])
## [1] -2

Berdasarkan hasil analisis Box-Cox transformation, diperoleh nilai optimal lambda sebesar -2. Nilai lambda ini menunjukkan bentuk transformasi terbaik terhadap variabel respons (IPM) untuk mendekati asumsi normalitas residual. Karena nilai lamda sebesar -2 maka Ybaru = 1/Y^2. Artinya, transfomasi IPM dilakukan dengan mengambil invers kuadrat dari IPM.

Transformasi Invers Kuadrat

lambda_opt <- -2
dt_numeric$IPM_bc <- (dt_numeric$IPM^lambda_opt - 1) / lambda_opt
model_bc <- lm(IPM_bc ~ PPM + SAML, data = dt_numeric)
summary(model_bc)
## 
## Call:
## lm(formula = IPM_bc ~ PPM + SAML, data = dt_numeric)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -1.131e-05 -2.259e-06 -1.067e-06  6.711e-07  1.668e-05 
## 
## Coefficients:
##               Estimate Std. Error   t value Pr(>|t|)    
## (Intercept)  4.999e-01  1.336e-05 37421.863  < 2e-16 ***
## PPM         -1.225e-06  2.120e-07    -5.776 2.32e-06 ***
## SAML         5.946e-07  1.427e-07     4.168 0.000229 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.149e-06 on 31 degrees of freedom
## Multiple R-squared:  0.6968, Adjusted R-squared:  0.6772 
## F-statistic: 35.62 on 2 and 31 DF,  p-value: 9.27e-09
shapiro.test(residuals(model_bc))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(model_bc)
## W = 0.90373, p-value = 0.005765

interpretasi

Cek Ulang Kehomogenan

Uji Formal

bptest(model2)
## 
##  studentized Breusch-Pagan test
## 
## data:  model2
## BP = 0.65066, df = 2, p-value = 0.7223

  • Hipotesis nol (H₀): Tidak ada heteroskedastisitas (varian residual konstan).

  • Hipotesis alternatif (H₁): Ada heteroskedastisitas (varian residual tidak konstan).

Berdasarkan hasil analisis uji Breusch-Pagan, diperoleh nilai p-value sebesar 0.9944. Nilai p-value yang lebih besar dari 0.05 ini menunjukkan bahwa tidak terdapat bukti adanya heteroskedastisitas pada model (model3). Dengan demikian, asumsi homoskedastisitas (varian residual konstan) pada model telah terpenuhi.

Uji Grafik

plot(model2,1)

Berdasarkan hasil plot Residuals vs Fitted pada model (model3), terlihat bahwa sebaran residual berada secara acak di sekitar garis horizontal nol dan tidak membentuk pola tertentu. Hal ini menunjukkan bahwa asumsi homogenitas varians (homoskedastisitas) dan linearitas hubungan antara variabel bebas dan variabel respon telah terpenuhi.