PML Prak 5

Soal 1

A

y <- c(7, 5, 8, 9)
x <- c(4, 3, 6, 7)

X1 <- cbind(1, x)
X1 

##        x
## [1,] 1 4
## [2,] 1 3
## [3,] 1 6
## [4,] 1 7

XtX1 <- t(X1) %*% X1
XtX1

##        x
##    4  20
## x 20 110

XtX_inv1 <- solve(XtX1) #invers
XtX_inv1

##            x
##    2.75 -0.5
## x -0.50  0.1

XtY1 <- t(X1) %*% y
XtY1

##   [,1]
##     29
## x  154

B

betaduga1 <- XtX_inv1 %*% XtY1
betaduga1

##   [,1]
##   2.75
## x 0.90

C

model1 <- lm(y ~ x) # modreg
model1

## 
## Call:
## lm(formula = y ~ x)
## 
## Coefficients:
## (Intercept)            x  
##        2.75         0.90

sigma1sq <- summary(model1)$sigma^2
sigma1sq

## [1] 0.325

D

yduga1 <- predict(model1) #prediksi y utk tiap x
yduga1

##    1    2    3    4 
## 6.35 5.45 8.15 9.05

sigma1sq <- summary(model1)$sigma^2
sigma1sq

## [1] 0.325

Soal 2

A

y <- c(3, 4, 5, 6)
x <- c(2, 3, 4, 5)

X2 <- cbind(1, x)
X2

##        x
## [1,] 1 2
## [2,] 1 3
## [3,] 1 4
## [4,] 1 5

XtX2 <- t(X2) %*% X2
XtX2

##       x
##    4 14
## x 14 54

XtX_inv2 <- solve(XtX2) #invers
XtX_inv2

##           x
##    2.7 -0.7
## x -0.7  0.2

XtY2 <- t(X2) %*% y
XtY2

##   [,1]
##     18
## x   68

B

betaduga2 <- XtX_inv2 %*% XtY2
betaduga2

##   [,1]
##      1
## x    1

C

model2 <- lm(y ~ x) # modreg
model2

## 
## Call:
## lm(formula = y ~ x)
## 
## Coefficients:
## (Intercept)            x  
##           1            1

sigma2sq <- summary(model2)$sigma^2

## Warning in summary.lm(model2): essentially perfect fit: summary may be
## unreliable

sigma2sq

## [1] 0

D

yduga2 <- predict(model2) #prediksi y utk tiap x
yduga2

## 1 2 3 4 
## 3 4 5 6

sigma2sq <- summary(model2)$sigma^2

## Warning in summary.lm(model2): essentially perfect fit: summary may be
## unreliable

sigma2sq

## [1] 0

Soal 3

A

y <- c(10, 12, 14, 16)
x <- c(1, 2, 3, 4)

X3 <- cbind(1, x)
X3

##        x
## [1,] 1 1
## [2,] 1 2
## [3,] 1 3
## [4,] 1 4

XtX3 <- t(X3) %*% X3
XtX3

##       x
##    4 10
## x 10 30

XtX_inv3 <- solve(XtX3) #invers
XtX_inv3

##           x
##    1.5 -0.5
## x -0.5  0.2

XtY3 <- t(X3) %*% y
XtY3

##   [,1]
##     52
## x  140

B

betaduga3 <- XtX_inv3 %*% XtY3
betaduga3

##   [,1]
##      8
## x    2

C

model3 <- lm(y ~ x) # modreg
model3

## 
## Call:
## lm(formula = y ~ x)
## 
## Coefficients:
## (Intercept)            x  
##           8            2

sigma3sq <- summary(model3)$sigma^2

## Warning in summary.lm(model3): essentially perfect fit: summary may be
## unreliable

sigma3sq

## [1] 0

D

yduga3 <- predict(model3) #prediksi y utk tiap x
yduga3

##  1  2  3  4 
## 10 12 14 16

sigma3sq <- summary(model3)$sigma^2

## Warning in summary.lm(model3): essentially perfect fit: summary may be
## unreliable

sigma3sq

## [1] 0

Soal 4

A

y <- c(2, 3, 5, 7)
x <- c(1, 2, 4, 6)

X4 <- cbind(1, x)
X4

##        x
## [1,] 1 1
## [2,] 1 2
## [3,] 1 4
## [4,] 1 6

XtX4 <- t(X4) %*% X4
XtX4

##       x
##    4 13
## x 13 57

XtX_inv4 <- solve(XtX4) #invers
XtX_inv4

##                        x
##    0.9661017 -0.22033898
## x -0.2203390  0.06779661

XtY4 <- t(X4) %*% y
XtY4

##   [,1]
##     17
## x   70

B

betaduga4 <- XtX_inv4 %*% XtY4
betaduga4

##   [,1]
##      1
## x    1

C

model4 <- lm(y ~ x) # modreg
model4

## 
## Call:
## lm(formula = y ~ x)
## 
## Coefficients:
## (Intercept)            x  
##           1            1

sigma4sq <- summary(model4)$sigma^2

## Warning in summary.lm(model4): essentially perfect fit: summary may be
## unreliable

sigma4sq

## [1] 4.930381e-31

D

yduga4 <- predict(model4) #prediksi y utk tiap x
yduga4

## 1 2 3 4 
## 2 3 5 7

sigma4sq <- summary(model4)$sigma^2

## Warning in summary.lm(model4): essentially perfect fit: summary may be
## unreliable

sigma4sq

## [1] 4.930381e-31

Soal 5

A

y <- c(15, 10, 20, 25)
x <- c(5, 3, 6, 8)

X5 <- cbind(1, x)
X5

##        x
## [1,] 1 5
## [2,] 1 3
## [3,] 1 6
## [4,] 1 8

XtX5 <- t(X5) %*% X5
XtX5

##        x
##    4  22
## x 22 134

XtX_inv5 <- solve(XtX5) #invers
XtX_inv5

##                        x
##    2.5769231 -0.42307692
## x -0.4230769  0.07692308

XtY5 <- t(X5) %*% y
XtY5

##   [,1]
##     70
## x  425

B

betaduga5 <- XtX_inv5 %*% XtY5
betaduga5

##        [,1]
##   0.5769231
## x 3.0769231

C

model5 <- lm(y ~ x) # modreg
model5

## 
## Call:
## lm(formula = y ~ x)
## 
## Coefficients:
## (Intercept)            x  
##      0.5769       3.0769

sigma5sq <- summary(model5)$sigma^2
sigma5sq

## [1] 0.9615385

D

yduga2 <- predict(model2) #prediksi y utk tiap x
yduga2

## 1 2 3 4 
## 3 4 5 6

sigma2sq <- summary(model2)$sigma^2

## Warning in summary.lm(model2): essentially perfect fit: summary may be
## unreliable

sigma2sq

## [1] 0

Soal 6

A

x <- c(2000, 2500, 1800, 2200, 2400)
y <- c(60, 65, 58, 62, 64)            

X6 <- cbind(1, x)
X6

##           x
## [1,] 1 2000
## [2,] 1 2500
## [3,] 1 1800
## [4,] 1 2200
## [5,] 1 2400

B

XtX6 <- t(X6) %*% X6
XtX6

##                x
##       5    10900
## x 10900 24090000

XtX_inv6 <- solve(XtX6) #invers
XtX_inv6

##                            x
##   14.689024390 -6.646341e-03
## x -0.006646341  3.048780e-06

XtY6 <- t(X6) %*% y
XtY6

##     [,1]
##      309
## x 676900

C

beta6 <- XtX_inv6 %*% XtY6 #least sq
beta6

##    [,1]
##   40.00
## x  0.01

D

yduga6 <- X6 %*% beta6 # prediksi yduga

n <- length(y)
s2 <- sum((y - yduga6)^2) / (n - 2)  # ragam galat
s2

## [1] 3.674114e-25

E

xisq6 <- x*x
xisq6

## [1] 4000000 6250000 3240000 4840000 5760000

varbeta6 <- s2/sum(xisq6)
varbeta6

## [1] 1.525162e-32

Soal 7

A

x <- c(1, 2, 3, 4, 5)
y <- c(55, 60, 70, 75, 80)            

X7 <- cbind(1, x)
X7

##        x
## [1,] 1 1
## [2,] 1 2
## [3,] 1 3
## [4,] 1 4
## [5,] 1 5

B

XtX7 <- t(X7) %*% X7
XtX7

##       x
##    5 15
## x 15 55

XtX_inv7 <- solve(XtX7) #invers
XtX_inv7

##           x
##    1.1 -0.3
## x -0.3  0.1

XtY7 <- t(X7) %*% y
XtY7

##   [,1]
##    340
## x 1085

C

beta7 <- XtX_inv7 %*% XtY7 #least sq
beta7

##   [,1]
##   48.5
## x  6.5

D

yduga7 <- X7 %*% beta7 # prediksi yduga

n <- length(y)
s2 <- sum((y - yduga7)^2) / (n - 2)  # ragam galat
s2

## [1] 2.5

E

xisq7 <- x*x
xisq7

## [1]  1  4  9 16 25

varbeta7 <- s2/sum(xisq7)
varbeta7

## [1] 0.04545455

Soal 8

A

x <- c(20, 30, 40, 50, 60)
y <- c(1, 2, 3, 4, 5)            

X8 <- cbind(1, x)
X8

##         x
## [1,] 1 20
## [2,] 1 30
## [3,] 1 40
## [4,] 1 50
## [5,] 1 60

B

XtX8 <- t(X8) %*% X8
XtX8

##          x
##     5  200
## x 200 9000

XtX_inv8 <- solve(XtX8) #invers
XtX_inv8

##              x
##    1.80 -0.040
## x -0.04  0.001

XtY8 <- t(X8) %*% y
XtY8

##   [,1]
##     15
## x  700

C

beta8 <- XtX_inv8 %*% XtY8 #least sq
beta8

##   [,1]
##   -1.0
## x  0.1

D

yduga8 <- X8 %*% beta8 # prediksi yduga

n <- length(y)
s2 <- sum((y - yduga8)^2) / (n - 2)  # ragam galat
s2

## [1] 2.366583e-29

E

xisq8 <- x*x
xisq8

## [1]  400  900 1600 2500 3600

varbeta8 <- s2/sum(xisq8)
varbeta8

## [1] 2.629536e-33

Soal 9

A

x <- c(100, 200, 300, 400, 500)
y <- c(0.5, 1.0, 1.5, 2.0, 2.5)            

X9 <- cbind(1, x)
X9

##          x
## [1,] 1 100
## [2,] 1 200
## [3,] 1 300
## [4,] 1 400
## [5,] 1 500

B

XtX9 <- t(X9) %*% X9
XtX9

##             x
##      5   1500
## x 1500 550000

XtX_inv9 <- solve(XtX9) #invers
XtX_inv9

##               x
##    1.100 -3e-03
## x -0.003  1e-05

XtY9 <- t(X9) %*% y
XtY9

##     [,1]
##      7.5
## x 2750.0

C

beta9 <- XtX_inv9 %*% XtY9 #least sq
beta9

##            [,1]
##   -1.776357e-15
## x  5.000000e-03

D

yduga9 <- X9 %*% beta9 # prediksi yduga

n <- length(y)
s2 <- sum((y - yduga9)^2) / (n - 2)  # ragam galat
s2

## [1] 3.718329e-30

E

xisq9 <- x*x
xisq9

## [1]  10000  40000  90000 160000 250000

varbeta9 <- s2/sum(xisq9)
varbeta9

## [1] 6.760598e-36

Soal 10

A

x <- c(1, 2, 3, 4, 5)
y <- c(3, 5, 7, 9, 11)            

X10 <- cbind(1, x)
X10

##        x
## [1,] 1 1
## [2,] 1 2
## [3,] 1 3
## [4,] 1 4
## [5,] 1 5

B

XtX10 <- t(X10) %*% X10
XtX10

##       x
##    5 15
## x 15 55

XtX_inv10 <- solve(XtX10) #invers
XtX_inv10

##           x
##    1.1 -0.3
## x -0.3  0.1

XtY10 <- t(X10) %*% y
XtY10

##   [,1]
##     35
## x  125

C

beta10 <- XtX_inv10 %*% XtY10 #least sq
beta10

##   [,1]
##      1
## x    2

D

yduga10 <- X10 %*% beta10 # prediksi yduga

n <- length(y)
s2 <- sum((y - yduga10)^2) / (n - 2)  # ragam galat
s2

## [1] 8.414516e-29

E

xisq10 <- x*x
xisq10

## [1]  1  4  9 16 25

varbeta10 <- s2/sum(xisq10)
varbeta10

## [1] 1.529912e-30