สัมประสิทธิ์เท่าเดิม แต่ standard error ต่างกัน: คำนวณ sandwich estimator ด้วยมือ

Clinical Epidemiology ResearchMethodology and Research Design THUniqcret doctor knowledges TH
สัมประสิทธิ์เท่าเดิม แต่ standard error ต่างกัน: คำนวณ sandwich estimator ด้วยมือ
On this page

Read the English version

บทคัดย่อ

การประมาณ regression coefficient กับการประมาณความไม่แน่นอนของ coefficient เป็นคนละกระบวนการ Robust standard error ไม่ใช่ model ใหม่ แต่เป็นการประมาณ variance-covariance matrix อีกแบบหนึ่งภายใต้ Poisson model เดิม ในการศึกษาสมมติเรื่อง treatment กับจำนวนครั้งที่ผู้ป่วยเข้ารับการรักษาในโรงพยาบาล coefficient ยังคงเป็น 0.40 ขณะที่ standard error เพิ่มจาก 0.10 เป็น 0.19 ข้อมูลสองชุด ชุดละผู้ป่วยสี่คน ซึ่งมี mean เท่ากันที่ 5 ครั้ง อธิบายว่าทำไมจึงเป็นเช่นนั้น Poisson regression ให้ทั้งสองชุดมี coefficient log 5 = 1.609 และ conventional standard error 0.224 เท่ากัน เพราะสมมติว่า variance เท่ากับ mean ส่วน robust (HC0) standard error อ่านจาก squared residuals ได้ 0.071 สำหรับชุดที่ค่าอยู่ใกล้กัน และ 0.400 สำหรับชุดที่กระจายมาก สูตร sandwich ซึ่งประกอบด้วย bread และ meat ขยายแนวคิดนี้ไปยัง model ที่มี predictors หลายตัว และ modified Poisson regression นำแนวคิดเดียวกันไปใช้กับ binary outcome บทความนี้คำนวณทีละขั้นด้วยมือ แล้วให้คำสั่ง Stata ขอบเขตของ robust standard error และตัวอย่างการเขียน Methods


ภาพสรุป ตัวอย่างคำนวณด้วยมือ: coefficient เท่ากันที่ 1.609 โดยมี conventional standard error 0.224 และ robust (HC0) standard error 0.071 และ 0.400

Coefficient เท่าเดิม แต่ standard error ต่างกัน

ทีมวิจัยศึกษาความสัมพันธ์ระหว่างการได้รับ treatment กับจำนวนครั้งที่ผู้ป่วยเข้ารับการรักษาในโรงพยาบาล outcome เป็นจำนวนนับ (count) ทีมจึงใช้ Poisson regression ซึ่งเป็น regression มาตรฐานสำหรับข้อมูลจำนวนนับ โดย model ลอการิทึมของ mean count ทีม run model สองครั้ง ครั้งแรกแบบ conventional และอีกครั้งด้วย robust standard error หลังวิเคราะห์ข้อมูล ได้ผลลัพธ์ดังนี้

Model เดียว standard error สองค่า

ตัวเลขสมมติเพื่อประกอบคำอธิบาย ไม่ใช่ผลการวิเคราะห์จากข้อมูลผู้ป่วยจริง SE คือ standard error และ CI คือ confidence interval
ผลลัพธ์Conventional PoissonPoisson with robust SE
Coefficient (β̂)0.400.40
Standard error0.100.19
95% CI ของ β̂0.204 ถึง 0.5960.028 ถึง 0.772
Exponentiated coefficient, exp(β̂)1.491.49
95% CI ของ exp(β̂)1.23 ถึง 1.811.03 ถึง 2.16

Variance และ standard error: ความแปรปรวนสองระดับ

ก่อนจะเข้าใจ robust standard error ต้องแยกความแปรปรวนออกเป็นสองระดับ

ระดับที่ 1: variance ของ outcome บอกว่าค่าของ outcome แตกต่างกันระหว่างผู้ป่วยมากเพียงใด ตัวอย่างเช่น ผู้ป่วยบางคนเข้ารับการรักษาในโรงพยาบาลเพียง 1 ครั้ง ขณะที่บางคนเข้า 9 ครั้ง

ระดับที่ 2: variance ของ estimator บอกว่าหากเก็บข้อมูลจากประชากรเดิมซ้ำหลายครั้ง แล้วประมาณ coefficient ใหม่ในแต่ละครั้ง ค่าประมาณเหล่านั้นจะเปลี่ยนแปลงมากเพียงใด standard error คือรากที่สองของ variance นี้

$$\mathrm{SE}(\hat\beta) = \sqrt{\mathrm{Var}(\hat\beta)}$$

ลองจินตนาการว่าเรามีประชากรผู้ป่วยจำนวนมาก และทำการศึกษาเดิมซ้ำ 1,000 ครั้ง โดยสุ่มผู้ป่วยกลุ่มใหม่ในแต่ละครั้ง ถ้า regression coefficient จากแต่ละการศึกษามีค่าใกล้เคียงกัน แสดงว่า estimator มีความแม่นยำสูงและ standard error ต่ำ แต่ถ้า coefficient แกว่งมากระหว่างการศึกษา แสดงว่า standard error สูง นี่คือความหมายของ sampling uncertainty

ประเด็นสำคัญคือ variance ของ outcome กับ standard error ของ coefficient ไม่ใช่สิ่งเดียวกัน แต่ความผันผวนของ outcome มีผลต่อความไม่แน่นอนในการประมาณ coefficient

ผู้ป่วยสองกลุ่มที่มี mean เท่ากัน

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

ตัวอย่างคำนวณด้วยมือ: จำนวนครั้งที่เข้ารับการรักษาในโรงพยาบาลของผู้ป่วยแต่ละคน
ผู้ป่วยDataset ADataset B
คนที่ 141
คนที่ 251
คนที่ 359
คนที่ 469
Mean55

Mean เท่ากัน แต่การกระจายต่างกัน

ทั้งสองชุดมี mean เท่ากัน

$$\text{Dataset A: } (4 + 5 + 5 + 6)/4 = 5$$$$\text{Dataset B: } (1 + 1 + 9 + 9)/4 = 5$$

แต่ลองสังเกตการกระจายของข้อมูล ใน Dataset A ผู้ป่วยมีจำนวนครั้งที่เข้ารับการรักษาใกล้เคียงกัน ในทางกลับกัน Dataset B มีผู้ป่วยบางคนที่เข้ารับการรักษาน้อยมาก และบางคนมากเป็นพิเศษ ดังนั้น แม้ mean จะเท่ากัน แต่ variance ของ Dataset B สูงกว่ามาก

คำนวณ variance ที่สังเกตได้

เพื่อให้สอดคล้องกับ sandwich estimator แบบ HC0 ที่ใช้ในหัวข้อถัดไป เราจะใช้ empirical variance ซึ่งหารด้วย $n$ ไม่ใช่ sample variance ที่หารด้วย $n - 1$ สูตรคือ

$$\text{Variance} = \sum_i (Y_i - \text{Mean})^2 / n$$

สำหรับ Dataset A ระยะห่างกำลังสองจาก mean คือ 1, 0, 0 และ 1

$$\text{Variance} = (1 + 0 + 0 + 1)/4 = 0.5$$

สำหรับ Dataset B ทุกค่าอยู่ห่างจาก mean 4 หน่วย ระยะห่างกำลังสองแต่ละค่าจึงเป็น 16

$$\text{Variance} = (16 + 16 + 16 + 16)/4 = 16$$
นี่คือจุดเริ่มต้นของความแตกต่างระหว่าง conventional standard error กับ robust standard error
ปริมาณDataset ADataset B
Mean55
Empirical variance0.516
ลักษณะข้อมูลกระจายน้อยกระจายมาก

Poisson regression มองข้อมูลทั้งสองชุดอย่างไร

เริ่มจาก Poisson regression ที่มีเฉพาะ intercept model คือ

$$\log(\mu_i) = \beta_0$$

โดยที่ $\mu_i = E(Y_i)$ คือจำนวนนับที่คาดไว้ของผู้ป่วยคนที่ $i$ และ log คือ natural logarithm ใน model นี้ไม่มี predictors อื่นนอกจาก intercept เมื่อ fit model จะได้

$$\hat\beta_0 = \log(\text{Mean})$$

เนื่องจากทั้งสองชุดมี mean เท่ากับ 5 จึงได้

$$\hat\beta_0 = \log(5) = 1.609$$

ดังนั้น regression coefficient ของ Dataset A และ Dataset B เท่ากัน

แต่ conventional Poisson model ยังมีสมมติฐานเกี่ยวกับการแจกแจงของข้อมูลว่า

$$\mathrm{Var}(Y_i) = E(Y_i) = \mu_i$$

นั่นหมายความว่า เมื่อ model คาดว่า mean เท่ากับ 5 ก็จะใช้ variance เท่ากับ 5 ด้วย ลองเปรียบเทียบกับข้อมูลที่พบจริง

ตัวอย่างคำนวณด้วยมือ: variance ที่ Poisson model สมมติ เทียบกับ variance ในข้อมูล
ปริมาณDataset ADataset B
Mean55
Poisson assumed variance55
Observed empirical variance0.516

Underdispersion และ overdispersion

เมื่อเทียบกับ Poisson variance assumption แล้ว Dataset A มีลักษณะ underdispersion คือกระจายน้อยกว่าที่ model สมมติ ส่วน Dataset B มีลักษณะ overdispersion คือกระจายมากกว่าที่ model สมมติ อย่างไรก็ตาม conventional Poisson model คำนวณ model-based standard error ภายใต้สมมติฐานเดียวกันทั้งสองชุด คือ variance เท่ากับ mean

คำนวณ conventional Poisson standard error ด้วยตัวเลขจริง

conventional standard error คือ model-based standard error ซึ่งได้จาก variance ที่ model สมมติไว้ สำหรับ Poisson regression ที่มีเฉพาะ intercept สามารถคำนวณ model-based standard error ของ estimated log mean ได้ดังนี้

$$\mathrm{SE}_{\text{Poisson}} = \sqrt{1/(n \times \hat\mu)}$$

ในที่นี้ $n = 4$ คน และ $\hat\mu = 5$

  1. แทนค่า

    \[ \mathrm{SE}_{\text{Poisson}} = \sqrt{1/(4 \times 5)} \]

    ผู้ป่วยสี่คน และ fitted mean เท่ากับ 5

  2. คำนวณค่าในราก

    \[ \mathrm{SE}_{\text{Poisson}} = \sqrt{0.05} \]

    1 หารด้วย 20 ได้ 0.05

  3. ถอดรากที่สอง

    \[ \mathrm{SE}_{\text{Poisson}} = 0.224 \]

    นี่คือ standard error ของ estimated log mean

ผลลัพธ์: conventional standard error เท่ากับ 0.224 ทั้ง Dataset A และ Dataset B เพราะ mean ของทั้งสองชุดเท่ากัน

ทำไมทั้งสองชุดได้ standard error เท่ากัน?

ตรงนี้เกิดคำถามสำคัญ Dataset A มีข้อมูลค่อนข้างใกล้กัน แต่ Dataset B มีข้อมูลกระจายมากกว่าอย่างชัดเจน เหตุใดทั้งสองชุดจึงมี model-based standard error เท่ากัน?

คำตอบคือ conventional Poisson ไม่ได้ใช้ observed squared residuals เพื่อประมาณ standard error แบบที่ robust sandwich estimator ใช้ แต่ใช้ variance structure ที่กำหนดโดย Poisson model ทั้งนี้ residual คือจำนวนนับที่สังเกตได้ของผู้ป่วยลบด้วย fitted mean คือ $Y_i - \hat\mu$

Robust standard error แก้ปัญหานี้อย่างไร

Robust standard error ยังคงใช้ Poisson regression เพื่อประมาณ coefficient เหมือนเดิม สิ่งที่เปลี่ยนคือวิธีคำนวณ variance ของ coefficient โดยอาศัย score contributions ที่สังเกตได้ score contribution คือส่วนของผู้ป่วยหนึ่งคนในสมการที่ซอฟต์แวร์แก้เพื่อหา coefficient ใน model ที่มีเฉพาะ intercept ค่านี้ก็คือ residual ของผู้ป่วย $Y_i - \hat\mu$

สำหรับ intercept-only Poisson model สามารถเขียน robust standard error ในรูปพื้นฐาน ซึ่งเรียกว่า HC0 ได้เป็น

$$\mathrm{SE}_{\text{Robust}} = \sqrt{\sum_i (Y_i - \hat\mu)^2} \,/\, (n \times \hat\mu)$$

ลองคำนวณสำหรับข้อมูลทั้งสองชุด

Robust (HC0) standard error ของข้อมูลทั้งสองชุด

ตัวอย่างคำนวณด้วยมือ Dataset A คือ 4, 5, 5 และ 6 ส่วน Dataset B คือ 1, 1, 9 และ 9 แต่ละชุดมีผู้ป่วยสี่คน และ fitted mean เท่ากับ 5

  1. Dataset A: ผลรวมของ squared residuals

    \[ 1 + 0 + 0 + 1 = 2 \]

    residuals คือ -1, 0, 0 และ +1

  2. Dataset A: robust standard error

    \[ \mathrm{SE}_{\text{Robust}} = \sqrt{2}/(4 \times 5) = 0.071 \]

    ค่านี้ต่ำกว่า conventional standard error 0.224

  3. Dataset B: ผลรวมของ squared residuals

    \[ 16 + 16 + 16 + 16 = 64 \]

    residuals คือ -4, -4, +4 และ +4

  4. Dataset B: robust standard error

    \[ \mathrm{SE}_{\text{Robust}} = \sqrt{64}/(4 \times 5) = 0.400 \]

    ค่านี้สูงกว่า conventional standard error 0.224

ผลลัพธ์: robust (HC0) standard error ของ Dataset A คือ 0.071 และของ Dataset B คือ 0.400

หัวใจของ robust standard error

ตัวอย่างคำนวณด้วยมือ ชุดละผู้ป่วยสี่คน SE คือ standard error
ผลลัพธ์Dataset ADataset B
Poisson coefficient1.6091.609
Conventional SE0.2240.224
Robust SE (HC0)0.0710.400

Sandwich variance estimator คืออะไร?

ตอนนี้เราสามารถกลับมาทำความเข้าใจสูตรของ sandwich estimator ได้แล้ว [1, 2] ใน regression model ที่มี predictors หลายตัว variance ของ coefficient ไม่ได้เป็นตัวเลขเดียวอีกต่อไป แต่เป็น variance-covariance matrix สูตร sandwich คือ

$$\hat V_{\text{Robust}}(\hat\beta) = A^{-1} B A^{-1}$$

สูตรนี้มีสองส่วนสำคัญ คือ bread $A$ และ meat $B$ ตัวอักษรในที่นี้เป็นชื่อส่วนของสูตร ไม่ใช่ Dataset A และ Dataset B

Bread: A

สำหรับ Poisson regression

$$A = \sum_i \hat\mu_i X_i X_i^\top$$

ในที่นี้ $X_i$ คือ column ของค่า predictors ของผู้ป่วยคนที่ $i$ โดยมีค่า 1 สำหรับ intercept และ $X_i^\top$ คือค่าชุดเดียวกันที่เรียงเป็นแถว Bread อธิบายความไวของ estimating equations ซึ่งคือสมการที่ซอฟต์แวร์แก้เพื่อหา coefficients ต่อการเปลี่ยนแปลงของ coefficients พูดง่าย ๆ คือ บอกว่าเมื่อ coefficient เปลี่ยนเพียงเล็กน้อย จะส่งผลต่อ model มากเพียงใด

Meat: B

เริ่มจาก score contribution ของ observation ที่ $i$

$$s_i = X_i (Y_i - \hat\mu_i)$$

จากนั้นคำนวณ

$$B = \sum_i s_i s_i^\top$$

หรือเขียนอีกแบบได้ว่า

$$B = \sum_i X_i X_i^\top (Y_i - \hat\mu_i)^2$$

Meat สะท้อนความผันผวนของ score contributions ที่สังเกตได้จากข้อมูล

รวมเป็น sandwich

นำ bread และ meat มาคำนวณร่วมกัน

$$\hat V_{\text{Robust}} = A^{-1} B A^{-1}$$

จากนั้นนำค่าบนแนวทแยงของ matrix มาถอดรากที่สอง เพื่อหา standard error ของ coefficient แต่ละตัว

$$\mathrm{SE}(\hat\beta_j) = \sqrt{[\hat V_{\text{Robust}}]_{jj}}$$

ย้อนกลับไปยัง Dataset B

ตัวอย่างคำนวณด้วยมือ เนื่องจากมีเพียง intercept จึงมี $X_i = 1$ สำหรับผู้ป่วยทุกคน $A$ และ $B$ จึงเป็นตัวเลขเดียวที่คำนวณได้โดยตรง

  1. Bread

    \[ A = 4 \times 5 = 20 \]

    จำนวนผู้ป่วยคูณ fitted mean

  2. Meat

    \[ B = 16 + 16 + 16 + 16 = 64 \]

    ผลรวมของ squared residuals

  3. Sandwich variance

    \[ \hat V_{\text{Robust}} = 64/20^2 = 0.16 \]

    meat หารด้วย bread ยกกำลังสอง

  4. Robust standard error

    \[ \mathrm{SE}_{\text{Robust}} = \sqrt{0.16} = 0.40 \]

    ได้ค่าเดียวกับสูตรในหัวข้อก่อนหน้า

ผลลัพธ์: sandwich ให้ robust standard error ของ Dataset B เท่ากับ 0.40 ซึ่งเป็นค่า HC0 เดียวกับที่คำนวณไว้ก่อนหน้า

ทำไมจึงเรียกว่า sandwich

$A^{-1}$ เป็น bread สองด้าน ส่วน $B$ เป็น meat ที่อยู่ตรงกลาง ใน regression ที่มี predictors หลายตัว วิธีการเหมือนเดิม แต่เปลี่ยนจากตัวเลขธรรมดาเป็น matrix operations

สิ่งสำคัญคือ robust standard error ของ coefficient แต่ละตัวไม่จำเป็นต้องเปลี่ยนไปด้วยอัตราส่วนเดียวกัน เพราะ score contributions และ predictor structure แตกต่างกัน

Robust standard error ส่งผลต่อ confidence interval อย่างไร

กลับมายัง regression coefficient เท่ากับ 0.40 จากตารางตอนต้น สูตร 95% confidence interval บน log scale โดยใช้ large-sample normal approximation คือ

$$\text{95\% CI} = \hat\beta \pm 1.96 \times \mathrm{SE}(\hat\beta)$$
  1. Conventional Poisson: lower limit

    \[ 0.40 - 1.96 \times 0.10 = 0.204 \]

    ในที่นี้ standard error เท่ากับ 0.10

  2. Conventional Poisson: upper limit

    \[ 0.40 + 1.96 \times 0.10 = 0.596 \]

    ช่วงนี้อยู่ระหว่าง 0.204 ถึง 0.596

  3. Robust Poisson: lower limit

    \[ 0.40 - 1.96 \times 0.19 = 0.028 \]

    ในที่นี้ standard error เท่ากับ 0.19

  4. Robust Poisson: upper limit

    \[ 0.40 + 1.96 \times 0.19 = 0.772 \]

    ช่วงนี้อยู่ระหว่าง 0.028 ถึง 0.772

  5. Exponentiate coefficient

    \[ \exp(0.40) \approx 1.49 \]

    ได้ค่าเดียวกันทั้งสองวิธี

ผลลัพธ์: ทั้งสองวิธีมี point estimate เท่ากับ 0.40 หรือ 1.49 หลัง exponentiate แต่ robust standard error ที่สูงกว่าทำให้ confidence interval กว้างขึ้น

เชื่อมโยงกับ modified Poisson regression สำหรับ binary outcome

ใน clinical epidemiology เรามักสนใจ outcome แบบ binary เช่น neonatal morbidity ภายใน 28 วัน โดย $Y = 1$ เมื่อเกิด outcome และ $Y = 0$ เมื่อไม่เกิด สำหรับ binary outcome ที่มีความน่าจะเป็นเกิดเหตุการณ์เท่ากับ $p$ expected value คือ

$$E(Y) = p$$

แต่ variance ภายใต้ Bernoulli distribution หรือ Bernoulli variance คือ

$$\mathrm{Var}(Y) = p(1 - p)$$

ในขณะที่ Poisson variance assumption คือ

$$\mathrm{Var}(Y) = p$$

ลองสมมติว่า risk เท่ากับ 0.40 Poisson จะใช้ variance เท่ากับ 0.40 แต่ Bernoulli variance คือ

$$0.40 \times (1 - 0.40) = 0.24$$

นี่หมายความว่า Poisson variance assumption ไม่ตรงกับ Bernoulli variance ที่สำคัญ ความแตกต่างนี้ไม่ได้เป็น overdispersion เหมือนตัวอย่าง count outcome ก่อนหน้า แต่ Bernoulli variance ต่ำกว่า Poisson variance เมื่อ $0 < p < 1$

Model และสิ่งที่ประมาณได้

Model คือ

$$\log[P(Y = 1 \mid X)] = \beta_0 + \beta_1 X$$

เมื่อ $X$ เป็น binary exposure และ model เหมาะสม

$$\mathrm{RR} = \exp(\beta_1)$$

ดังนั้น modified Poisson regression ไม่ได้เปลี่ยน outcome ให้เป็น count data จริง ๆ แต่ใช้ Poisson estimating framework กับ binary outcome

วิธีนี้เป็นอีกทางเลือกหนึ่งของ log-binomial regression ซึ่ง model ลอการิทึมของ risk โดยตรง โดยเฉพาะเมื่อ log-binomial model มีปัญหา convergence อย่างไรก็ตาม robust standard error ไม่ได้แก้ปัญหา predicted risk ที่อาจมากกว่า 1 หรือ mean model ที่ระบุไม่ถูกต้อง บทความเรื่อง risk ratio เมื่อ outcome พบบ่อยเปรียบเทียบสอง model นี้อย่างละเอียดกว่า

Stata คำนวณ robust standard error อย่างไร

คำสั่งด้านล่างเขียนสำหรับตัวอย่างการเข้ารับการรักษาในโรงพยาบาลตอนต้น ได้แก่ จำนวนครั้งที่เข้ารับการรักษา ตัวแปร treatment และอายุ Stata กำหนดวิธีประมาณ variance ผ่าน option vce() [4]

Conventional Poisson และ Poisson with robust standard error

A. Conventional Poisson

poisson admission_count i.treatment c.age

B. Poisson with robust standard error

poisson admission_count i.treatment c.age, vce(robust)

เมื่อใช้ข้อมูลและ model specification ชุดเดียวกัน ทั้งสองคำสั่งจะให้ regression coefficients เท่ากัน แต่ standard errors, confidence intervals และ p-values อาจแตกต่างกัน คำสั่ง vce(robust) ขอให้ Stata ใช้ robust variance-covariance estimator

Modified Poisson สำหรับ binary outcome

C. ตัวอย่างการศึกษาความสัมพันธ์ระหว่าง inadequate antenatal care (ANC) กับ neonatal morbidity

poisson morbidity28 i.inadequate_anc i.age_group i.low_ses i.coverage_gap, irr vce(robust)

หากเป็น binary outcome ที่ติดตามทุกคนครบช่วงเวลาเดียวกัน และ log-risk model เหมาะสม ค่าที่ได้จากการ exponentiate coefficient สามารถตีความเป็น adjusted risk ratio แม้ Stata จะแสดงหัวคอลัมน์ว่า IRR ตาม option irr ซึ่งย่อมาจาก incidence rate ratio แต่ในบริบทนี้ effect measure คือ risk ratio ไม่ใช่ incidence rate ratio

ในทางกลับกัน ถ้า outcome เป็นจำนวนเหตุการณ์ต่อ person-time และใส่ exposure time เข้าไปใน model การ exponentiate coefficient จะมีความหมายเป็น rate ratio

ดังนั้นชื่อ effect measure ต้องมาจาก outcome, time denominator และ estimand ไม่ใช่เพียงชื่อ option ใน Stata ทั้งนี้ estimand คือปริมาณที่การศึกษาต้องการประมาณอย่างเจาะจง

Cluster-robust standard error

D. กรณีมี repeated measurements หลายครั้งต่อผู้ป่วย

poisson admission_count i.treatment c.age, vce(cluster patient_id)

คำสั่งนี้คำนวณ cluster-robust standard error โดยรวม score contributions ภายในผู้ป่วยคนเดียวกันก่อนสร้างส่วน meat ของ sandwich estimator การอนุมานแบบ large-sample (asymptotic) ต้องมีจำนวน independent clusters เพียงพอ และเมื่อจำนวน clusters น้อย ควรพิจารณาวิธีแก้ small-sample bias แนวคิดเดียวกันนี้ขยาย modified Poisson regression ไปยัง correlated binary data [5] ตอนที่ 2 ของชุดบทความนี้อธิบายว่าเมื่อใดควรใช้ option vce() แบบใด

Robust standard error ทำอะไรได้ และทำอะไรไม่ได้

Robust standard error ช่วยได้Robust standard error ไม่ได้แก้ไข
ลดการพึ่งพา Poisson variance assumptionMean model ที่ระบุผิด
ประมาณความไม่แน่นอนของ coefficient ภายใต้เงื่อนไขที่เหมาะสมConfounding
ใช้กับ modified Poisson สำหรับ binary outcomeSelection bias
ใช้แบบ cluster-robust เพื่อรองรับ within-cluster correlationCorrelation ที่ไม่ได้ระบุโครงสร้าง cluster อย่างเหมาะสม
ช่วยให้ large-sample inference เหมาะสมขึ้นปัญหาการอนุมานจาก sample size ที่เล็กมาก
ใช้กับ count mean model ได้ไม่ได้สร้าง negative binomial distribution

Count outcome ที่มี overdispersion: วิเคราะห์เพื่ออะไร?

โดยเฉพาะเมื่อ count outcome มี overdispersion ผู้วิจัยต้องพิจารณาว่าเป้าหมายของการวิเคราะห์คืออะไร ถ้าสนใจ conditional mean association และ mean model ถูกต้อง Poisson regression with robust standard error อาจเพียงพอสำหรับ inference แต่ถ้าสนใจ count distribution, predictive distribution หรือการอธิบาย variance structure อาจพิจารณา negative binomial regression ซึ่งเป็น count model ที่ยอมให้ variance มากกว่า mean หรือ model อื่นที่เหมาะสม

Robust standard error ไม่ได้ทำให้ variance ของข้อมูลลดลง และไม่ได้ทำให้ regression model กลายเป็น model ที่ถูกต้องโดยอัตโนมัติ [6]

ตัวอย่างการเขียน Methods ใน manuscript

สำหรับ binary outcome ใน cohort study สามารถเขียนได้ว่า

"We estimated adjusted risk ratios using modified Poisson regression with a log link and robust sandwich variance estimation to obtain 95% confidence intervals."

ควรระบุ adjustment variables, การจัดการ missing data, follow-up และ clustering เพิ่มเติมตาม study design

สำหรับ count outcome อาจเขียนว่า

"We fitted Poisson regression models with robust sandwich standard errors to estimate associations with the expected event count."

หากใส่ person-time offset และสนใจ incidence rate ratio ต้องระบุด้วยว่าใช้ log person-time เป็น offset ทั้งนี้ offset คือพจน์ที่ใส่ไว้ใน model โดยกำหนด coefficient ไว้ที่ 1 เพื่อให้ model อธิบายจำนวนเหตุการณ์ต่อหน่วยเวลา

บทสรุป: point estimation และ variance estimation

สิ่งสำคัญที่สุดของ robust standard errors ไม่ใช่การจำสูตร sandwich แต่เป็นการเข้าใจว่า statistical model มีสองส่วนที่ต้องแยกออกจากกัน

ส่วนแรก: point estimation model ใช้ข้อมูลเพื่อประมาณ regression coefficient ซึ่งบอกขนาดและทิศทางของความสัมพันธ์

ส่วนที่สอง: variance estimation เราใช้วิธีทางสถิติเพื่อประเมินว่า regression coefficient ที่ประมาณได้นั้นมีความไม่แน่นอนมากเพียงใด

Conventional Poisson ใช้ variance assumption ของ Poisson model ส่วน robust Poisson ใช้ sandwich estimator เพื่อประมาณความไม่แน่นอนจาก observed score contributions โดยยังคง mean model เดิมไว้

ดังนั้น robust standard error จึงไม่ใช่เทคนิคที่ทำให้ข้อมูลดีขึ้น แต่เป็นเครื่องมือที่ช่วยให้การอนุมานเกี่ยวกับ coefficient เหมาะสมขึ้นเมื่อ variance structure ของ model ไม่ตรงกับข้อมูล ภายใต้สมมติฐานและเงื่อนไขที่จำเป็น

ความเข้าใจผิดที่พบบ่อยและวิธีแก้

  • "Robust standard error คือ model ใหม่"

    coefficient ในตารางตอนต้นยังคงเป็น 0.40 ไม่ว่าจะใช้ standard error แบบใด และข้อมูลทั้งสองชุดในตัวอย่างคำนวณด้วยมือก็ยังคงมี coefficient 1.609

    วิธีแก้: Robust standard error ใช้ model เดิมและ coefficient เดิม เปลี่ยนเฉพาะวิธีประมาณ variance-covariance matrix

  • "Robust standard error สูงกว่าเสมอ"

    ใน Dataset A robust standard error เท่ากับ 0.071 ต่ำกว่า conventional standard error 0.224

    วิธีแก้: Robust standard error ตามการกระจายที่สังเกตได้ จะเล็กลงเมื่อข้อมูลแปรปรวนน้อยกว่าที่ model สมมติ และใหญ่ขึ้นเมื่อแปรปรวนมากกว่า

  • "Poisson regression ของ binary outcome มีปัญหา overdispersion"

    ที่ risk เท่ากับ 0.40 Bernoulli variance เท่ากับ 0.24 ต่ำกว่า Poisson variance ที่ 0.40

    วิธีแก้: ความไม่ตรงกันเป็นไปในทิศตรงข้าม Bernoulli variance ต่ำกว่า Poisson variance เสมอเมื่อ risk อยู่ระหว่าง 0 ถึง 1 และ robust standard error เข้ามาแทนที่ในการอนุมาน

  • "Stata แสดง IRR ผลจึงเป็น incidence rate ratio"

    option irr กำหนดเพียงชื่อหัวคอลัมน์ เมื่อ outcome เป็น binary และทุกคนติดตามเท่ากัน ค่าที่ได้จากการ exponentiate coefficient คือ risk ratio

    วิธีแก้: ตั้งชื่อ effect measure จาก outcome, time denominator และ estimand ไม่ใช่จาก option ของ Stata

  • "Robust standard error ทำให้ model ถูกต้อง"

    Robust standard error ไม่ได้แก้ mean model ที่ระบุผิด confounding, selection bias, โครงสร้าง cluster ที่ไม่ได้ระบุ หรือ sample ที่เล็กมาก และไม่ได้เปลี่ยน Poisson model ให้เป็น negative binomial model

    วิธีแก้: ใช้ robust standard error สำหรับความไม่แน่นอนของ coefficient ส่วน mean model, confounding และ study design ต้องจัดการแยกต่างหาก

อภิธานศัพท์

Poisson regression
regression สำหรับข้อมูลจำนวนนับ ซึ่ง model ลอการิทึมของ mean count และสมมติว่า variance เท่ากับ mean
variance-covariance matrix (เมทริกซ์ความแปรปรวนและความแปรปรวนร่วม)
matrix ที่เก็บ variance ที่ประมาณได้ของ coefficient แต่ละตัวไว้บนแนวทแยง และเก็บ covariance ระหว่าง coefficients ไว้นอกแนวทแยง standard errors คือรากที่สองของค่าบนแนวทแยง
sampling uncertainty
ค่าประมาณจะแปรผันมากเพียงใดระหว่างการศึกษาซ้ำบน sample ใหม่จากประชากรเดิม standard error เป็นตัววัดสิ่งนี้
empirical variance
ค่าเฉลี่ยของระยะห่างกำลังสองของค่าสังเกตจาก mean โดยหารด้วย n ไม่ใช่ n - 1
underdispersion (การกระจายน้อยกว่าที่แบบจำลองกำหนด)
count outcome กระจายน้อยกว่าที่ Poisson model สมมติ คือ variance ต่ำกว่า mean
overdispersion (การกระจายเกินกว่าที่แบบจำลองกำหนด)
count outcome กระจายมากกว่าที่ Poisson model สมมติ คือ variance มากกว่า mean
model-based standard error (standard error แบบอิงแบบจำลอง)
standard error ที่ได้จาก variance ซึ่ง model สมมติไว้ เช่น variance เท่ากับ mean ใน Poisson regression
score contribution (score contribution ของแต่ละหน่วยสังเกต)
ส่วนของ observation หนึ่งใน estimating equations คือค่า predictors คูณ residual เมื่อมีเฉพาะ intercept ก็คือ residual นั้นเอง
HC0
robust (sandwich) standard error ในรูปพื้นฐาน ที่ยังไม่มี finite-sample correction ซึ่งซอฟต์แวร์บางชนิดใช้
sandwich estimator
การประมาณ variance ในรูป inverse bread, meat, inverse bread ซึ่งคง coefficient เดิมไว้ และอ่านการกระจายจาก observed score contributions
bread (ส่วน bread ของ sandwich estimator)
ตัวประกอบด้านนอกของ sandwich คือ A บอกความไวของ estimating equations ต่อ coefficients และปรากฏในรูป inverse ทั้งสองข้างของ meat
meat (ส่วน meat ของ sandwich estimator)
ตัวประกอบตรงกลางของ sandwich คือ B ผลรวมของ squared score contributions ซึ่งอ่านจากข้อมูล ไม่ใช่จาก variance ที่ model สมมติ
Bernoulli variance (ความแปรปรวนแบบ Bernoulli)
variance p(1 - p) ของ outcome แบบใช่หรือไม่ใช่ที่มีความน่าจะเป็น p ซึ่งต่ำกว่าค่า Poisson คือ p เสมอเมื่อ p อยู่ระหว่าง 0 ถึง 1
modified Poisson regression (Poisson ดัดแปลง)
Poisson log-link mean model สำหรับ binary outcome ที่ใช้ robust sandwich standard error เพื่อประมาณ risk ratio
log-binomial regression (แบบจำลอง log-binomial)
regression สำหรับ outcome แบบใช่หรือไม่ใช่ที่ model ลอการิทึมของ risk โดยตรง และอาจไม่ converge
rate ratio (อัตราส่วนอัตราอุบัติการณ์)
อัตราส่วนของจำนวนเหตุการณ์ต่อ person-time ระหว่างสองกลุ่ม Stata ติดป้ายว่า IRR ซึ่งย่อมาจาก incidence rate ratio
estimand (ปริมาณเป้าหมายของการประมาณ)
ปริมาณที่การศึกษาต้องการประมาณอย่างเจาะจง รวมถึงประชากร การเปรียบเทียบ และ summary measure
cluster-robust standard error (cluster-robust standard error)
sandwich standard error ที่ส่วน meat รวม score contributions ภายในแต่ละ cluster เช่นผู้ป่วยหนึ่งคน ก่อนยกกำลังสอง และต้องมี independent clusters มากพอ
negative binomial regression (negative binomial regression)
count model ที่มี mean structure แบบ Poisson และมี dispersion parameter เพิ่ม ทำให้ variance มากกว่า mean ได้
offset
พจน์ที่ใส่ไว้ใน model โดยกำหนด coefficient ไว้ที่ 1 เช่น log ของ person-time เพื่อให้ model อธิบายจำนวนเหตุการณ์ต่อหน่วยเวลา

เอกสารอ้างอิง

  1. Huber PJ. The behavior of maximum likelihood estimates under nonstandard conditions. In: Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, Vol. 1. University of California Press; 1967. p. 221-233. https://projecteuclid.org/euclid.bsmsp/1200512988
  2. 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
  3. Zou G. A modified Poisson regression approach to prospective studies with binary data. Am J Epidemiol. 2004;159(7):702-706. doi:10.1093/aje/kwh090 https://doi.org/10.1093/aje/kwh090
  4. StataCorp. Stata Base Reference Manual: [R] vce_option. Stata documentation at stata.com (PDF). Accessed 2026-10-11. https://www.stata.com/manuals/rvce_option.pdf
  5. Zou GY, Donner A. Extension of the modified Poisson regression model to prospective studies with correlated binary data. Stat Methods Med Res. 2013;22(6):661-670. doi:10.1177/0962280211427759 https://doi.org/10.1177/0962280211427759
  6. Freedman DA. On the so-called 'Huber sandwich estimator' and 'robust standard errors'. Am Stat. 2006;60(4):299-302. doi:10.1198/000313006X152207 https://doi.org/10.1198/000313006X152207

ประเด็นสำคัญ

  • การประมาณ coefficient กับการประมาณความไม่แน่นอนของ coefficient เป็นคนละกระบวนการ robust standard error คง model และ coefficient เดิมไว้ และเปลี่ยนเฉพาะการประมาณ variance-covariance
  • เมื่อ mean เท่ากัน Dataset A และ B มี coefficient 1.609 และ conventional standard error 0.224 เท่ากัน แต่ robust (HC0) standard error เป็น 0.071 และ 0.400 เพราะการกระจายที่สังเกตได้ต่างกัน
  • Sandwich คือ A⁻¹BA⁻¹ ส่วน bread A มาจาก model ส่วน meat B มาจาก observed score contributions และ standard error แต่ละค่าคือรากที่สองของค่าบนแนวทแยง
  • Modified Poisson regression ใช้ Poisson log-link model สำหรับ risk ratio และใช้ robust sandwich variance สำหรับ inference ชื่อ effect measure ต้องมาจาก outcome, time denominator และ estimand ไม่ใช่จาก option ของ Stata
  • Robust standard error เปลี่ยนวิธีประมาณความไม่แน่นอน ไม่ได้เปลี่ยน coefficient ที่ประมาณได้ และไม่ได้แก้ mean model ที่ผิด confounding หรือ selection bias

อ่านต่อในวิกิ: [[heteroskedasticity-robust-standard-errors-th]] [[risk-regression-models-epidemiology-th]] [[choosing-vce-robust-cluster-standard-errors-th]]

0
ถึงนักอ่านชาวไทยและต่างชาติทำความเข้าใจบริบททางการแพทย์ของผมอ่านต่อ →ถึงนักอ่านชาวไทยและต่างชาติทำความเข้าใจเนื้อหาของผมที่นอกเหนือจากการแพทย์อ่านต่อ →

ความคิดเห็น

ยังไม่มีความคิดเห็น มาเป็นคนแรกกันเลย

เข้าสู่ระบบเพื่อแสดงความคิดเห็น