Packages

library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.3.2
## Warning: package 'ggplot2' was built under R version 4.3.2
## Warning: package 'tidyr' was built under R version 4.3.3
## Warning: package 'readr' was built under R version 4.3.2
## Warning: package 'stringr' was built under R version 4.3.2
## Warning: package 'forcats' was built under R version 4.3.2
## Warning: package 'lubridate' was built under R version 4.3.2
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.3     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.4.4     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(ggridges)
## Warning: package 'ggridges' was built under R version 4.3.2
library(treemap)
## Warning: package 'treemap' was built under R version 4.3.2
library(treemapify)
## Warning: package 'treemapify' was built under R version 4.3.3
library(GGally)
## Warning: package 'GGally' was built under R version 4.3.2
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
library(plotly)
## Warning: package 'plotly' was built under R version 4.3.3
## 
## Attaching package: 'plotly'
## 
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## 
## The following object is masked from 'package:stats':
## 
##     filter
## 
## The following object is masked from 'package:graphics':
## 
##     layout
library(ggplot2)
library(readxl)
## Warning: package 'readxl' was built under R version 4.3.2
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## 
## The following object is masked from 'package:dplyr':
## 
##     recode
## 
## The following object is masked from 'package:purrr':
## 
##     some
library(dplyr)
library(randtests)
## Warning: package 'randtests' was built under R version 4.3.3
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.3.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.3.2
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(olsrr)
## Warning: package 'olsrr' was built under R version 4.3.3
## 
## Attaching package: 'olsrr'
## 
## The following object is masked from 'package:datasets':
## 
##     rivers
library(corrplot)
## Warning: package 'corrplot' was built under R version 4.3.2
## corrplot 0.92 loaded
library(nortest)
library(MASS)
## 
## Attaching package: 'MASS'
## 
## The following object is masked from 'package:olsrr':
## 
##     cement
## 
## The following object is masked from 'package:plotly':
## 
##     select
## 
## The following object is masked from 'package:dplyr':
## 
##     select

Data

library(readxl)
umk <- read_xlsx("C:/## SEMESTER 4 IPB UNIVERSITY/#MATA KULIAH/ANALISIS REGRESI/UAS/TUGAS AKHIR/DATA123.xlsx")
umk
## # A tibble: 27 × 12
##    `Kabupaten/Kota`      `Upah Minimum (Juta Rupiah)` KHL Per Orang (Juta Rupi…¹
##    <chr>                                        <dbl>                      <dbl>
##  1 KABUPATEN BOGOR                               4.52                      0.683
##  2 KABUPATEN SUKABUMI                            3.35                      0.653
##  3 KABUPATEN CIANJUR                             2.89                      0.575
##  4 KABUPATEN BANDUNG                             3.49                      0.631
##  5 KABUPATEN GARUT                               2.12                      0.615
##  6 KABUPATEN TASIKMALAYA                         2.50                      0.63 
##  7 KABUPATEN CIAMIS                              2.02                      0.59 
##  8 KABUPATEN KUNINGAN                            2.01                      0.626
##  9 KABUPATEN CIREBON                             2.43                      0.63 
## 10 KABUPATEN MAJALENGKA                          2.18                      0.585
## # ℹ 17 more rows
## # ℹ abbreviated name: ¹​`KHL Per Orang (Juta Rupiah)`
## # ℹ 9 more variables: `Angka Inflasi (Persen)` <dbl>,
## #   `Tingkat Pengangguran Terbuka (Persen)` <dbl>,
## #   `Pertumbuhan Ekonomi (Persen)` <dbl>, `Total APBD (Miliar)` <dbl>,
## #   `Penanaman Modal Dalam Negeri (Triliun)` <dbl>,
## #   `Tingkat Partisipasi Angkatan Kerja (Persen)` <dbl>, …
Y <- umk$`Upah Minimum (Juta Rupiah)`
X1 <- umk$`KHL Per Orang (Juta Rupiah)`
X2 <-umk$`Angka Inflasi (Persen)`
X3<-umk$`Tingkat Pengangguran Terbuka (Persen)`
X4 <- umk$`Pertumbuhan Ekonomi (Persen)`
X5<-umk$`Total APBD (Miliar)`
X6 <-umk$`Penanaman Modal Dalam Negeri (Triliun)`
X7 <- umk$`Tingkat Partisipasi Angkatan Kerja (Persen)`
X8<- umk$`Indeks Harga Konsumen (poin)`
X9 <- umk$`PDRB (Miliar)`
X10 <-umk$`Pendapatan Perkapita`

umk1 <- data.frame(Y,X1,X2,X3,X4,X5,X6,X7,X8,X9,X10)
umk1
##           Y    X1   X2    X3   X4    X5      X6    X7     X8      X9     X10
## 1  4.520212 0.683 3.36  7.14 5.51 12250 12.4720 64.22 118.92  98.700  55.206
## 2  3.351883 0.653 2.78  7.46 5.12  7320  4.6540 67.75 117.05 567.000  48.011
## 3  2.893229 0.575 3.14  7.27 4.97  4910  3.2730 72.31 117.40 235.000  44.576
## 4  3.492466 0.631 4.11  8.22 5.18 11310 18.3140 66.26 118.29 103.000  57.317
## 5  2.117318 0.615 3.02  7.13 5.03  5640  3.2340 70.10 112.44  30.300  43.433
## 6  2.499954 0.630 2.96  8.94 6.55  3624 22.1000 68.37 112.70  46.100  24.170
## 7  2.021657 0.590 2.89  7.64 5.07  5180  3.3480 66.26 112.02  27.035  45.452
## 8  2.010734 0.626 3.29  7.23 5.00  5630  3.9450 61.95 114.19  27.032  46.650
## 9  2.430781 0.630 3.21  7.75 5.36  7340  5.9040 66.16 112.38  33.500  66.439
## 10 2.180603 0.585 3.05 11.21 4.74  5110  3.4430 68.50 116.97  34.500  47.322
## 11 3.471134 0.611 3.22  6.80 5.24  5940  3.5330 67.76 115.80  42.400  36.733
## 12 2.541997 0.623 3.70  9.41 4.54  5940  3.6880 63.88 114.96  75.500  48.715
## 13 3.273811 0.625 4.90  7.94 5.31  6240  3.9150 70.03 116.32  49.900  52.454
## 14 4.464675 0.646 2.87  6.24 5.40  6290  4.1190 66.73 117.46  51.500  64.340
## 15 5.176179 0.743 3.10  5.98 5.74  8230  1.0692 63.40 119.30 112.100  88.314
## 16 5.137575 0.771 3.32  8.73 4.81 14960 16.9370 65.00 120.77 147.600  82.991
## 17 3.480795 0.685 3.74  6.63 5.00  6690 25.4600 67.01 112.48  56.900  53.886
## 18 2.018389 0.694 3.00  7.22 5.71  2330  2.4780 80.15 111.14  10.010  46.840
## 19 4.639429 0.732 3.36  8.00 4.91  3080 28.3000 64.81 117.96  57.000  74.977
## 20 2.893229 0.746 2.57  7.32 4.84  1156  0.3610 62.57 116.40  15.349  42.561
## 21 4.048463 0.740 0.63  6.34 4.72  7202 29.0600 66.97 116.16 351.280 140.140
## 22 2.456517 0.673 3.07  7.16 5.13  2710  0.8000 68.71 113.87  28.772  66.439
## 23 5.158248 0.787 3.32  9.02 4.79  7118  5.6800 64.65 119.30 118.963  45.280
## 24 4.694494 0.759 2.34  5.81 4.77  4372 15.7900 62.57 117.76  87.568  40.816
## 25 3.514093 0.710 3.08  9.74 5.05  1624  2.3000 68.43 118.57  40.499  68.550
## 26 2.533341 0.662 3.13  9.00 4.78  1156  1.1000 65.44 114.96  27.411  36.950
## 27 1.998119 0.651 2.81  7.84 4.89  2161  1.1200 67.44 116.94  43.327  25.280
umk1[] <- lapply(umk1, as.numeric)
head(umk1)
##          Y    X1   X2   X3   X4    X5     X6    X7     X8    X9    X10
## 1 4.520212 0.683 3.36 7.14 5.51 12250 12.472 64.22 118.92  98.7 55.206
## 2 3.351883 0.653 2.78 7.46 5.12  7320  4.654 67.75 117.05 567.0 48.011
## 3 2.893229 0.575 3.14 7.27 4.97  4910  3.273 72.31 117.40 235.0 44.576
## 4 3.492466 0.631 4.11 8.22 5.18 11310 18.314 66.26 118.29 103.0 57.317
## 5 2.117318 0.615 3.02 7.13 5.03  5640  3.234 70.10 112.44  30.3 43.433
## 6 2.499954 0.630 2.96 8.94 6.55  3624 22.100 68.37 112.70  46.1 24.170
str(umk1)
## 'data.frame':    27 obs. of  11 variables:
##  $ Y  : num  4.52 3.35 2.89 3.49 2.12 ...
##  $ X1 : num  0.683 0.653 0.575 0.631 0.615 0.63 0.59 0.626 0.63 0.585 ...
##  $ X2 : num  3.36 2.78 3.14 4.11 3.02 2.96 2.89 3.29 3.21 3.05 ...
##  $ X3 : num  7.14 7.46 7.27 8.22 7.13 ...
##  $ X4 : num  5.51 5.12 4.97 5.18 5.03 6.55 5.07 5 5.36 4.74 ...
##  $ X5 : num  12250 7320 4910 11310 5640 ...
##  $ X6 : num  12.47 4.65 3.27 18.31 3.23 ...
##  $ X7 : num  64.2 67.8 72.3 66.3 70.1 ...
##  $ X8 : num  119 117 117 118 112 ...
##  $ X9 : num  98.7 567 235 103 30.3 ...
##  $ X10: num  55.2 48 44.6 57.3 43.4 ...

Model Regresi Linear Berganda Awal

model1<- lm(Y~., data=umk1)
summary(model1)
## 
## Call:
## lm(formula = Y ~ ., data = umk1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.85924 -0.15173  0.05552  0.20403  0.71234 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.750e+01  6.393e+00  -4.301 0.000549 ***
## X1           6.377e+00  2.120e+00   3.008 0.008344 ** 
## X2           2.130e-01  1.890e-01   1.127 0.276433    
## X3          -1.896e-01  8.615e-02  -2.201 0.042768 *  
## X4           2.805e-01  2.602e-01   1.078 0.297027    
## X5           3.458e-05  4.046e-05   0.855 0.405360    
## X6           1.977e-02  1.244e-02   1.589 0.131572    
## X7          -1.060e-02  3.162e-02  -0.335 0.741827    
## X8           2.244e-01  5.288e-02   4.244 0.000619 ***
## X9           2.723e-05  1.012e-03   0.027 0.978872    
## X10          3.613e-03  5.727e-03   0.631 0.537045    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.475 on 16 degrees of freedom
## Multiple R-squared:  0.8817, Adjusted R-squared:  0.8078 
## F-statistic: 11.93 on 10 and 16 DF,  p-value: 1.216e-05

Didapatkan model regresi sebagai berikut: \[ \hat Y=-27.50+6.377X_1+0.213X_2-0.189X_3-0.280X_4+0.00003X_5+0.0197X_6-0.0106X_7+0.224X_8+0.00002X_9+0.0036_{10} \] Hasil tersebut belum bisa dipastikan menjadi model terbaik karena belum melalui serangkaian uji, sehingga diperlukan pengecekan nilai korelasi, multikolinearitas antarvariabel X dengan pembandingan nilai VIF yang lebih dari 10, dan kriteria model terbaik.

Eksplorasi Data

ggpairs(umk1,
        upper = list(continuous = wrap('cor', size = 2)),
        title = "Matriks Scatterplot dan Korelasi")

Pada matriks diatas tanda *(bintang) artinya terdapat nilai korelasi antara 2 peubah. Pada matriks tersebut, peubah-peubah yang memiliki korelasi terhadap respon Y adalah X1, X5, X6, X7, X8, X10. Sedangkan peubah X2, X3, X4,dan X9 memiliki korelasi yang kecil dan telah terwakili oleh variabel X yang lainnya sehingga dapat dihilangkan untuk kemudian diperiksa model dan multikolinearitasnya.

Model Regresi Linear Berganda Awal ( Setelah Reduksi Peubah Penjelas)

model2<- lm(Y~X1+X5+X6+X7+X8+X10, data=umk1)
summary(model2)
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X7 + X8 + X10, data = umk1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.90794 -0.26521  0.00891  0.31101  1.09809 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)   
## (Intercept) -2.385e+01  6.505e+00  -3.667  0.00153 **
## X1           7.259e+00  2.218e+00   3.273  0.00381 **
## X5           7.202e-05  3.850e-05   1.871  0.07610 . 
## X6           1.853e-02  1.341e-02   1.381  0.18240   
## X7          -1.386e-03  3.177e-02  -0.044  0.96565   
## X8           1.874e-01  5.210e-02   3.598  0.00180 **
## X10          1.163e-03  5.519e-03   0.211  0.83525   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5226 on 20 degrees of freedom
## Multiple R-squared:  0.8211, Adjusted R-squared:  0.7674 
## F-statistic:  15.3 on 6 and 20 DF,  p-value: 1.556e-06

Pengecekan Multikolinearitas

vif(model2)
##       X1       X5       X6       X7       X8      X10 
## 1.721655 1.534201 1.405791 1.271787 1.766033 1.521672

Berdasarkan Uji VIF diatas, nilai VIF dari tiap peubah penjelas (X1,X5,X6,X7,X8,X10) kurang dari 10. Hal ini, menandakan tidak terdapat multikolinearitas antar peubah X

Eksplorasi Kondisi Gauss Markov

Plot Sisaan vs Y duga

plot(model2,1)

  1. Sisaan menyebar di sekitar 0, sehingga nilai harapan galat sama dengan nol
  2. Lebar pita sama untuk setiap nilai dugaan, sehingga ragam homogen

Plot Sisaan vs Urutan

plot(x = 1:dim(umk1)[1],
     y = model2$residuals,
     type = 'b', 
     ylab = "Residuals",
     xlab = "Observation")

Tebaran tidak berpola, sehingga sisaan saling bebas, dan model pas

Normalitas Sisaan dengan QQ-Plot

plot(model2,2)

Titik-titik data pada QQ-plot cukup dekat dengan garis.Hal ini menunjukkan bahwa distribusi sisaan mirip dengan Sebaran normal.

Uji Formal Kondisi Gauss-Markov

1. Nilai harapan sisaan sama dengan nol

H0: Nilai harapan sisaan sama dengan nol
H1: Nilai harapan sisaan tidak sama dengan nol

t.test(model2$residuals,mu = 0,conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  model2$residuals
## t = -1.0193e-16, df = 26, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.1813143  0.1813143
## sample estimates:
##     mean of x 
## -8.990847e-18

Uji t.tes tersebut menunjukkan hasil p-value > alpha = 0.05, maka tak tolak H0,sehingg nilai harapan sisaan sama dengan nol pada taraf nyata 5%. Asumsi terpenuhi.

2. Ragam Sisaan Homogen

H0:ragam sisan homogen  H1:ragam siaan tidak homogen

bptest(model2)
## 
##  studentized Breusch-Pagan test
## 
## data:  model2
## BP = 4.2861, df = 6, p-value = 0.638

Uji ini sering disebut dengan uji homokesdasitas.Karena p-value = 0.638> alpha = 0.05, maka tak tolak H0, ragam sisaan homogen pada taraf nyata 5%. Asumsi terpenuhi.

3. Sisaan Saling Bebas

H0:sisaan saling bebas/tidak ada autokorelasi
H1:sisaan tidak saling bebas/ada autokorelasi

dwtest(model2)
## 
##  Durbin-Watson test
## 
## data:  model2
## DW = 1.9803, p-value = 0.3608
## alternative hypothesis: true autocorrelation is greater than 0

Uji ini sering disebut dengan uji autokorelasi yang dilakukan dengan Durbin watson. Berdasarakan Uji Durbin-Watson yang dilakukan, menghasilkan nilai p-value = 0.3608 > alpha = 0.05, maka tak tolak H0, sisaan saling bebas pada taraf nyata 5%, sehingga asumsi terpenuhi.

Uji Formal Normalitas Sisaan

H0:sisaan menyebar Normal
H1:sisaan tidak menyebar Normal

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

Berdasarkan Uji Shapiro-test yang dilakukan mengahasilkan nilai p-value = 0.8346 > alpha = 0.05, maka tak tolak H0, sehingga sisaan menyebar normal pada taraf nyata 5%.

Pencilan, Titik Leverage, dan Amatan Berpengaruh

Perhitungan ri, ci, Hi, DFFITSi

index <- c(1:27)
ri <- studres(model2)
ci <- cooks.distance(model2)
DFFITSi <- dffits(model2)
hi <- ols_hadi(model2)
hii <- hatvalues(model2)
hasil <- data.frame(index,ri,ci,DFFITSi,hi,hii); round(hasil,4)
##    index      ri     ci DFFITSi   hadi potential residual    hii
## 1      1  0.0727 0.0002  0.0384 0.2801    0.2781   0.0019 0.2176
## 2      2 -0.1015 0.0001 -0.0285 0.0829    0.0791   0.0038 0.0733
## 3      3  0.4542 0.0147  0.3148 0.5560    0.4803   0.0757 0.3244
## 4      4 -1.1958 0.0643 -0.6780 0.8388    0.3214   0.5174 0.2432
## 5      5  0.0180 0.0000  0.0074 0.1691    0.1690   0.0001 0.1446
## 6      6  0.1147 0.0011  0.0844 0.5465    0.5416   0.0048 0.3513
## 7      7  0.4124 0.0064  0.2067 0.3138    0.2512   0.0626 0.2008
## 8      8 -1.1627 0.0526 -0.6118 0.7674    0.2769   0.4905 0.2169
## 9      9  0.0440 0.0001  0.0244 0.3092    0.3085   0.0007 0.2357
## 10    10 -1.1353 0.0414 -0.5421 0.6969    0.2280   0.4690 0.1857
## 11    11  1.5323 0.0311  0.4820 0.9545    0.0990   0.8556 0.0900
## 12    12 -0.2834 0.0018 -0.1083 0.1757    0.1462   0.0296 0.1275
## 13    13  0.5933 0.0057  0.1969 0.2396    0.1102   0.1294 0.0993
## 14    14  2.4734 0.0719  0.7949 2.2913    0.1033   2.1880 0.0936
## 15    15  1.4984 0.1172  0.9335 1.1888    0.3881   0.8007 0.2796
## 16    16 -1.7132 0.2938 -1.5018 1.7818    0.7684   1.0134 0.4345
## 17    17  0.8016 0.0407  0.5290 0.6698    0.4355   0.2343 0.3034
## 18    18 -0.5776 0.1381 -0.9667 2.9226    2.8013   0.1213 0.7369
## 19    19  0.7906 0.0609  0.6465 0.8960    0.6687   0.2273 0.4007
## 20    20 -1.2510 0.0906 -0.8077 0.9798    0.4169   0.5630 0.2942
## 21    21 -1.2417 0.4191 -1.7359 2.4936    1.9546   0.5391 0.6615
## 22    22 -0.2311 0.0015 -0.1001 0.2073    0.1876   0.0197 0.1580
## 23    23  0.8029 0.0404  0.5271 0.6661    0.4309   0.2351 0.3012
## 24    24  0.8699 0.0374  0.5088 0.6181    0.3421   0.2760 0.2549
## 25    25 -0.3433 0.0062 -0.2031 0.3933    0.3500   0.0433 0.2593
## 26    26 -0.0505 0.0001 -0.0208 0.1709    0.1700   0.0009 0.1453
## 27    27 -2.0496 0.1033 -0.9159 1.6924    0.1997   1.4927 0.1664

Plot ri, ci, Hi, DFFITSi, dan Potential Residual

#Plot ri
ggplot(hasil) +
  geom_point(aes(y = `ri`, x = `index`),
             color="blue",size=5) +
  ylab("ri") +
  xlab("index") +
  ggtitle("Plot ri") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))

#Plot ci
ggplot(hasil) +
  geom_point(aes(y = `ci`, x = `index`),
             color="red",size=4) +
  ylab("ci") +
  xlab("index") +
  ggtitle("Plot ci") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))

#plot DFFITSi
ggplot(hasil) +
  geom_point(aes(y = `DFFITSi`, x = `index`),
             color="green",size=3) +
  ylab("DFFITSi") +
  xlab("index") +
  ggtitle("Plot DFFITSi") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))

#plot Hi
ggplot(hasil) +
  geom_point(aes(y = `hadi`, x = `index`),
             color="navyblue",size=2.5) +
  ylab("Hi") +
  xlab("index") +
  ggtitle("Plot Hi") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))

#plot Potensial residual
ggplot(hasil) +
  geom_point(aes(y = `potential`, x = `residual`),
             color="yellow",size=5.5) +
  ylab("potential") +
  xlab("residual") +
  ggtitle("P-R Plot") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))

Pencilan atau Outlier

for (i in 1:dim(hasil)[1]){
  absri <- abs(hasil$ri)
  pencilan <- which(absri > 2)
}
pencilan
## [1] 14 27

Berdasarkan data tersebut,terdapat dua titik observasi yang terdeteksi sebagai sebuah pencilan, yaitu amatan ke-14 dengan nilai rinya 2.4734 dan amatan ke-27 dengan nilai rinya -2.0496.

Titik Leverage

titik_leverage <- vector("list", dim(hasil)[1])
for (i in 1:dim(hasil)[1]) {
  cutoff <- 2 * 7/ 27
  titik_leverage[[i]] <- which(hii > cutoff)
}
leverages <- unlist(titik_leverage)
titik_leverage <- sort(unique(leverages))
titik_leverage
## [1] 18 21

Berdasarkan data tersebut,terdapat dua titik observasi yang termasuk sebagai titik leverage pada data , yaitu amatan ke-18 dengan nilai hii 0.7369 dan amatan ke-21 dengan nilai hii 0.1580.

Amatan Berpengaruh

amatan_berpengaruh <- vector("list", dim(hasil)[1])
for (i in 1:dim(hasil)[1]) {
  cutoff <- 2 * sqrt((7 / 22))
  amatan_berpengaruh[[i]] <- which(abs(DFFITSi) > cutoff)
}
berpengaruh <- unlist(amatan_berpengaruh)
amatan_berpengaruh <- sort(unique(berpengaruh))
amatan_berpengaruh
## [1] 16 21

Berdasarkan data tersebut menunjukan bahwa amatan ke-16 dan ke-21 termasuk sebagai amatan berpengaruh dengan nilai DFFITSI-nya -1.5018 dan -1.7359.

ols_plot_resid_lev(model2)

#pembentukan data baru tanpa amatan 14,18,dan 27
umk2 <- subset(umk1, select = -c(X2,X3,X4,X9))
umk2 <- umk2[-c(14,27,18),]
#pemodelan rlb terbaru
model3 <- lm(Y~.,data= umk2)
summary(model3)
## 
## Call:
## lm(formula = Y ~ ., data = umk2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.63474 -0.21423 -0.04211  0.22633  0.81547 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.557e+01  5.447e+00  -4.694 0.000209 ***
## X1           9.466e+00  2.384e+00   3.971 0.000988 ***
## X5           6.942e-05  3.282e-05   2.115 0.049452 *  
## X6           1.865e-02  1.123e-02   1.660 0.115178    
## X7           3.717e-02  4.236e-02   0.878 0.392391    
## X8           1.699e-01  4.866e-02   3.492 0.002793 ** 
## X10         -3.645e-03  4.713e-03  -0.773 0.449868    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4236 on 17 degrees of freedom
## Multiple R-squared:  0.8816, Adjusted R-squared:  0.8398 
## F-statistic:  21.1 on 6 and 17 DF,  p-value: 5.302e-07
anova(model3)
## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value   Pr(>F)    
## X1         1 16.3992 16.3992 91.4002 2.98e-08 ***
## X5         1  3.5049  3.5049 19.5343 0.000375 ***
## X6         1  0.0775  0.0775  0.4318 0.519925    
## X7         1  0.3789  0.3789  2.1120 0.164357    
## X8         1  2.2417  2.2417 12.4942 0.002545 ** 
## X10        1  0.1073  0.1073  0.5983 0.449868    
## Residuals 17  3.0502  0.1794                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pemilihan Model Terbaik

Forward

step(model3,direction="forward")
## Start:  AIC=-35.51
## Y ~ X1 + X5 + X6 + X7 + X8 + X10
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X7 + X8 + X10, data = umk2)
## 
## Coefficients:
## (Intercept)           X1           X5           X6           X7           X8  
##  -2.557e+01    9.466e+00    6.942e-05    1.865e-02    3.717e-02    1.699e-01  
##         X10  
##  -3.645e-03

Backward

step(model3,direction="backward")
## Start:  AIC=-35.51
## Y ~ X1 + X5 + X6 + X7 + X8 + X10
## 
##        Df Sum of Sq    RSS     AIC
## - X10   1   0.10734 3.1575 -36.678
## - X7    1   0.13819 3.1884 -36.445
## <none>              3.0502 -35.509
## - X6    1   0.49462 3.5448 -33.902
## - X5    1   0.80297 3.8531 -31.900
## - X8    1   2.18770 5.2379 -24.531
## - X1    1   2.82861 5.8788 -21.761
## 
## Step:  AIC=-36.68
## Y ~ X1 + X5 + X6 + X7 + X8
## 
##        Df Sum of Sq    RSS     AIC
## - X7    1   0.09057 3.2481 -38.000
## <none>              3.1575 -36.678
## - X6    1   0.42332 3.5808 -35.659
## - X5    1   0.70762 3.8651 -33.825
## - X8    1   2.24173 5.3992 -25.803
## - X1    1   2.83231 5.9898 -23.312
## 
## Step:  AIC=-38
## Y ~ X1 + X5 + X6 + X8
## 
##        Df Sum of Sq    RSS     AIC
## <none>              3.2481 -38.000
## - X6    1    0.5198 3.7679 -36.437
## - X5    1    0.6373 3.8854 -35.700
## - X8    1    2.5301 5.7782 -26.175
## - X1    1    3.2949 6.5430 -23.192
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X8, data = umk2)
## 
## Coefficients:
## (Intercept)           X1           X5           X6           X8  
##  -2.322e+01    7.909e+00    5.881e-05    1.838e-02    1.788e-01

Stepwise

step(model3,direction="both")
## Start:  AIC=-35.51
## Y ~ X1 + X5 + X6 + X7 + X8 + X10
## 
##        Df Sum of Sq    RSS     AIC
## - X10   1   0.10734 3.1575 -36.678
## - X7    1   0.13819 3.1884 -36.445
## <none>              3.0502 -35.509
## - X6    1   0.49462 3.5448 -33.902
## - X5    1   0.80297 3.8531 -31.900
## - X8    1   2.18770 5.2379 -24.531
## - X1    1   2.82861 5.8788 -21.761
## 
## Step:  AIC=-36.68
## Y ~ X1 + X5 + X6 + X7 + X8
## 
##        Df Sum of Sq    RSS     AIC
## - X7    1   0.09057 3.2481 -38.000
## <none>              3.1575 -36.678
## - X6    1   0.42332 3.5808 -35.659
## + X10   1   0.10734 3.0502 -35.509
## - X5    1   0.70762 3.8651 -33.825
## - X8    1   2.24173 5.3992 -25.803
## - X1    1   2.83231 5.9898 -23.312
## 
## Step:  AIC=-38
## Y ~ X1 + X5 + X6 + X8
## 
##        Df Sum of Sq    RSS     AIC
## <none>              3.2481 -38.000
## + X7    1    0.0906 3.1575 -36.678
## + X10   1    0.0597 3.1884 -36.445
## - X6    1    0.5198 3.7679 -36.437
## - X5    1    0.6373 3.8854 -35.700
## - X8    1    2.5301 5.7782 -26.175
## - X1    1    3.2949 6.5430 -23.192
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X8, data = umk2)
## 
## Coefficients:
## (Intercept)           X1           X5           X6           X8  
##  -2.322e+01    7.909e+00    5.881e-05    1.838e-02    1.788e-01

Model terbaik

best <- ols_step_best_subset(model3)
best
##      Best Subsets Regression     
## ---------------------------------
## Model Index    Predictors
## ---------------------------------
##      1         X1                 
##      2         X1 X8              
##      3         X1 X5 X8           
##      4         X1 X5 X6 X8        
##      5         X1 X5 X6 X7 X8     
##      6         X1 X5 X6 X7 X8 X10 
## ---------------------------------
## 
##                                                   Subsets Regression Summary                                                   
## -------------------------------------------------------------------------------------------------------------------------------
##                        Adj.        Pred                                                                                         
## Model    R-Square    R-Square    R-Square     C(p)        AIC        SBIC        SBC       MSEP       FPE       HSP       APC  
## -------------------------------------------------------------------------------------------------------------------------------
##   1        0.6366      0.6201      0.5692    32.1704    51.5118    -19.3403    55.0460    10.2149    0.4609    0.0203    0.4294 
##   2        0.8047      0.7861      0.7467    10.0342    38.6056    -30.4083    43.3178     5.7635    0.2695    0.0120    0.2511 
##   3        0.8537      0.8318       0.772     5.0000    33.6719    -33.3351    39.5622     4.5446    0.2198    0.0099    0.2048 
##   4        0.8739      0.8474       0.785     4.1030    32.1093    -32.9545    39.1776     4.1353    0.2066    0.0095    0.1925 
##   5        0.8774      0.8434      0.7678     5.5983    33.4306    -30.5779    41.6770     4.2565    0.2193    0.0103    0.2043 
##   6        0.8816      0.8398      0.7492     7.0000    34.6005    -28.0829    44.0250     4.3687    0.2318    0.0112    0.2159 
## -------------------------------------------------------------------------------------------------------------------------------
## 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

Model Regresi Linear Berganda Terbaik

model4<- lm(Y~X1+X5+X6+X8, data=umk2)
summary(model4)
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X8, data = umk2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.66376 -0.23229 -0.01766  0.29315  0.74941 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.322e+01  4.699e+00  -4.943 9.03e-05 ***
## X1           7.909e+00  1.801e+00   4.390 0.000315 ***
## X5           5.881e-05  3.046e-05   1.931 0.068572 .  
## X6           1.838e-02  1.054e-02   1.744 0.097372 .  
## X8           1.788e-01  4.647e-02   3.847 0.001086 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4135 on 19 degrees of freedom
## Multiple R-squared:  0.8739, Adjusted R-squared:  0.8474 
## F-statistic: 32.92 on 4 and 19 DF,  p-value: 2.661e-08

Eksplorasi Kondisi Gauss Markov

Plot Sisaan vs Y duga

plot(model4,1)

  1. Sisaan menyebar di sekitar 0, sehingga nilai harapan galat sama dengan nol
  2. Lebar pita sama untuk setiap nilai dugaan, sehingga ragam homogen

Plot Sisaan vs Urutan

plot(x = 1:dim(umk2)[1],
     y = model4$residuals,
     type = 'b', 
     ylab = "Residuals",
     xlab = "Observation")

Tebaran tidak berpola, sehingga sisaan saling bebas, dan model pas

Normalitas Sisaan dengan QQ-Plot

plot(model4,2)

Titik-titik data pada QQ-plot cukup dekat dengan garis.Hal ini menunjukkan bahwa distribusi sisaan mirip dengan Sebaran normal.

Uji Formal Kondisi Gauss-Markov

1. Nilai harapan sisaan sama dengan nol

H0: Nilai harapan sisaan sama dengan nol
H1: Nilai harapan sisaan tidak sama dengan nol

t.test(model4$residuals,mu = 0,conf.level = 0.95)
## 
##  One Sample t-test
## 
## data:  model4$residuals
## t = 1.4887e-16, df = 23, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.1586838  0.1586838
## sample estimates:
##    mean of x 
## 1.141942e-17

Uji t-tes tersebut menunjukkan hasil p-value > alpha = 0.05, maka tak tolak H0,sehingg nilai harapan sisaan sama dengan nol pada taraf nyata 5%. Asumsi terpenuhi.

2. Ragam Sisaan Homogen

H0:ragam sisan homogen
H1:ragam siaan tidak homogen

bptest(model4)
## 
##  studentized Breusch-Pagan test
## 
## data:  model4
## BP = 3.622, df = 4, p-value = 0.4596

Uji ini sering disebut dengan uji homokesdasitas.Karena p-value = 0.4596> alpha = 0.05, maka tak tolak H0, ragam sisaan homogen pada taraf nyata 5%. Asumsi terpenuhi.

3. Sisaan Saling Bebas

H0:sisaan saling bebas/tidak ada autokorelasi
H1:sisaan tidak saling bebas/ada autokorelasi

dwtest(model4)
## 
##  Durbin-Watson test
## 
## data:  model4
## DW = 2.5461, p-value = 0.8733
## alternative hypothesis: true autocorrelation is greater than 0

Uji ini sering disebut dengan uji autokorelasi yang dilakukan dengan Durbin watson. Berdasarakan Uji Durbin-Watson yang dilakukan, menghasilkan nilai p-value = 0.8733 > alpha = 0.05, maka tak tolak H0, sisaan saling bebas pada taraf nyata 5%, sehingga asumsi terpenuhi.

Uji Formal Normalitas Sisaan

H0:sisaan menyebar Normal
H1:sisaan tidak menyebar Normal

shapiro.test(model4$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  model4$residuals
## W = 0.96793, p-value = 0.6163

Berdasarkan Uji Shapiro-test yang dilakukan mengahasilkan nilai p-value = 0.6163 > alpha = 0.05, maka tak tolak H0, sehingga sisaan menyebar normal pada taraf nyata 5%.

Uji Kelayakan Model

summary(model4)
## 
## Call:
## lm(formula = Y ~ X1 + X5 + X6 + X8, data = umk2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.66376 -0.23229 -0.01766  0.29315  0.74941 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.322e+01  4.699e+00  -4.943 9.03e-05 ***
## X1           7.909e+00  1.801e+00   4.390 0.000315 ***
## X5           5.881e-05  3.046e-05   1.931 0.068572 .  
## X6           1.838e-02  1.054e-02   1.744 0.097372 .  
## X8           1.788e-01  4.647e-02   3.847 0.001086 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4135 on 19 degrees of freedom
## Multiple R-squared:  0.8739, Adjusted R-squared:  0.8474 
## F-statistic: 32.92 on 4 and 19 DF,  p-value: 2.661e-08
anova(model4)
## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value    Pr(>F)    
## X1         1 16.3992 16.3992 95.9288 7.348e-09 ***
## X5         1  3.5049  3.5049 20.5021 0.0002301 ***
## X6         1  0.0775  0.0775  0.4532 0.5089486    
## X8         1  2.5301  2.5301 14.8001 0.0010865 ** 
## Residuals 19  3.2481  0.1710                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1