วิธีการคำนวณข้อผิดพลาดมาตรฐานที่แข็งแกร่งใน r

ข้อสันนิษฐานประการหนึ่งของการถดถอยเชิงเส้น คือ ส่วนที่เหลือของ แบบจำลองจะกระจายเท่ากันในแต่ละระดับของตัวแปรทำนาย

เมื่อไม่เป็นไปตามสมมติฐานนี้ กล่าวกันว่าความ มีความต่างแบบเฮเทอโรสเกดาติสติก ปรากฏอยู่ในแบบจำลองการถดถอย

เมื่อสิ่งนี้เกิดขึ้น ข้อผิดพลาดมาตรฐานของสัมประสิทธิ์การถดถอยของแบบจำลองจะไม่น่าเชื่อถือ

เพื่ออธิบายสิ่งนี้ เราสามารถคำนวณ ข้อผิดพลาดมาตรฐานที่แข็งแกร่ง ซึ่ง “แข็งแกร่ง” เมื่อเทียบกับค่าความคลาดเคลื่อนแบบเฮเทอโรสเคดาสติก และช่วยให้เราเข้าใจค่าความผิดพลาดมาตรฐานที่แท้จริงสำหรับค่าสัมประสิทธิ์การถดถอยได้ดีขึ้น

ตัวอย่างต่อไปนี้แสดงวิธีการคำนวณข้อผิดพลาดมาตรฐานที่แข็งแกร่งสำหรับแบบจำลองการถดถอยใน R

ตัวอย่าง: การคำนวณข้อผิดพลาดมาตรฐานที่แข็งแกร่งใน R

สมมติว่าเรามีกรอบข้อมูลต่อไปนี้ใน R ซึ่งมีข้อมูลเกี่ยวกับชั่วโมงเรียนและคะแนนสอบที่นักเรียน 20 คนในชั้นเรียนได้รับ:

 #create data frame
df <- data. frame (hours=c(1, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4,
                         4, 5, 5, 5, 6, 6, 7, 7, 8),
                 score=c(67, 68, 74, 70, 71, 75, 80, 70, 84, 72,
                         88, 75, 95, 75, 99, 78, 99, 65, 96, 70))

#view head of data frame
head(df)

  hours score
1 1 67
2 1 68
3 1 74
4 1 70
5 2 71
6 2 75

เราสามารถใช้ฟังก์ชัน lm() เพื่อให้พอดีกับโมเดลการถดถอยใน R ที่ใช้ ชั่วโมง เป็นตัวแปรทำนายและ ใช้คะแนน เป็นตัวแปรตอบสนอง:

 #fit regression model
fit <- lm(score ~ hours, data=df)

#view summary of model
summary(fit)

Call:
lm(formula = score ~ hours, data = df)

Residuals:
    Min 1Q Median 3Q Max 
-19,775 -5,298 -3,521 7,520 18,116 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 71.158 4.708 15.11 1.14e-11 ***
hours 1.945 1.075 1.81 0.087 .  
---
Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 10.48 on 18 degrees of freedom
Multiple R-squared: 0.154, Adjusted R-squared: 0.107 
F-statistic: 3.278 on 1 and 18 DF, p-value: 0.08696

วิธีที่ง่ายที่สุดในการตรวจสอบด้วยสายตาว่าปัญหาแบบเฮเทอโรสเคดาสติติสเป็นปัญหาในแบบจำลองการถดถอยหรือไม่คือการสร้างพล็อตส่วนที่เหลือ:

 #create residual vs. fitted plot
plot(fitted(fit), reside(fit))

#add a horizontal line at y=0 
abline(0,0) 

แกน x แสดงค่าที่พอดีของตัวแปรตอบสนอง และแกน y แสดงค่าคงเหลือที่สอดคล้องกัน

จากกราฟเราจะเห็นว่าความแปรปรวนของค่าคงเหลือเพิ่มขึ้นเมื่อค่าที่ติดตั้งเพิ่มขึ้น

สิ่งนี้บ่งชี้ว่าความต่างกันมีแนวโน้มที่จะเป็นปัญหาในแบบจำลองการถดถอย และข้อผิดพลาดมาตรฐานของสรุปแบบจำลองนั้นไม่น่าเชื่อถือ

ในการคำนวณข้อผิดพลาดมาตรฐานที่เข้มงวด เราสามารถใช้ฟังก์ชัน coeftest() จากแพ็คเกจ lmtest และฟังก์ชัน vcovHC() จากแพ็คเกจ แซนวิช ได้ดังนี้:

 library (lmtest)
library (sandwich)

#calculate robust standard errors for model coefficients
coeftest(fit, vcov = vcovHC(fit, type = ' HC0 '))

t test of coefficients:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 71.1576 3.3072 21.5160 2.719e-14 ***
hours 1.9454 1.2072 1.6115 0.1245    
---
Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

โปรดทราบว่าข้อผิดพลาดมาตรฐานสำหรับตัวแปรตัวทำนาย ชั่วโมง เพิ่มขึ้นจาก 1.075 ในสรุปโมเดลก่อนหน้าเป็น 1.2072 ในสรุปโมเดลนี้

เนื่องจากความต่างศักย์ต่างกันมีอยู่ในแบบจำลองการถดถอยดั้งเดิม การประมาณค่าข้อผิดพลาดมาตรฐานนี้จึงมีความน่าเชื่อถือมากกว่า และควรใช้เมื่อคำนวณช่วงความเชื่อมั่นสำหรับตัวแปรทำนาย ชั่วโมง

หมายเหตุ : ประเภทการประมาณค่าที่พบบ่อยที่สุดในการคำนวณในฟังก์ชัน vcovHC() คือ ‘HC0’ แต่คุณสามารถดูการประมาณค่าประเภทอื่นได้จาก เอกสารประกอบ

แหล่งข้อมูลเพิ่มเติม

บทช่วยสอนต่อไปนี้จะอธิบายวิธีดำเนินการงานทั่วไปอื่นๆ ใน R:

วิธีทำการทดสอบของ White สำหรับความต่างศักย์ไฟฟ้าใน R
วิธีการตีความเอาต์พุตการถดถอยเชิงเส้นใน R
วิธีสร้างพล็อตที่เหลือใน R