Weighted Least Square

Data

Contoh. Seorang peneliti kesehatan, tertarik untuk mempelajari hubungan antara tekanan darah diastolik dan usia di antara wanita dewasa yang sehat berusia 20 hingga 60 tahun, mengumpulkan data dari 54 subjek

library(readr)
blood<-read.table("D:\\bloodpressure.txt",head=T)
str(blood) #dbp = diastolic blood pressure

## 'data.frame':    54 obs. of  3 variables:
##  $ subject: int  1 2 3 4 5 6 7 8 9 10 ...
##  $ age    : int  27 21 22 24 25 23 20 20 29 24 ...
##  $ dbp    : int  73 66 63 75 71 70 65 70 79 72 ...

Model Regresi

mod.reg <- lm(dbp~age,data=blood)
summary(mod.reg)

## 
## Call:
## lm(formula = dbp ~ age, data = blood)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -16.4786  -5.7877  -0.0784   5.6117  19.7813 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 56.15693    3.99367  14.061  < 2e-16 ***
## age          0.58003    0.09695   5.983 2.05e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.146 on 52 degrees of freedom
## Multiple R-squared:  0.4077, Adjusted R-squared:  0.3963 
## F-statistic: 35.79 on 1 and 52 DF,  p-value: 2.05e-07

Scatter plot

plot(blood$age,blood$dbp,xlab="Age",ylab="Diastolic blood pressure")
abline(mod.reg)

Plot residuals vs fitted value

plot(fitted(mod.reg), resid(mod.reg), col = "grey", pch = 20,
xlab = "Fitted", ylab = "Residuals", main = "Blood Data")
abline(h = 0, col = "darkorange", lwd = 2)

Plot standardized residuals vs `Age`

std.res <- rstandard(mod.reg)
plot(blood$age,std.res,xlab="Age",
ylab="Standardized Residuals",ylim=c(-4,4))
abline(h=0) # Terlihat bahwa variasi galat/error konstan

Berdasarkan hasil di atas terlihat bahwa terdapat pola khusus (tidak acak) pada kedua plot sehingga menyebabkan asumsi kehomogenan ragam tidak terpenuhi.

Uji Homogenitas

library(lmtest)
bptest(mod.reg) #studentized Breusch-Pagan test

## 
##  studentized Breusch-Pagan test
## 
## data:  mod.reg
## BP = 12.541, df = 1, p-value = 0.0003981

Berdasarkan hasil pengujian menggunakan uji Breusch Pagan diperoleh p-value < 0.05 sehingga dapat disimpulkan bahwa terjadi pelanggaran asumsi kehomogenan ragam pada model regresi. Salah satu metode yang dapat digunakan untuk menangani permasalahan ragam yang tidak homogen adalah Metode Weighted Least Square (WLS).

Model Regresi dengan WLS

Koefisien regresi dugaan dengan metode WLS dapat disusun dengan terlebih dahulu menyusun bobot dari masing-masing amatan sebagai berikut,

res <- residuals(mod.reg) #residual
abs_res <- abs(res) #residual absolut
sqr_res <- res^2 #square residual
si_hat1 <- fitted(lm(abs_res~blood$age))
wi_hat1 <- 1/si_hat1^2
si_hat2 <- fitted(lm(sqr_res~blood$age))
wi_hat2 <- 1/si_hat2^2
data <-cbind(blood,res,abs_res,si_hat1,si_hat2,wi_hat1,wi_hat2)
head(data)

##   subject age dbp        res   abs_res  si_hat1    si_hat2    wi_hat1
## 1       1  27  73  1.1822391 1.1822391 3.801175 19.9862708 0.06920928
## 2       2  21  66 -2.3375761 2.3375761 2.612141 -0.9660684 0.14655708
## 3       3  22  63 -5.9176069 5.9176069 2.810313  2.5259882 0.12661657
## 4       4  24  75  4.9223315 4.9223315 3.206658  9.5101012 0.09725115
## 5       5  25  71  0.3423007 0.3423007 3.404830 13.0021578 0.08625993
## 6       6  23  70  0.5023623 0.5023623 3.008485  6.0180447 0.11048521
##       wi_hat2
## 1 0.002503436
## 2 1.071480527
## 3 0.156724674
## 4 0.011056807
## 5 0.005915196
## 6 0.027611448

Koefisien Regresi WLS dengan Matriks

X <- as.matrix(cbind(1,blood$age))
Y <- as.matrix(blood$dbp)
W <- diag(wi_hat1)
b_w <- solve(t(X)%*%W%*%X)%*%t(X)%*%W%*%Y #keofisien WLS
b_w

##            [,1]
## [1,] 55.5657664
## [2,]  0.5963417

n <- nrow(blood) #n=54
p <- ncol(X) #p =2
MSEw <- sum(wi_hat1*res^2)/(n-p)
s2_bw <- MSEw*solve(t(X)%*%W%*%X)
s_bw <- sqrt(diag(s2_bw)) # sb0 dan sb1 dari WLS
s_bw

## [1] 2.52231503 0.07928196

Koefisien Regresi WLS dengan fungsi `lm`

Dengan menggunakan fungsi lm untuk model linier, persamaan model regresi dugaan dengan metode WLS juga dapat diperoleh dengan memasukkan bobot yang digunakan.

model.WLS <- lm(dbp~age,data=blood,weights=wi_hat1)
summary(model.WLS)

## 
## Call:
## lm(formula = dbp ~ age, data = blood, weights = wi_hat1)
## 
## Weighted Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.0230 -0.9939 -0.0327  0.9250  2.2008 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 55.56577    2.52092  22.042  < 2e-16 ***
## age          0.59634    0.07924   7.526 7.19e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.213 on 52 degrees of freedom
## Multiple R-squared:  0.5214, Adjusted R-squared:  0.5122 
## F-statistic: 56.64 on 1 and 52 DF,  p-value: 7.187e-10

Plot residuals vs fitted value

Pemeriksaan kembali apakah variasi galat/error sudah konstan

plot(fitted(model.WLS), resid(model.WLS), col = "grey", pch = 20,
xlab = "Fitted", ylab = "Residuals", main = "Blood Data")
abline(h = 0, col = "darkorange", lwd = 2)

Plot standardized residuals vs `Age`

std.res <- rstandard(model.WLS)
plot(blood$age,std.res,xlab="Age",
ylab="Standardized Residuals",ylim=c(-4,4))
abline(h=0) # Terlihat bahwa variasi galat/error konstan

Berdasarkan kedua plot di atas terlihat bahwa keduanya memiliki pola yang acak sehingga asumsi kehomogenan ragam terpenuhi setelah di atasi dengan WLS.

Uji Homogenitas

library(lmtest)
bptest(model.WLS) #studentized Breusch-Pagan test

## 
##  studentized Breusch-Pagan test
## 
## data:  model.WLS
## BP = 0.43608, df = 1, p-value = 0.509

Berdasarkan hasil pengujian menggunakan uji Breusch Pagan diperoleh p-value > 0.05 sehingga dapat disimpulkan bahwa tidak terjadi pelanggaran asumsi kehomogenan ragam pada model regresi WLS.

Weighted Least Square