library(readxl)
## Warning: package 'readxl' was built under R version 4.2.3
data <- read_xlsx("C:\\Users\\Sherinda Engelina\\Downloads\\datafix.xlsx")
y<-data$y
xpadi<-data$xpadi
xpdrb<-data$xpdrb
xmiskin<-data$xmiskin
xpenduduk<-data$xpenduduk
xpengeluaran<-data$xpengeluaran
xuhh<-data$xuhh
xtpt<-data$xtpt
xair<-data$xair
datareg<-data.frame(cbind(y,xpadi,xpdrb,xmiskin,xpenduduk,xpengeluaran,xuhh,xtpt,xair))
head(datareg[24,])
## y xpadi xpdrb xmiskin xpenduduk xpengeluaran xuhh xtpt xair
## 24 35.63 4.12 1.110182 26.03 68 86.60588 66.44 2.67 66.49
library(PerformanceAnalytics)
## Warning: package 'PerformanceAnalytics' was built under R version 4.2.3
## Loading required package: xts
## Warning: package 'xts' was built under R version 4.2.3
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.2.3
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
##
## legend
chart.Correlation(datareg,
histogram = TRUE, # tambahkan histogram pada diagonal
pch = 1) # jenis point = 1
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
y <- data.matrix(data[, 1])
xpadi <- data.matrix(data[, 2])
xpdrb <- data.matrix(data[, 3])
xmiskin <- data.matrix(data[, 4])
xpenduduk <- data.matrix(data[,5])
xpengeluaran <- data.matrix(data[,6])
xuhh <- data.matrix(data[,7])
xtpt <- data.matrix(data[, 8])
xair <- data.matrix(data[, 9])
x <- cbind(xpadi,xpdrb,xmiskin,xpenduduk,xpengeluaran,xuhh,xtpt,xair)
modelawal <- lm(y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +xtpt + xair, data = datareg)
summary(modelawal)
##
## Call:
## lm(formula = y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran +
## xuhh + xtpt + xair, data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.4303 -2.6843 0.0079 1.9858 15.9264
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 27.4657911 41.4553636 0.663 0.513688
## xpadi -5.1614106 1.4457335 -3.570 0.001480 **
## xpdrb 0.0611397 0.1177828 0.519 0.608264
## xmiskin 0.9316159 0.2488774 3.743 0.000955 ***
## xpenduduk -0.0002914 0.0004253 -0.685 0.499481
## xpengeluaran 0.0789229 0.1257942 0.627 0.536088
## xuhh -0.1218253 0.6145677 -0.198 0.844469
## xtpt 0.0220401 0.9251916 0.024 0.981184
## xair 0.0117572 0.1660051 0.071 0.944101
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.462 on 25 degrees of freedom
## Multiple R-squared: 0.6294, Adjusted R-squared: 0.5108
## F-statistic: 5.307 on 8 and 25 DF, p-value: 0.0005887
model.ols<-modelawal
model.ols.summary<-summary(modelawal)
anova(modelawal)
## Analysis of Variance Table
##
## Response: y
## Df Sum Sq Mean Sq F value Pr(>F)
## xpadi 1 501.63 501.63 16.8141 0.0003826 ***
## xpdrb 1 4.69 4.69 0.1572 0.6950823
## xmiskin 1 739.44 739.44 24.7852 3.941e-05 ***
## xpenduduk 1 4.67 4.67 0.1566 0.6956432
## xpengeluaran 1 14.83 14.83 0.4971 0.4872805
## xuhh 1 1.08 1.08 0.0363 0.8504669
## xtpt 1 0.07 0.07 0.0024 0.9615686
## xair 1 0.15 0.15 0.0050 0.9441009
## Residuals 25 745.85 29.83
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Berdasarkan pemodelan dan uji anova diperoleh peubah yang signifikan adalah peubah produktivitas padi dan peubah kemiskinan penduduk
plot(modelawal,1)
plot(x = 1:dim(data)[1],
y = modelawal$residuals,
type = 'b',
ylab = "Residuals",
xlab = "Observation")
plot(modelawal,2)
residu_data=resid(modelawal)
library(car)
## Warning: package 'car' was built under R version 4.2.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.2.3
shapiro.test(modelawal$residuals)
##
## Shapiro-Wilk normality test
##
## data: modelawal$residuals
## W = 0.94307, p-value = 0.07612
H0: menyebar normal ; H1: Tidak menyebar normal
P-value (0.076) > Alpha (0.05),
terima H0, berarti data berdistribusi normal
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.2.3
Homoskedasticity_test <- bptest(modelawal)
print(Homoskedasticity_test$p.value)
## BP
## 0.2433766
H0 : Homoskedastik ; H1 : Heterokedastik
P-value (0.243) >Alpha (0.05), terima H0, data homoskedastik
library(randtests)
runs.test(modelawal$residuals)
##
## Runs Test
##
## data: modelawal$residuals
## statistic = 0.69663, runs = 20, n1 = 17, n2 = 17, n = 34, p-value =
## 0.486
## alternative hypothesis: nonrandomness
H0 : Tidak ada serial korelasi (autokorelasi) ; H1 : Terdapat serial korelasi (autokorelasi)
P-Value (0.486) > Alpha (0.05), sehingga terima H0 saling bebas (tidak autokorelasi)
t.test(modelawal$residuals,mu = 0,conf.level = 0.95)
##
## One Sample t-test
##
## data: modelawal$residuals
## t = -8.01e-18, df = 33, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -1.658781 1.658781
## sample estimates:
## mean of x
## -6.530724e-18
library(olsrr)
## Warning: package 'olsrr' was built under R version 4.2.3
##
## Attaching package: 'olsrr'
## The following object is masked from 'package:datasets':
##
## rivers
ols_vif_tol(modelawal)
## Variables Tolerance VIF
## 1 xpadi 0.6118596 1.634362
## 2 xpdrb 0.4739731 2.109824
## 3 xmiskin 0.5432208 1.840872
## 4 xpenduduk 0.5844683 1.710957
## 5 xpengeluaran 0.4211194 2.374624
## 6 xuhh 0.4138111 2.416562
## 7 xtpt 0.5244779 1.906658
## 8 xair 0.5457675 1.832282
VIF<10 tidak terdapat multikolinearitas antar peubah
step(lm(y~.,data=datareg),direction="backward")
## Start: AIC=123
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xtpt + xair
##
## Df Sum of Sq RSS AIC
## - xtpt 1 0.02 745.86 121.00
## - xair 1 0.15 746.00 121.00
## - xuhh 1 1.17 747.02 121.05
## - xpdrb 1 8.04 753.88 121.36
## - xpengeluaran 1 11.74 757.59 121.53
## - xpenduduk 1 14.01 759.86 121.63
## <none> 745.85 123.00
## - xpadi 1 380.25 1126.10 135.00
## - xmiskin 1 418.03 1163.88 136.13
##
## Step: AIC=121
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xair
##
## Df Sum of Sq RSS AIC
## - xair 1 0.20 746.07 119.01
## - xuhh 1 1.27 747.14 119.06
## - xpdrb 1 11.13 756.99 119.50
## - xpenduduk 1 14.19 760.06 119.64
## - xpengeluaran 1 15.02 760.88 119.68
## <none> 745.86 121.00
## - xpadi 1 411.71 1157.57 133.94
## - xmiskin 1 427.82 1173.68 134.41
##
## Step: AIC=119.01
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh
##
## Df Sum of Sq RSS AIC
## - xuhh 1 1.08 747.15 117.06
## - xpdrb 1 11.21 757.27 117.51
## - xpenduduk 1 14.87 760.93 117.68
## - xpengeluaran 1 15.81 761.88 117.72
## <none> 746.07 119.01
## - xpadi 1 419.32 1165.39 132.17
## - xmiskin 1 432.36 1178.43 132.55
##
## Step: AIC=117.06
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpdrb 1 10.15 757.30 115.52
## - xpengeluaran 1 14.83 761.98 115.72
## - xpenduduk 1 14.91 762.05 115.73
## <none> 747.15 117.06
## - xpadi 1 433.41 1180.56 130.61
## - xmiskin 1 682.03 1429.17 137.11
##
## Step: AIC=115.52
## y ~ xpadi + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpenduduk 1 13.21 770.50 114.10
## - xpengeluaran 1 13.56 770.85 114.12
## <none> 757.30 115.52
## - xpadi 1 485.73 1243.03 130.36
## - xmiskin 1 673.92 1431.22 135.16
##
## Step: AIC=114.1
## y ~ xpadi + xmiskin + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpengeluaran 1 4.37 774.87 112.30
## <none> 770.50 114.10
## - xpadi 1 569.66 1340.16 130.92
## - xmiskin 1 683.74 1454.25 133.70
##
## Step: AIC=112.3
## y ~ xpadi + xmiskin
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## - xpadi 1 598.40 1373.27 129.75
## - xmiskin 1 735.91 1510.78 133.00
##
## Call:
## lm(formula = y ~ xpadi + xmiskin, data = datareg)
##
## Coefficients:
## (Intercept) xpadi xmiskin
## 25.7220 -5.0804 0.9139
bestmodel<-step(lm(y~.,data=datareg),direction="backward")
## Start: AIC=123
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xtpt + xair
##
## Df Sum of Sq RSS AIC
## - xtpt 1 0.02 745.86 121.00
## - xair 1 0.15 746.00 121.00
## - xuhh 1 1.17 747.02 121.05
## - xpdrb 1 8.04 753.88 121.36
## - xpengeluaran 1 11.74 757.59 121.53
## - xpenduduk 1 14.01 759.86 121.63
## <none> 745.85 123.00
## - xpadi 1 380.25 1126.10 135.00
## - xmiskin 1 418.03 1163.88 136.13
##
## Step: AIC=121
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xair
##
## Df Sum of Sq RSS AIC
## - xair 1 0.20 746.07 119.01
## - xuhh 1 1.27 747.14 119.06
## - xpdrb 1 11.13 756.99 119.50
## - xpenduduk 1 14.19 760.06 119.64
## - xpengeluaran 1 15.02 760.88 119.68
## <none> 745.86 121.00
## - xpadi 1 411.71 1157.57 133.94
## - xmiskin 1 427.82 1173.68 134.41
##
## Step: AIC=119.01
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh
##
## Df Sum of Sq RSS AIC
## - xuhh 1 1.08 747.15 117.06
## - xpdrb 1 11.21 757.27 117.51
## - xpenduduk 1 14.87 760.93 117.68
## - xpengeluaran 1 15.81 761.88 117.72
## <none> 746.07 119.01
## - xpadi 1 419.32 1165.39 132.17
## - xmiskin 1 432.36 1178.43 132.55
##
## Step: AIC=117.06
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpdrb 1 10.15 757.30 115.52
## - xpengeluaran 1 14.83 761.98 115.72
## - xpenduduk 1 14.91 762.05 115.73
## <none> 747.15 117.06
## - xpadi 1 433.41 1180.56 130.61
## - xmiskin 1 682.03 1429.17 137.11
##
## Step: AIC=115.52
## y ~ xpadi + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpenduduk 1 13.21 770.50 114.10
## - xpengeluaran 1 13.56 770.85 114.12
## <none> 757.30 115.52
## - xpadi 1 485.73 1243.03 130.36
## - xmiskin 1 673.92 1431.22 135.16
##
## Step: AIC=114.1
## y ~ xpadi + xmiskin + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpengeluaran 1 4.37 774.87 112.30
## <none> 770.50 114.10
## - xpadi 1 569.66 1340.16 130.92
## - xmiskin 1 683.74 1454.25 133.70
##
## Step: AIC=112.3
## y ~ xpadi + xmiskin
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## - xpadi 1 598.40 1373.27 129.75
## - xmiskin 1 735.91 1510.78 133.00
summary(bestmodel)
##
## Call:
## lm(formula = y ~ xpadi + xmiskin, data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.0285 -3.2965 -0.4622 2.2292 15.5129
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 25.7220 5.0318 5.112 1.56e-05 ***
## xpadi -5.0804 1.0383 -4.893 2.92e-05 ***
## xmiskin 0.9139 0.1684 5.426 6.33e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5 on 31 degrees of freedom
## Multiple R-squared: 0.615, Adjusted R-squared: 0.5901
## F-statistic: 24.75 on 2 and 31 DF, p-value: 3.763e-07
step(lm(y~1,data=datareg),direction="forward",scope=~xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +xtpt + xair)
## Start: AIC=140.74
## y ~ 1
##
## Df Sum of Sq RSS AIC
## + xmiskin 1 639.14 1373.3 129.75
## + xpadi 1 501.63 1510.8 133.00
## + xuhh 1 450.89 1561.5 134.12
## + xair 1 286.79 1725.6 137.52
## + xpdrb 1 172.24 1840.2 139.70
## + xpenduduk 1 123.83 1888.6 140.59
## <none> 2012.4 140.75
## + xtpt 1 43.57 1968.8 142.00
## + xpengeluaran 1 13.24 1999.2 142.52
##
## Step: AIC=129.75
## y ~ xmiskin
##
## Df Sum of Sq RSS AIC
## + xpadi 1 598.40 774.87 112.30
## + xpdrb 1 103.87 1269.40 129.08
## + xair 1 103.85 1269.42 129.08
## <none> 1373.27 129.75
## + xuhh 1 53.40 1319.87 130.40
## + xpenduduk 1 33.72 1339.55 130.91
## + xpengeluaran 1 33.11 1340.16 130.92
## + xtpt 1 5.69 1367.58 131.61
##
## Step: AIC=112.3
## y ~ xmiskin + xpadi
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## + xpdrb 1 8.2223 766.65 113.93
## + xtpt 1 5.0614 769.81 114.07
## + xpengeluaran 1 4.3704 770.50 114.10
## + xpenduduk 1 4.0199 770.85 114.12
## + xair 1 1.0376 773.84 114.25
## + xuhh 1 0.2075 774.67 114.29
##
## Call:
## lm(formula = y ~ xmiskin + xpadi, data = datareg)
##
## Coefficients:
## (Intercept) xmiskin xpadi
## 25.7220 0.9139 -5.0804
bestmodel<-step(lm(y~1,data=datareg),direction="forward",scope=~xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +xtpt + xair)
## Start: AIC=140.74
## y ~ 1
##
## Df Sum of Sq RSS AIC
## + xmiskin 1 639.14 1373.3 129.75
## + xpadi 1 501.63 1510.8 133.00
## + xuhh 1 450.89 1561.5 134.12
## + xair 1 286.79 1725.6 137.52
## + xpdrb 1 172.24 1840.2 139.70
## + xpenduduk 1 123.83 1888.6 140.59
## <none> 2012.4 140.75
## + xtpt 1 43.57 1968.8 142.00
## + xpengeluaran 1 13.24 1999.2 142.52
##
## Step: AIC=129.75
## y ~ xmiskin
##
## Df Sum of Sq RSS AIC
## + xpadi 1 598.40 774.87 112.30
## + xpdrb 1 103.87 1269.40 129.08
## + xair 1 103.85 1269.42 129.08
## <none> 1373.27 129.75
## + xuhh 1 53.40 1319.87 130.40
## + xpenduduk 1 33.72 1339.55 130.91
## + xpengeluaran 1 33.11 1340.16 130.92
## + xtpt 1 5.69 1367.58 131.61
##
## Step: AIC=112.3
## y ~ xmiskin + xpadi
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## + xpdrb 1 8.2223 766.65 113.93
## + xtpt 1 5.0614 769.81 114.07
## + xpengeluaran 1 4.3704 770.50 114.10
## + xpenduduk 1 4.0199 770.85 114.12
## + xair 1 1.0376 773.84 114.25
## + xuhh 1 0.2075 774.67 114.29
summary(bestmodel)
##
## Call:
## lm(formula = y ~ xmiskin + xpadi, data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.0285 -3.2965 -0.4622 2.2292 15.5129
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 25.7220 5.0318 5.112 1.56e-05 ***
## xmiskin 0.9139 0.1684 5.426 6.33e-06 ***
## xpadi -5.0804 1.0383 -4.893 2.92e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5 on 31 degrees of freedom
## Multiple R-squared: 0.615, Adjusted R-squared: 0.5901
## F-statistic: 24.75 on 2 and 31 DF, p-value: 3.763e-07
step(lm(y~.,data=datareg),direction="both")
## Start: AIC=123
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xtpt + xair
##
## Df Sum of Sq RSS AIC
## - xtpt 1 0.02 745.86 121.00
## - xair 1 0.15 746.00 121.00
## - xuhh 1 1.17 747.02 121.05
## - xpdrb 1 8.04 753.88 121.36
## - xpengeluaran 1 11.74 757.59 121.53
## - xpenduduk 1 14.01 759.86 121.63
## <none> 745.85 123.00
## - xpadi 1 380.25 1126.10 135.00
## - xmiskin 1 418.03 1163.88 136.13
##
## Step: AIC=121
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xair
##
## Df Sum of Sq RSS AIC
## - xair 1 0.20 746.07 119.01
## - xuhh 1 1.27 747.14 119.06
## - xpdrb 1 11.13 756.99 119.50
## - xpenduduk 1 14.19 760.06 119.64
## - xpengeluaran 1 15.02 760.88 119.68
## <none> 745.86 121.00
## + xtpt 1 0.02 745.85 123.00
## - xpadi 1 411.71 1157.57 133.94
## - xmiskin 1 427.82 1173.68 134.41
##
## Step: AIC=119.01
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh
##
## Df Sum of Sq RSS AIC
## - xuhh 1 1.08 747.15 117.06
## - xpdrb 1 11.21 757.27 117.51
## - xpenduduk 1 14.87 760.93 117.68
## - xpengeluaran 1 15.81 761.88 117.72
## <none> 746.07 119.01
## + xair 1 0.20 745.86 121.00
## + xtpt 1 0.07 746.00 121.00
## - xpadi 1 419.32 1165.39 132.17
## - xmiskin 1 432.36 1178.43 132.55
##
## Step: AIC=117.06
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpdrb 1 10.15 757.30 115.52
## - xpengeluaran 1 14.83 761.98 115.72
## - xpenduduk 1 14.91 762.05 115.73
## <none> 747.15 117.06
## + xuhh 1 1.08 746.07 119.01
## + xtpt 1 0.13 747.02 119.05
## + xair 1 0.01 747.14 119.06
## - xpadi 1 433.41 1180.56 130.61
## - xmiskin 1 682.03 1429.17 137.11
##
## Step: AIC=115.52
## y ~ xpadi + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpenduduk 1 13.21 770.50 114.10
## - xpengeluaran 1 13.56 770.85 114.12
## <none> 757.30 115.52
## + xpdrb 1 10.15 747.15 117.06
## + xtpt 1 3.41 753.88 117.36
## + xair 1 0.30 756.99 117.50
## + xuhh 1 0.02 757.27 117.51
## - xpadi 1 485.73 1243.03 130.36
## - xmiskin 1 673.92 1431.22 135.16
##
## Step: AIC=114.1
## y ~ xpadi + xmiskin + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpengeluaran 1 4.37 774.87 112.30
## <none> 770.50 114.10
## + xpenduduk 1 13.21 757.30 115.52
## + xpdrb 1 8.45 762.05 115.73
## + xtpt 1 2.60 767.90 115.99
## + xair 1 0.54 769.96 116.08
## + xuhh 1 0.00 770.50 116.10
## - xpadi 1 569.66 1340.16 130.92
## - xmiskin 1 683.74 1454.25 133.70
##
## Step: AIC=112.3
## y ~ xpadi + xmiskin
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## + xpdrb 1 8.22 766.65 113.93
## + xtpt 1 5.06 769.81 114.07
## + xpengeluaran 1 4.37 770.50 114.10
## + xpenduduk 1 4.02 770.85 114.12
## + xair 1 1.04 773.84 114.25
## + xuhh 1 0.21 774.67 114.29
## - xpadi 1 598.40 1373.27 129.75
## - xmiskin 1 735.91 1510.78 133.00
##
## Call:
## lm(formula = y ~ xpadi + xmiskin, data = datareg)
##
## Coefficients:
## (Intercept) xpadi xmiskin
## 25.7220 -5.0804 0.9139
bestmodel<-step(lm(y~.,data=datareg),direction="both")
## Start: AIC=123
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xtpt + xair
##
## Df Sum of Sq RSS AIC
## - xtpt 1 0.02 745.86 121.00
## - xair 1 0.15 746.00 121.00
## - xuhh 1 1.17 747.02 121.05
## - xpdrb 1 8.04 753.88 121.36
## - xpengeluaran 1 11.74 757.59 121.53
## - xpenduduk 1 14.01 759.86 121.63
## <none> 745.85 123.00
## - xpadi 1 380.25 1126.10 135.00
## - xmiskin 1 418.03 1163.88 136.13
##
## Step: AIC=121
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh +
## xair
##
## Df Sum of Sq RSS AIC
## - xair 1 0.20 746.07 119.01
## - xuhh 1 1.27 747.14 119.06
## - xpdrb 1 11.13 756.99 119.50
## - xpenduduk 1 14.19 760.06 119.64
## - xpengeluaran 1 15.02 760.88 119.68
## <none> 745.86 121.00
## + xtpt 1 0.02 745.85 123.00
## - xpadi 1 411.71 1157.57 133.94
## - xmiskin 1 427.82 1173.68 134.41
##
## Step: AIC=119.01
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran + xuhh
##
## Df Sum of Sq RSS AIC
## - xuhh 1 1.08 747.15 117.06
## - xpdrb 1 11.21 757.27 117.51
## - xpenduduk 1 14.87 760.93 117.68
## - xpengeluaran 1 15.81 761.88 117.72
## <none> 746.07 119.01
## + xair 1 0.20 745.86 121.00
## + xtpt 1 0.07 746.00 121.00
## - xpadi 1 419.32 1165.39 132.17
## - xmiskin 1 432.36 1178.43 132.55
##
## Step: AIC=117.06
## y ~ xpadi + xpdrb + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpdrb 1 10.15 757.30 115.52
## - xpengeluaran 1 14.83 761.98 115.72
## - xpenduduk 1 14.91 762.05 115.73
## <none> 747.15 117.06
## + xuhh 1 1.08 746.07 119.01
## + xtpt 1 0.13 747.02 119.05
## + xair 1 0.01 747.14 119.06
## - xpadi 1 433.41 1180.56 130.61
## - xmiskin 1 682.03 1429.17 137.11
##
## Step: AIC=115.52
## y ~ xpadi + xmiskin + xpenduduk + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpenduduk 1 13.21 770.50 114.10
## - xpengeluaran 1 13.56 770.85 114.12
## <none> 757.30 115.52
## + xpdrb 1 10.15 747.15 117.06
## + xtpt 1 3.41 753.88 117.36
## + xair 1 0.30 756.99 117.50
## + xuhh 1 0.02 757.27 117.51
## - xpadi 1 485.73 1243.03 130.36
## - xmiskin 1 673.92 1431.22 135.16
##
## Step: AIC=114.1
## y ~ xpadi + xmiskin + xpengeluaran
##
## Df Sum of Sq RSS AIC
## - xpengeluaran 1 4.37 774.87 112.30
## <none> 770.50 114.10
## + xpenduduk 1 13.21 757.30 115.52
## + xpdrb 1 8.45 762.05 115.73
## + xtpt 1 2.60 767.90 115.99
## + xair 1 0.54 769.96 116.08
## + xuhh 1 0.00 770.50 116.10
## - xpadi 1 569.66 1340.16 130.92
## - xmiskin 1 683.74 1454.25 133.70
##
## Step: AIC=112.3
## y ~ xpadi + xmiskin
##
## Df Sum of Sq RSS AIC
## <none> 774.87 112.30
## + xpdrb 1 8.22 766.65 113.93
## + xtpt 1 5.06 769.81 114.07
## + xpengeluaran 1 4.37 770.50 114.10
## + xpenduduk 1 4.02 770.85 114.12
## + xair 1 1.04 773.84 114.25
## + xuhh 1 0.21 774.67 114.29
## - xpadi 1 598.40 1373.27 129.75
## - xmiskin 1 735.91 1510.78 133.00
summary(bestmodel)
##
## Call:
## lm(formula = y ~ xpadi + xmiskin, data = datareg)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.0285 -3.2965 -0.4622 2.2292 15.5129
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 25.7220 5.0318 5.112 1.56e-05 ***
## xpadi -5.0804 1.0383 -4.893 2.92e-05 ***
## xmiskin 0.9139 0.1684 5.426 6.33e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5 on 31 degrees of freedom
## Multiple R-squared: 0.615, Adjusted R-squared: 0.5901
## F-statistic: 24.75 on 2 and 31 DF, p-value: 3.763e-07
model.best.summary<-summary(bestmodel)
which(abs(rstandard(bestmodel))>2)
## 21
## 21
leverage <- hatvalues(bestmodel)
batas_leverage <- 2 * (ncol(datareg) + 1) / nrow(datareg)
observasi_high_leverage <- which(leverage > batas_leverage)
datareg[observasi_high_leverage, ]
## [1] y xpadi xpdrb xmiskin xpenduduk
## [6] xpengeluaran xuhh xtpt xair
## <0 rows> (or 0-length row.names)
leverage <- hatvalues(bestmodel)
print(leverage)
## 1 2 3 4 5 6 7
## 0.08186813 0.18787981 0.07815914 0.05047351 0.03809087 0.07148918 0.06482301
## 8 9 10 11 12 13 14
## 0.03675119 0.09686113 0.06524964 0.08169921 0.13339702 0.07956892 0.12947670
## 15 16 17 18 19 20 21
## 0.06206502 0.08246267 0.06621075 0.18187806 0.04462592 0.13338498 0.09742539
## 22 23 24 25 26 27 28
## 0.06563466 0.15443942 0.33713350 0.15477973 0.05820114 0.03648457 0.03961020
## 29 30 31 32 33 34
## 0.03551610 0.03917871 0.03981297 0.05452306 0.07422683 0.04661888
Tidak ada amatan yang merupakan titik leverage
ri = rstandard(bestmodel)
hii = hatvalues(bestmodel)
Obs = c(1:34)
summ <- cbind.data.frame(Obs, hii, ri)
library(olsrr)
ols_plot_resid_lev(bestmodel)
data1 <- datareg[-c(0,34),]
data2 <- datareg[-c(0,34),]
model1 <- lm(y~., data=data1)
model2 <- lm(y~., data=data2)
summary(model1)$r.squared
## [1] 0.6282588
summary(model2)$r.squared
## [1] 0.6282588
ri = rstandard(bestmodel)
hii = hatvalues(bestmodel)
Obs = c(1:34)
summ <- cbind.data.frame(Obs, hii, ri)
plot(bestmodel, 5)
Dengan metode jarak Cook (Cook’s distance) secara eksploratif, terlihat bahwa tidak ada amatan berpengaruh pada taraf nyata 5% karena tidak ada titik plot yang berada di luar garis putus-putus skala 0.5.
influence.measures(bestmodel)
## Influence measures of
## lm(formula = y ~ xpadi + xmiskin, data = datareg) :
##
## dfb.1_ dfb.xpad dfb.xmsk dffit cov.r cook.d hat inf
## 1 0.056728 -0.05454 -0.040928 -0.08862 1.191 2.70e-03 0.0819
## 2 -0.320714 0.49253 -0.316680 0.61624 1.158 1.24e-01 0.1879
## 3 0.018476 -0.03614 0.029590 -0.05699 1.192 1.12e-03 0.0782
## 4 0.060889 -0.05122 -0.111197 -0.19571 1.082 1.29e-02 0.0505
## 5 0.000956 -0.00134 -0.000381 -0.00298 1.147 3.06e-06 0.0381
## 6 -0.009813 -0.03672 0.091913 -0.12586 1.164 5.42e-03 0.0715
## 7 -0.027955 0.02197 0.044347 0.06924 1.172 1.65e-03 0.0648
## 8 0.045740 -0.00983 -0.052658 0.12194 1.102 5.06e-03 0.0368
## 9 -0.079460 0.11451 -0.047470 0.14480 1.198 7.18e-03 0.0969
## 10 -0.095952 0.11823 0.005039 0.16071 1.138 8.79e-03 0.0652
## 11 -0.078436 0.09670 -0.003912 0.12097 1.182 5.01e-03 0.0817
## 12 -0.076700 0.06673 0.019759 -0.08074 1.268 2.24e-03 0.1334
## 13 -0.205246 0.11510 0.193015 -0.29359 1.087 2.87e-02 0.0796
## 14 -0.341936 0.27475 0.148941 -0.36817 1.159 4.53e-02 0.1295
## 15 -0.074091 0.04823 0.049026 -0.09878 1.159 3.34e-03 0.0621
## 16 0.148679 -0.11557 -0.063027 0.16999 1.165 9.85e-03 0.0825
## 17 0.157295 -0.05095 -0.218686 0.30734 1.037 3.12e-02 0.0662
## 18 -0.698765 0.61052 0.200733 -0.72036 1.077 1.66e-01 0.1819
## 19 0.012261 -0.01578 -0.003311 -0.02814 1.153 2.73e-04 0.0446
## 20 0.322942 -0.42980 0.409588 0.64750 0.981 1.32e-01 0.1334
## 21 1.184965 -0.96277 -0.431592 1.30337 0.345 3.84e-01 0.0974 *
## 22 0.292948 -0.28735 -0.231511 -0.51707 0.825 8.17e-02 0.0656
## 23 -0.061008 0.23776 -0.657281 -0.75939 0.967 1.80e-01 0.1544
## 24 -0.036090 -0.32017 1.200662 1.27888 1.226 5.09e-01 0.3371 *
## 25 -0.018332 -0.03067 0.186552 0.20805 1.275 1.48e-02 0.1548
## 26 0.025837 -0.01761 -0.014664 0.03392 1.169 3.96e-04 0.0582
## 27 0.028167 -0.03837 -0.022796 -0.10514 1.112 3.77e-03 0.0365
## 28 0.002164 -0.00510 0.003049 -0.01133 1.149 4.42e-05 0.0396
## 29 -0.005763 0.00390 -0.056269 -0.13573 1.089 6.24e-03 0.0355
## 30 -0.084959 0.07824 -0.046078 -0.17603 1.065 1.04e-02 0.0392
## 31 -0.086167 0.03873 0.073175 -0.16769 1.075 9.47e-03 0.0398
## 32 0.001758 0.02406 -0.045716 0.07384 1.156 1.87e-03 0.0545
## 33 0.053456 -0.06099 -0.012130 -0.08147 1.182 2.28e-03 0.0742
## 34 -0.017328 0.03633 -0.023574 0.06890 1.146 1.63e-03 0.0466
Berdasarkan ringkasan di atas dengan pendekatan Jarak Cook, tidak ada amatan berpengaruh.
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:olsrr':
##
## cement
model.huber<- rlm(y~xpadi+xmiskin, data = datareg, psi=psi.huber,k=10)
summary(model.huber)
##
## Call: rlm(formula = y ~ xpadi + xmiskin, data = datareg, psi = psi.huber,
## k = 10)
## Residuals:
## Min 1Q Median 3Q Max
## -9.0285 -3.2965 -0.4622 2.2292 15.5129
##
## Coefficients:
## Value Std. Error t value
## (Intercept) 25.7220 5.0318 5.1119
## xpadi -5.0804 1.0383 -4.8928
## xmiskin 0.9139 0.1684 5.4260
##
## Residual standard error: 4.024 on 31 degrees of freedom
model.huber.summary <-summary(model.huber)
library(MASS)
model.tukbis <- rlm(y~xpadi+xmiskin, data = datareg, psi = psi.bisquare, )
summary(model.tukbis)
##
## Call: rlm(formula = y ~ xpadi + xmiskin, data = datareg, psi = psi.bisquare)
## Residuals:
## Min 1Q Median 3Q Max
## -9.5044 -2.4878 -0.1062 2.3528 17.6150
##
## Coefficients:
## Value Std. Error t value
## (Intercept) 19.2147 4.3196 4.4483
## xpadi -3.8904 0.8914 -4.3644
## xmiskin 0.9584 0.1446 6.6283
##
## Residual standard error: 3.674 on 31 degrees of freedom
model.tukbis.summary <-summary(model.tukbis)
rss <- sum(resid(model.huber )^2)
tss <- sum((datareg$y - mean(datareg$y))^2)
r_squared <- 1 - (rss / tss)
n <- nrow(datareg)
p <- length(model.huber$coefficients)
adjusted_r_squared <- 1 - ((1 - r_squared) * (n - 1) / (n - p))
print(adjusted_r_squared)
## [1] 0.59011
Model OLS
AIC.OLS<-AIC(model.ols)
AIC.OLS
## [1] 221.4852
Model Stepwise
AIC.best<-AIC(bestmodel)
AIC.best
## [1] 210.7834
Model Regresi Robust: Pembobot Huber
AIC.Huber<-AIC(model.huber)
AIC.Huber
## [1] 210.7834
Model Regresi Robust: Pembobot Tukey Bisquare
AIC.Tukbis<-AIC(model.tukbis)
AIC.Tukbis
## [1] 212.7782
predictions <- predict(model.huber)
residuals <- datareg$y - predictions
tss <- sum((datareg$y - mean(datareg$y))^2)
rss <- sum(residuals^2)
r_squared.Huber <- 1 - (rss / tss)
print(r_squared.Huber)
## [1] 0.6149518
predictions <- predict(model.tukbis)
residuals <- datareg$y - predictions
tss <- sum((datareg$y - mean(datareg$y))^2)
rss <- sum(residuals^2)
r_squared.Tukbis <- 1 - (rss / tss)
print(r_squared.Tukbis)
## [1] 0.5916847
rss <- sum(resid(model.huber )^2)
tss <- sum((datareg$y - mean(datareg$y))^2)
r_squared <- 1 - (rss / tss)
n <- nrow(datareg)
p <- length(model.huber$coefficients)
adjusted_r_squared <- 1 - ((1 - r_squared) * (n - 1) / (n - p))
print(adjusted_r_squared)
## [1] 0.59011
residual_standar <- model.huber.summary$sigma
residual_standar
## [1] 4.023752
residual_standar <- model.tukbis.summary$sigma
residual_standar
## [1] 3.673946
dt<-data.frame(
model=c( "Huber", "Tukbis"),
AIC= c(AIC.Huber, AIC.Tukbis),
R2= c(r_squared.Huber, r_squared.Tukbis)
)
print(dt)
## model AIC R2
## 1 Huber 210.7834 0.6149518
## 2 Tukbis 212.7782 0.5916847
library(MASS)
library(sandwich)
## Warning: package 'sandwich' was built under R version 4.2.3
library(lmtest)
library(car)
hypothesis_test <- linearHypothesis(model.huber,c("xpadi = 0", "xmiskin = 0"))
hypothesis_test
## Linear hypothesis test
##
## Hypothesis:
## xpadi = 0
## xmiskin = 0
##
## Model 1: restricted model
## Model 2: y ~ xpadi + xmiskin
##
## Res.Df Df F Pr(>F)
## 1 33
## 2 31 2 24.755 3.763e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1