Heteroskedasticity: residual กระจายไม่เท่ากัน สิ่งที่เสี่ยงคือ Standard Error ไม่ใช่ Coefficient
On this page
Read the English version
บทคัดย่อ
Skedasticity คือความกว้างของ residual ใน regression คือดูว่าระยะห่างระหว่างค่าจริงกับค่าที่โมเดลทำนาย กระจายกว้างเท่ากันไหมในทุกระดับของค่าที่ทำนาย ถ้าความกว้างใกล้เคียงกันตลอด เช่น ทำนาย systolic blood pressure จากอายุแล้วคลาดไปประมาณ 10 mmHg ในทุกช่วงอายุ เรียกว่า error เป็น homoskedastic ถ้าความกว้างเพิ่มขึ้น เช่น ค่ารักษาประมาณ 500,000 บาท ทำนายคลาดได้ประมาณ 200,000 บาท แต่ค่ารักษาประมาณ 5,000 บาท คลาดแค่ 1,000 บาท residual plot จะบานออกเป็นรูปพัด และ error เป็น heteroskedastic Ordinary least squares มี assumption ว่า variance ของ error คงที่ การผิด assumption นี้ไม่ได้ทำให้ coefficient ผิดทันที แต่ทำให้ standard error แบบปกติผิด และพาให้ 95% confidence interval, P value และ hypothesis test ผิดตามไปด้วย เช่น slope เท่ากับ 2 อาจได้ standard error 0.5 ทั้งที่ค่าที่ถูกคือ 1.0 บทความนี้แสดงวิธีจำรูปแบบนี้ให้ได้ และเหตุผลที่คำตอบคือ heteroskedasticity-robust standard error
โมเดลค่ารักษาที่ residual บานออกเป็นรูปพัด
ทีมโรงพยาบาลทีมหนึ่งสร้าง linear regression ทำนายค่ารักษาต่อการนอนโรงพยาบาลหนึ่งครั้ง จากอายุของผู้ป่วย จำนวนวันนอน และจำนวนหัตถการ coefficient ทุกตัวดูสมเหตุสมผลทางคลินิก จากนั้นมีคน plot residual เทียบกับ fitted value สำหรับ admission ที่ค่ารักษาถูก จุดเกาะอยู่ใกล้ศูนย์ แต่สำหรับ admission ที่แพงที่สุด จุดกระจายออกไปไกลทั้งด้านบนและด้านล่าง
คำถามสองเรื่องตามมาทันที coefficient ของจำนวนวันนอนผิดไปแล้วหรือเปล่า และทีมยังเชื่อ confidence interval กับ P value ของมันได้อยู่ไหม คำตอบของทั้งสองเรื่องเริ่มจากคำคำเดียว: skedasticity
Skedasticity: residual กระจายกว้างแค่ไหน
คำว่า skedasticity หมายถึงเรื่องของความกระจายของ error หรือ residual ใน regression พูดให้ง่ายที่สุด regression model แบ่งค่าผลลัพธ์ทุกค่าออกเป็นสองส่วน:
$$Y = \text{the model's predicted value} + \text{error}$$ค่าที่โมเดลทำนายเรียกอีกชื่อว่า fitted value ส่วน residual คือระยะห่างระหว่างค่าจริงที่สังเกตได้กับค่าที่โมเดลทำนาย คือค่าจริงลบค่าที่ทำนาย
สิ่งที่ skedasticity สนใจคือคำถามเดียวเกี่ยวกับ residual เหล่านั้น:
residual กระจายกว้างเท่ากันไหม ในทุกระดับของค่าที่ทำนาย
Homoskedasticity: กระจายพอ ๆ กันทุกช่วง
สมมติโมเดลทำนาย systolic blood pressure (SBP) จากอายุ ที่อายุ 20 ปี ค่าที่ทำนายคลาดไปประมาณ ±10 mmHg ที่อายุ 40 ปี คลาดประมาณ ±10 mmHg และที่อายุ 70 ปี ก็คลาดประมาณ ±10 mmHg อีก variance ของ residual ค่อนข้างคงที่
แบบนี้เรียกว่า homoskedasticity ภาพในหัวคือ plot ของ residual เทียบกับ fitted value ที่เป็นแถบกว้างใกล้เคียงกันตลอดแนว เหมือนแผงซ้ายของภาพด้านล่าง
| อายุ | ความคลาดเคลื่อนโดยทั่วไปของการทำนาย SBP |
|---|---|
| 20 ปี | ประมาณ ±10 mmHg |
| 40 ปี | ประมาณ ±10 mmHg |
| 70 ปี | ประมาณ ±10 mmHg |
Heteroskedasticity: ความกระจายเปลี่ยนไป
คราวนี้สมมติโมเดลทำนาย hospital cost คนที่ค่ารักษาประมาณ 5,000 บาท ค่าที่ทำนายอาจคลาดแค่ ประมาณ ±1,000 บาท แต่คนที่ค่ารักษาประมาณ 500,000 บาท อาจคลาดถึงประมาณ ±200,000 บาท
พอ fitted value สูงขึ้น residual ก็ยิ่งกระจายกว้าง บน residual plot จุดจะบานออกเป็นรูปกรวย หรือ fan shape เหมือนแผงขวาด้านล่าง แบบนี้เรียกว่า heteroskedasticity
| ค่ารักษาที่ทำนาย | ความคลาดเคลื่อนโดยทั่วไปของการทำนาย |
|---|---|
| ประมาณ 5,000 บาท | ประมาณ ±1,000 บาท |
| ประมาณ 500,000 บาท | ประมาณ ±200,000 บาท |
ใช้แนวคิดนี้ที่ไหน: linear regression และ OLS
แนวคิดนี้ใช้หลัก ๆ ใน linear regression ที่ fit ด้วย ordinary least squares (OLS) ซึ่งเป็นวิธี fit ที่ใช้กันทั่วไป คือเลือก coefficient ที่ทำให้ผลรวมกำลังสองของ residual น้อยที่สุด สมมติว่าเราวิเคราะห์
$$\text{Blood pressure} = \beta_0 + \beta_1(\text{age}) + \beta_2(\text{BMI}) + \varepsilon$$โดย BMI คือ body mass index และ $\varepsilon$ คือ error OLS มี assumption สำคัญเกี่ยวกับ error ตัวนี้ว่า
$$\mathrm{Var}(\varepsilon \mid X) = \sigma^2$$ในที่นี้ $X$ แทน covariate คืออายุและ BMI ส่วน $\sigma^2$ คือ variance ค่าเดียวที่คงที่ แปลว่า เมื่อ $X$ เปลี่ยนไป variance ของ error ควรใกล้เคียงกัน
ถ้า variance ไม่คงที่ เราเรียกว่า error มี heteroskedasticity
ทำไมต้องสนใจ: coefficient อาจยังดี แต่ standard error อาจไม่ดี
จุดสำคัญที่สุดคือ heteroskedasticity ไม่ได้แปลว่า coefficient ผิดทันที ในโมเดล blood pressure ข้างบน estimate อย่าง
$$\hat\beta_{\text{BMI}} = 2.5$$อาจยังเป็น estimate ที่ดีอยู่
แต่ปัญหาคือ standard error อาจคำนวณผิด สิ่งที่สร้างขึ้นบน standard error ก็อาจผิดตามไปด้วย:
- 95% confidence interval อาจผิด
- P value อาจผิด
- hypothesis test อาจผิด
assumption เรื่อง variance ของ error
ถ้าพูดให้เป็นทางการ skedasticity ใช้กับ assumption เรื่อง variance ของ residual หรือ error ใน regression โดยเฉพาะ linear regression ที่ fit ด้วย OLS assumption นั้นคือ
$$\mathrm{Var}(\varepsilon \mid X) = \sigma^2$$แปลว่า เมื่อค่า $X$ เปลี่ยนไป ความแปรปรวนของ error ควรคงที่ นี่คือ homoskedasticity
ถ้า $\mathrm{Var}(\varepsilon \mid X)$ เปลี่ยนตาม $X$ หรือตาม fitted value เราเรียกว่า heteroskedasticity
ทำไม assumption นี้จำเป็น: OLS ให้สองเรื่องที่แยกกัน
OLS ให้สองเรื่องที่แยกกัน:
- coefficient estimate เช่น $\hat\beta_1$
- uncertainty ของ coefficient คือ standard error, 95% confidence interval และ P value
ถ้า mean model ถูกต้องและมี exogeneity (error มีค่าเฉลี่ยเป็นศูนย์ไม่ว่า $X$ จะเป็นค่าใด covariate จึงไม่มีข้อมูลอะไรเกี่ยวกับ error) coefficient ก็ยัง estimate ได้ดีแม้มี heteroskedasticity
แต่ standard error แบบปกติของ OLS จะไม่ถูกต้อง เพราะสูตรนั้นอาศัย assumption ว่า variance ของ error คงที่ ผลกระทบหลักคือ:
- standard error ผิด
- 95% confidence interval ผิด
- P value ผิด
- hypothesis test อาจสรุปผิด
ตัวอย่างคำนวณด้วยมือ: slope เท่ากับ 2 ที่ดูแม่นกว่าความจริง
ตัวอย่างคำนวณด้วยมือ ใช้โมเดลที่ง่ายที่สุด มี covariate ตัวเดียวคือ $X$:
$$Y = \beta_0 + \beta_1 X + \varepsilon$$-
ค่า estimate
\[ \hat\beta_1 = 2 \]
สมมติ slope ที่ fit ได้คือ 2 heteroskedasticity โดยตัวมันเองไม่ได้เปลี่ยนตัวเลขนี้
-
ถ้า error เป็น homoskedastic
\[ SE = 0.5 \]
เมื่อ variance ของ error คงที่ standard error แบบปกติของ OLS ถูกต้อง คือ 0.5
-
ถ้า error เป็น heteroskedastic แต่ยังใช้ standard error แบบปกติ
\[ SE_{\text{ordinary}} = 0.5 \quad \text{but it should be} \quad SE = 1.0 \]
เมื่อความกระจายเป็นรูปพัด สูตรปกติอาจยังให้ 0.5 ทั้งที่ standard error จริง ๆ ควรเป็น 1.0
-
95% confidence interval จะอ้างว่าอะไร
\[ 2 \pm 1.96 \times 0.5 = 1.02 \text{ to } 2.98 \quad \text{versus} \quad 2 \pm 1.96 \times 1.0 = 0.04 \text{ to } 3.96 \]
ถ้าสร้างจาก standard error แบบปกติ interval จะกว้างแค่ครึ่งเดียวของที่ควรเป็น และ P value จะเล็กกว่าที่ควรเป็น
ผลลัพธ์: standard error แบบปกติบอก 0.5 ทั้งที่ควรเป็น 1.0 เราจึงอาจคิดว่า effect “แม่น” เกินจริง
รูปแบบเดียวกันใน cohort จำลอง
เรื่องเดียวกันเห็นได้ใน cohort จำลองของผู้ใหญ่ 2,000 คน อายุ 30 ถึง 80 ปี ที่ systolic blood pressure กระจายมากขึ้นตามอายุ: ความกระจายแบบสุ่มที่เติมให้ SBP ของแต่ละคนมี standard deviation 4 mmHg ที่อายุ 30 ปี 9 ที่อายุ 55 ปี และ 24 ที่อายุ 80 ปี residual จากเส้นตรงของ SBP ตามอายุจึงบานออกไปทางอายุมาก
regression เส้นตรงของ SBP ตามอายุให้ slope 0.66 mmHg ต่อปี ใกล้กับ slope เส้นตรง 0.67 ในประชากรจำลองขนาดใหญ่ที่ cohort นี้สุ่มมา standard error แบบปกติคือ 0.0215 เมื่อ fit ใหม่ด้วย heteroskedasticity-robust standard error คือ vce(robust) ของ Stata slope เท่าเดิมทุกประการ และ standard error เพิ่มเป็น 0.0243 คือ 1.13 เท่าของค่าปกติ
เหมือนในตัวอย่างคำนวณด้วยมือ coefficient อยู่ที่เดิม มีแค่ uncertainty ที่เปลี่ยน ในข้อมูลนี้ standard error แบบปกติเล็กเกินไป แต่มันไม่ได้เล็กเกินไปในทุกชุดข้อมูล ซึ่งเป็นอีกเหตุผลหนึ่งที่ไม่ควรพึ่งมันเมื่อความกระจายไม่เท่ากัน
Stata: SBP ตามอายุ พร้อม standard error แบบปกติ
* Specification A: SBP on age (the high-leverage ages are the high-variance ages)
regress sbp age. * Specification A: SBP on age (the high-leverage ages are the high-variance ages)
. regress sbp age
Source | SS df MS Number of obs = 2,000
-------------+---------------------------------- F(1, 1998) = 951.44
Model | 189678.648 1 189678.648 Prob > F = 0.0000
Residual | 398321.189 1,998 199.359955 R-squared = 0.3226
-------------+---------------------------------- Adj R-squared = 0.3222
Total | 587999.837 1,999 294.146992 Root MSE = 14.119
------------------------------------------------------------------------------
sbp | Coefficient Std. err. t P>|t| [95% conf. interval]
-------------+----------------------------------------------------------------
age | .6644293 .0215406 30.85 0.000 .6221848 .7066737
_cons | 88.26957 1.232598 71.61 0.000 85.85226 90.68689
------------------------------------------------------------------------------Stata: fit เดิมพร้อม heteroskedasticity-robust standard error
regress sbp age, vce(robust). regress sbp age, vce(robust)
Linear regression Number of obs = 2,000
F(1, 1998) = 750.63
Prob > F = 0.0000
R-squared = 0.3226
Root MSE = 14.119
------------------------------------------------------------------------------
| Robust
sbp | Coefficient std. err. t P>|t| [95% conf. interval]
-------------+----------------------------------------------------------------
age | .6644293 .0242514 27.40 0.000 .6168686 .7119899
_cons | 88.26957 1.203478 73.35 0.000 85.90937 90.62978
------------------------------------------------------------------------------vce(robust): slope ไม่เปลี่ยน ส่วน standard error ของมันเพิ่มจาก 0.0215 เป็น 0.0243อ่านรูปแบบนี้: ค่าเฉลี่ยโอเค แต่ variance ไม่โอเค
รูปแบบนี้มักถูกสรุปไว้ในประโยคเดียว ในการทบทวนระเบียบวิธีหรือในรายงานทางสถิติ:
The mean structure is appropriate, but the residual variance changes with the fitted values.
ประโยคนี้แปลว่า coefficient model ไม่ได้ผิด แต่ variance assumption ผิด coefficient จึงอยู่ที่เดิม ที่เปลี่ยนมีแค่วิธีคำนวณ uncertainty ของมัน:
Use heteroskedasticity-robust standard errors.
ความเข้าใจผิดที่พบบ่อยและวิธีแก้
-
"residual บานออกเป็นรูปพัด แปลว่า coefficient ผิด"
Heteroskedasticity เป็นเรื่องของ variance ของ error ไม่ใช่ค่าเฉลี่ยของมัน ถ้า mean model ถูกต้อง coefficient อย่าง estimate ของผลของ BMI ที่เท่ากับ 2.5 ก็ยังเป็น estimate ที่ดีได้
วิธีแก้: เก็บ coefficient ไว้ แล้วซ่อม uncertainty ของมัน: รายงาน heteroskedasticity-robust standard error
-
"heteroskedasticity แปลว่าเส้น regression โค้ง"
เส้นโค้งหรือไม่เป็นคำถามของ mean model ส่วน heteroskedasticity ถามอีกเรื่องหนึ่ง คือความกระจายรอบเส้นกว้างเท่ากันตลอดแนวหรือไม่
วิธีแก้: ดูความกว้างของแถบใน plot ของ residual เทียบกับ fitted value ไม่ใช่ดูรูปทรงของมัน
-
"coefficient ดูดี standard error ที่พิมพ์ออกมาก็ต้องดีด้วย"
standard error แบบปกติของ OLS สร้างบน assumption ว่า variance ของ error คงที่ เมื่อความกระจายเป็นรูปพัด มันอาจบอก 0.5 ทั้งที่ค่าที่ถูกคือ 1.0
วิธีแก้: เมื่อความกระจายไม่เท่ากัน ให้อ่าน robust standard error และ confidence interval กับ P value ที่สร้างจากมัน
-
"ความกระจายไม่เท่ากัน แปลว่าต้องสร้างโมเดลใหม่ทั้งหมด"
เมื่อ mean structure เหมาะสม และมีแค่ residual variance ที่เปลี่ยนตาม fitted value สมการค่าเฉลี่ยไม่ใช่ปัญหา ปัญหาคือ variance assumption
วิธีแก้: เก็บ mean model ไว้ และใช้ heteroskedasticity-robust standard errors
อภิธานศัพท์
- skedasticity
- ความกระจาย หรือความกว้าง ของ residual รอบค่าที่ regression ทำนาย และดูว่าความกว้างนั้นเท่ากันในทุกค่าที่ทำนายหรือไม่
- homoskedasticity
- residual กระจายกว้างพอ ๆ กันในทุกระดับของ fitted value คือ variance ของ error คงที่
- heteroskedasticity (ความแปรปรวนไม่คงที่)
- variance ของ residual ที่เปลี่ยนไปตาม covariate หรือ fitted value มักเห็นเป็นรูปพัดใน residual plot
- residual
- ค่าที่สังเกตได้ลบด้วยค่าที่โมเดลทำนาย
- fitted value
- ค่าที่ regression model ทำนายให้แต่ละคนจาก covariate ของคนนั้น
- ordinary least squares (OLS)
- วิธีปกติในการ fit linear regression คือเลือก coefficient ที่ทำให้ผลรวมกำลังสองของ residual น้อยที่สุด
- exogeneity
- error มีค่าเฉลี่ยเป็นศูนย์ไม่ว่า covariate จะเป็นค่าใด covariate จึงไม่มีข้อมูลเกี่ยวกับ error
- model-based standard error (standard error แบบอิงแบบจำลอง)
- standard error แบบปกติจากสูตรของโมเดลเอง ใน OLS มันสมมติว่า variance ของ error มีค่าเดียวคงที่
- robust (sandwich) standard error (ค่าคลาดเคลื่อนมาตรฐานแบบแซนด์วิช)
- standard error ที่สร้างจาก residual ยกกำลังสองของแต่ละคน จึงไม่พึ่ง assumption ว่า variance คงที่ ตัว coefficient ไม่เปลี่ยน
เอกสารอ้างอิง
- White H. A heteroskedasticity-consistent covariance matrix estimator and a direct test for heteroskedasticity. Econometrica. 1980;48(4):817-838. doi:10.2307/1912934 https://doi.org/10.2307/1912934
ประเด็นสำคัญ
- Skedasticity คือความกว้างของ residual: homo แปลว่า variance เท่ากันในทุก fitted value ส่วน hetero แปลว่า variance ต่างกัน
- Heteroskedasticity ไม่ได้ถามว่าเส้น regression โค้งไหม แต่ถามว่าความกระจายรอบเส้นกว้างเท่ากันตลอดหรือไม่ และลายเซ็นของมันคือรูปพัดใน plot ของ residual เทียบกับ fitted value
- OLS มี assumption ว่า variance ของ error คงที่ เมื่อ assumption นี้ผิด coefficient ไม่ได้ผิดโดยอัตโนมัติ
- สิ่งที่ผิดคือ standard error แบบปกติ และ 95% confidence interval, P value และ test ที่ตามมา: standard error 0.5 ทั้งที่ค่าที่ถูกคือ 1.0 ทำให้ effect ดูแม่นกว่าความจริง
- เมื่อ mean structure เหมาะสม แต่ residual variance เปลี่ยนตาม fitted value สิ่งที่ผิดคือ variance assumption ไม่ใช่ mean model: ใช้ heteroskedasticity-robust standard errors
อ่านต่อในวิกิ: [[risk-regression-models-epidemiology-th]] [[repeated-measures-modeling-guide-th]]