logit, ologit และ clogit: link function เดียว แต่ตอบคำถามต่างกันสามแบบ

On this page
Read the English version
บทคัดย่อ
logistic regression, ordinal logistic regression และ conditional logistic regression ใช้ link เดียวกันคือ logit สัมประสิทธิ์ที่ยกกำลังด้วยเลขชี้กำลังจึงเป็น odds ratio ชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่างเป็นตัวกำหนดแบบจำลอง ในทารกแรกเกิดจำลอง 500 คน การฝากครรภ์ไม่เพียงพอมี odds ratio เท่ากับ 2.01 สำหรับระดับภาวะแทรกซ้อนที่สูงขึ้น ค่านี้ปรับด้วยการคลอดก่อนกำหนด ซึ่งอาจเป็นตัวแปรสื่อกลาง (mediator) ดังนั้น 2.01 จึงไม่รวมผลใด ๆ ที่ส่งผ่านการคลอดก่อนกำหนด แบบจำลองแบบลำดับสมมติว่ามี odds ratio ค่าเดียวที่ทุกจุดตัดระหว่างระดับ Brant test พบหลักฐานของการเบี่ยงเบนจากข้อสมมตินี้สำหรับการคลอดก่อนกำหนด (P = 0.003) แต่ไม่พบสำหรับการฝากครรภ์ ในคู่ case-control ที่จับคู่ 20 คู่ logistic regression ทั่วไปให้ odds ratio 1.83 ส่วน conditional logistic regression ให้ 4.00 ทั้งสองค่าตอบคำถามต่างกัน คือเปรียบเทียบแบบรวมกลุ่มกับเปรียบเทียบภายในคู่ และแบบแผนการศึกษาเป็นตัวกำหนดว่าควรรายงานค่าใด มีเพียง 5 คู่ที่ต่างกันในสิ่งสัมผัส ค่า 4.00 จึงไม่แม่นยำ (ช่วงเชื่อมั่น 95% 0.45 ถึง 35.79) บทความนี้คำนวณค่าประมาณทั้งสองด้วยมือ แล้วแสดงโค้ด Stata และ R และตารางสำหรับตัดสินใจ
สองแบบแผน สามคำสั่ง และ odds ratio สองค่า
หน่วยทารกแรกเกิดแห่งหนึ่งต้องการทราบว่าการฝากครรภ์ไม่เพียงพอ (inadequate antenatal care, ANC) สัมพันธ์กับภาวะแทรกซ้อนในทารกแรกเกิดมากขึ้นหรือไม่ คำว่าไม่เพียงพอในที่นี้หมายถึงการตรวจระหว่างตั้งครรภ์ที่น้อยเกินไปหรือเริ่มช้าเกินไป ทีมศึกษาคำถามนี้ด้วยแบบแผนสองแบบพร้อมกัน หน่วยนี้ ทารก และตัวเลขทุกตัวด้านล่างเป็นข้อมูลจำลอง
ใน birth cohort ทารกแรกเกิด 500 คนแต่ละคนมีระดับภาวะแทรกซ้อนภายใน 28 วัน ได้แก่ ไม่มี เล็กน้อย ปานกลาง หรือรุนแรง ระดับเหล่านี้มีลำดับ แต่ช่วงห่างระหว่างระดับไม่จำเป็นต้องเท่ากัน อีกส่วนหนึ่งเป็นการศึกษาแบบ case-control ซึ่งจับคู่ทารกที่มีภาวะแทรกซ้อน 20 คน (case) กับกลุ่มควบคุม (control) คนละหนึ่งคน control แต่ละคนเป็นทารกที่ยังไม่มีภาวะแทรกซ้อนในวันที่ case ของมันมี
นักวิเคราะห์รันคำสั่ง Stata สามคำสั่ง คือ logistic, ologit และ clogit ในคู่ที่จับคู่ logistic ให้ odds ratio 1.83 และ clogit ให้ 4.00 ทารก สิ่งสัมผัส และผลลัพธ์เป็นชุดเดียวกัน ผู้ทบทวนจึงถามว่าตัวเลขไหนถูก
ตัวเลขทั้งสองคำนวณถูกต้อง แต่เปรียบเทียบคนละแบบ และมีเพียงการวิเคราะห์ที่จับคู่เท่านั้นที่สอดคล้องกับแบบแผนการศึกษา แบบจำลองทั้งสามใช้ส่วนประกอบเดียวกันคือ logit link ต่างกันที่ข้อสมมติเกี่ยวกับผลลัพธ์และเกี่ยวกับวิธีสุ่มตัวอย่างผู้คนในข้อมูล ก่อนถึงแบบจำลองแต่ละแบบ บทความนี้ระบุชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่างก่อน เพราะสองข้อเท็จจริงนี้เป็นตัวกำหนดคำสั่ง
logit link ที่แบบจำลองทั้งสามใช้ร่วมกัน
แบบจำลองทุกแบบในบทความนี้อธิบายความน่าจะเป็น $p$ เช่น ความน่าจะเป็นที่ทารกแรกเกิดมีภาวะแทรกซ้อน ความน่าจะเป็นถูกจำกัดอยู่ระหว่าง 0 ถึง 1 ขณะที่สมการถดถอยรับค่าใดก็ได้ ฟังก์ชันเชื่อมโยง (link function) เชื่อมสองสิ่งนี้เข้าด้วยกัน โดยแปลง $p$ ไปอยู่บนมาตรวัดที่ปรับแบบจำลองเชิงเส้นได้
logit link คือ log ของ odds:
$$\operatorname{logit}(p) = \ln\left(\frac{p}{1 - p}\right)$$ในที่นี้ $p/(1 - p)$ คือ odds และ $\ln$ คือลอการิทึมธรรมชาติ logit มีค่าตั้งแต่ลบอนันต์ถึงบวกอนันต์ ผลรวมของสัมประสิทธิ์คูณตัวแปรร่วม ซึ่งคือ ตัวทำนายเชิงเส้น (linear predictor) จึงรับค่าใดก็ได้ และยังแปลงกลับเป็นความน่าจะเป็นที่ถูกต้องได้
สำหรับสิ่งสัมผัสแบบสองค่า $x$ ที่รหัส 1 คือการฝากครรภ์ไม่เพียงพอ และ 0 คือเพียงพอ logistic regression ทั่วไปประมาณแบบจำลอง
$$\operatorname{logit} P(Y = 1 \mid x) = \beta_0 + \beta_1 x$$ในที่นี้ $Y = 1$ หมายถึงทารกที่มีผลลัพธ์ และ $\beta_0$ คือ log odds เมื่อ $x = 0$ ความชัน $\beta_1$ คือการเปลี่ยนแปลงของ log odds เมื่อ $x$ เปลี่ยนจาก 0 เป็น 1 ผลต่างของ log odds สองค่าคือ log ของอัตราส่วนของมัน ดังนั้น $\beta_1$ คือ log odds ratio และ $e^{\beta_1}$ คือ odds ratio [1]
คำสั่ง logit ของ Stata พิมพ์ $\beta_1$ ส่วนคำสั่ง logistic ประมาณแบบจำลองเดียวกันและพิมพ์ $e^{\beta_1}$ ฟังก์ชัน glm() ของ R ที่ใช้ family = binomial ใช้ logit link เป็นค่าเริ่มต้นและพิมพ์ $\beta_1$
logistic regression ทั่วไปสมมติสองอย่างเกี่ยวกับข้อมูล คือผลลัพธ์เป็นแบบสองค่า และหน่วยสังเกตเป็นอิสระต่อกัน หมายถึงผลลัพธ์ของคนหนึ่งไม่บอกอะไรเกี่ยวกับอีกคนเมื่อกำหนดตัวแปรร่วมแล้ว ระดับภาวะแทรกซ้อนที่มีลำดับฝืนข้อสมมติแรก และคู่ที่จับคู่ฝืนข้อสมมติหลัง การจับคู่ยังผูก control แต่ละคนไว้กับ case ของมัน จึงต้องเข้าไปอยู่ในการเปรียบเทียบเอง ไม่ใช่แค่ในค่าคลาดเคลื่อนมาตรฐาน
ระบุชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่างก่อน
การวิเคราะห์สามแบบในบทความนี้ใช้ logit link ร่วมกัน แต่ต้องใช้แบบจำลองต่างกัน การเลือกขึ้นกับสองข้อเท็จจริงที่ทราบก่อนวิเคราะห์ข้อมูลใด ๆ สิ่งแรกคือชนิดของผลลัพธ์ (outcome type) มีสองหมวด หรือหลายหมวดที่มีลำดับ สิ่งที่สองคือแบบแผนการสุ่มตัวอย่าง (sampling design) เป็นคนอิสระต่อกัน หรือเป็นชุดที่จับคู่ซึ่งประกอบด้วย case หนึ่งคนกับ control ของมัน
- ผลลัพธ์สองค่า คนอิสระต่อกัน: logistic regression ทั่วไป Stata ใช้
logisticR ใช้glm(..., family = binomial) - ผลลัพธ์หลายหมวดที่มีลำดับ คนอิสระต่อกัน: การถดถอยโลจิสติกสำหรับผลลัพธ์แบบลำดับ (ordinal logistic regression) Stata ใช้
ologitR ใช้MASS::polr() - ผลลัพธ์สองค่า (case หรือ control) ชุดที่จับคู่: conditional logistic regression Stata ใช้
clogitR ใช้survival::clogit()
logit ปรากฏในตำแหน่งต่างกันเล็กน้อยในแต่ละแบบจำลอง logistic regression สร้างแบบจำลอง logit ของความน่าจะเป็นหนึ่งค่า ordinal logistic regression สร้างแบบจำลอง logit หนึ่งค่าต่อจุดตัดระหว่างระดับแต่ละจุด ส่วน conditional logistic regression สร้างแบบจำลอง logit โดยให้ intercept แยกกันสำหรับชุดที่จับคู่แต่ละชุด แล้วกำจัด intercept เหล่านั้นออก
แต่ละหัวข้อของแบบจำลองด้านล่างเริ่มด้วยการระบุสองข้อเท็จจริงนี้ ตารางสำหรับตัดสินใจใกล้ท้ายบทความวางแบบจำลองทั้งสามเคียงกัน พร้อมความหมายของสัมประสิทธิ์ที่ยกกำลังแต่ละค่า
birth cohort จำลอง
cohort มีทารกแรกเกิด 500 คน โดย 200 คนเกิดจากมารดาที่ฝากครรภ์ไม่เพียงพอ และ 120 คนคลอดก่อนกำหนด (preterm birth) ซึ่งบางคนอยู่ในทั้งสองกลุ่ม ทารกแต่ละคนมีระดับ $Y$ ภายใน 28 วัน คือ 0 ไม่มี 1 เล็กน้อย 2 ปานกลาง 3 รุนแรง โดยรวมแล้วทารก 256 คนไม่มีภาวะแทรกซ้อน 101 คนมีภาวะเล็กน้อย 71 คนปานกลาง และ 72 คนรุนแรง
จำนวนเหล่านี้สร้างจากแบบจำลองที่ทราบค่า จึงตรวจสอบค่าประมาณทุกค่ากับแบบจำลองนั้นได้ สำหรับจุดตัดแต่ละจุดระหว่างระดับ แบบจำลองกำหนด log odds ของการอยู่ที่ระดับ $j$ หรือแย่กว่า:
$$\operatorname{logit} P(Y \ge j) = \alpha_j + \beta x + \gamma_j z$$ในที่นี้ $x$ เท่ากับ 1 เมื่อฝากครรภ์ไม่เพียงพอ $z$ เท่ากับ 1 เมื่อคลอดก่อนกำหนด และ $j = 1, 2, 3$ ระบุจุดตัดทั้งสามจุด intercept คือ $\alpha_1 = -0.4$, $\alpha_2 = -1.4$ และ $\alpha_3 = -2.5$ ความชันของ ANC คือ $\beta = \ln 2$ ซึ่งเป็น odds ratio 2.00 ที่ทุกจุดตัด ความชันของการคลอดก่อนกำหนดคือ $\gamma_1 = 0.3$, $\gamma_2 = 0.7$ และ $\gamma_3 = 1.2$ ซึ่งเป็น odds ratio 1.35, 2.01 และ 3.32 เพิ่มขึ้นตามความรุนแรง
ความจริงจำลองนี้ให้ ANC มี odds ratio ค่าเดียวที่ทุกจุดตัด และให้การคลอดก่อนกำหนดมี odds ratio ต่างกันที่แต่ละจุด จำนวนที่ใช้เป็นจำนวนคาดหวังภายใต้แบบจำลอง ปัดเป็นทารกจำนวนเต็ม ไม่ใช่การสุ่ม ค่า odds ratio ที่ประมาณได้ด้านล่างจึงใกล้ค่าที่ใช้สร้างข้อมูล
Stata: สร้างจำนวนของ cohort
* expected number at or above each cut, rounded to whole newborns
forvalues j = 1/3 {
generate int atleast`j' = floor(newborns*invlogit(`alpha`j'' + ln(2)*inadequate_anc + `gamma`j''*preterm) + 0.5)
}
generate int sev0 = newborns - atleast1
generate int sev1 = atleast1 - atleast2
generate int sev2 = atleast2 - atleast3
generate int sev3 = atleast3
list inadequate_anc preterm newborns sev0 sev1 sev2 sev3, noobs
. * expected number at or above each cut, rounded to whole newborns
. forvalues j = 1/3 {
2. generate int atleast`j' = floor(newborns*invlogit(`alpha`j'' + ln(2)*inadequate_anc + `gamma`j''*preterm) + 0.5)
3. }
. generate int sev0 = newborns - atleast1
. generate int sev1 = atleast1 - atleast2
. generate int sev2 = atleast2 - atleast3
. generate int sev3 = atleast3
. list inadequate_anc preterm newborns sev0 sev1 sev2 sev3, noobs
+-----------------------------------------------------------+
| inadeq~c preterm newborns sev0 sev1 sev2 sev3 |
|-----------------------------------------------------------|
| 0 0 240 144 49 29 18 |
| 0 1 60 31 9 7 13 |
| 1 0 140 60 34 26 20 |
| 1 1 60 21 9 9 21 |
+-----------------------------------------------------------+
แบบจำลองแรก: ผลลัพธ์สองค่าในทารกแรกเกิดที่เป็นอิสระต่อกัน
ชนิดของผลลัพธ์: สองค่า คือมีภาวะแทรกซ้อนใด ๆ (ระดับ 1 หรือแย่กว่า) เทียบกับไม่มี แบบแผนการสุ่มตัวอย่าง: ทารกแรกเกิดที่เป็นอิสระต่อกันจาก birth cohort เดียว แบบจำลอง: logistic regression ทั่วไป โดยมีการฝากครรภ์ไม่เพียงพอและการคลอดก่อนกำหนดเป็นตัวแปรร่วม การคลอดก่อนกำหนดอาจอยู่บนเส้นทางจาก ANC ไปสู่ภาวะแทรกซ้อนเอง ซึ่งจะกลับมากล่าวถึงใกล้ท้ายบทความ
จากทารก 500 คน มี 244 คนที่มีภาวะแทรกซ้อนระดับใดระดับหนึ่ง สัมประสิทธิ์ของ ANC ที่ประมาณได้เท่ากับ 0.691 บนมาตรวัด log odds เมื่อยกกำลังด้วยเลขชี้กำลังจะได้ odds ratio 2.00 (ช่วงเชื่อมั่น 95% 1.38 ถึง 2.88) ปรับด้วยการคลอดก่อนกำหนดแล้ว odds ratio ของการคลอดก่อนกำหนดคือ 1.40 (0.92 ถึง 2.13)
odds ratio ทั้งสองค่าอธิบายจุดตัดเพียงจุดเดียว คือไม่มีเทียบกับมีภาวะแทรกซ้อนใด ๆ แบบจำลองนี้แยกภาวะเล็กน้อยออกจากภาวะรุนแรงไม่ได้ การยุบสี่ระดับเหลือสองระดับทิ้งข้อมูลนั้นไป และแบบจำลองถัดไปเก็บข้อมูลนั้นไว้
Stata: logit และ logistic กับภาวะแทรกซ้อนใด ๆ
generate byte any_morbidity = severity >= 1
quietly count if any_morbidity == 1
canon bin.n_any 0 r(N)
logit any_morbidity i.inadequate_anc i.preterm, nolog
canon logit.b_anc 4 _b[1.inadequate_anc]
canon logit.se_anc 4 _se[1.inadequate_anc]
canon logit.b_preterm 4 _b[1.preterm]
canon logit.se_preterm 4 _se[1.preterm]
canon logit.b_cons 4 _b[_cons]
canon logit.ll 4 e(ll)
logistic any_morbidity i.inadequate_anc i.preterm, nolog
. generate byte any_morbidity = severity >= 1
. quietly count if any_morbidity == 1
. canon bin.n_any 0 r(N)
CANON c7.bin.n_any 244
. logit any_morbidity i.inadequate_anc i.preterm, nolog
Logistic regression Number of obs = 500
LR chi2(2) = 17.79
Prob > chi2 = 0.0001
Log likelihood = -337.53237 Pseudo R2 = 0.0257
----------------------------------------------------------------------------------
any_morbidity | Coefficient Std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | .6912976 .1868748 3.70 0.000 .3250298 1.057565
1.preterm | .3354302 .2149308 1.56 0.119 -.0858265 .756687
_cons | -.404775 .1254744 -3.23 0.001 -.6507003 -.1588497
----------------------------------------------------------------------------------
. canon logit.b_anc 4 _b[1.inadequate_anc]
CANON c7.logit.b_anc 0.6913
. canon logit.se_anc 4 _se[1.inadequate_anc]
CANON c7.logit.se_anc 0.1869
. canon logit.b_preterm 4 _b[1.preterm]
CANON c7.logit.b_preterm 0.3354
. canon logit.se_preterm 4 _se[1.preterm]
CANON c7.logit.se_preterm 0.2149
. canon logit.b_cons 4 _b[_cons]
CANON c7.logit.b_cons -0.4048
. canon logit.ll 4 e(ll)
CANON c7.logit.ll -337.5324
. logistic any_morbidity i.inadequate_anc i.preterm, nolog
Logistic regression Number of obs = 500
LR chi2(2) = 17.79
Prob > chi2 = 0.0001
Log likelihood = -337.53237 Pseudo R2 = 0.0257
----------------------------------------------------------------------------------
any_morbidity | Odds ratio Std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | 1.996304 .3730589 3.70 0.000 1.384072 2.879352
1.preterm | 1.398542 .3005898 1.56 0.119 .9177535 2.131204
_cons | .6671269 .0837074 -3.23 0.001 .5216803 .8531246
----------------------------------------------------------------------------------
Note: _cons estimates baseline odds.
R: glm() กับภาวะแทรกซ้อนใด ๆ
o$any_morbidity <- as.integer(o$severity >= 1)
canon("bin.n_any", sum(o$any_morbidity), 0)
fl <- glm(any_morbidity ~ inadequate_anc + preterm, family = binomial, data = o)
print(summary(fl)$coefficients)
se_l <- sqrt(diag(vcov(fl)))
canon("logit.b_anc", coef(fl)[["inadequate_anc"]])
canon("logit.se_anc", se_l[["inadequate_anc"]])
canon("logit.b_preterm", coef(fl)[["preterm"]])
canon("logit.se_preterm", se_l[["preterm"]])
canon("logit.b_cons", coef(fl)[["(Intercept)"]])
canon("logit.ll", as.numeric(logLik(fl)))
# the same fit shown as odds ratios
print(exp(cbind(OR = coef(fl), confint.default(fl))))
> o$any_morbidity <- as.integer(o$severity >= 1)
> canon("bin.n_any", sum(o$any_morbidity), 0)
CANON c7.bin.n_any 244
> fl <- glm(any_morbidity ~ inadequate_anc + preterm,
+ family = binomial, data = o)
> print(summary(fl)$coefficients)
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.4047750 0.1254744 -3.225956 0.0012555244
inadequate_anc 0.6912976 0.1868748 3.699256 0.0002162321
preterm 0.3354302 0.2149308 1.560643 0.1186080947
> se_l <- sqrt(diag(vcov(fl)))
> canon("logit.b_anc", coef(fl)[["inadequate_anc"]])
CANON c7.logit.b_anc 0.6913
> canon("logit.se_anc", se_l[["inadequate_anc"]])
CANON c7.logit.se_anc 0.1869
> canon("logit.b_preterm", coef(fl)[["preterm"]])
CANON c7.logit.b_preterm 0.3354
> canon("logit.se_preterm", se_l[["preterm"]])
CANON c7.logit.se_preterm 0.2149
> canon("logit.b_cons", coef(fl)[["(Intercept)"]])
CANON c7.logit.b_cons -0.4048
> canon("logit.ll", as.numeric(logLik(fl)))
CANON c7.logit.ll -337.5324
> print(exp(cbind(OR = coef(fl), confint.default(fl))))
OR 2.5 % 97.5 %
(Intercept) 0.6671269 0.5216803 0.8531246
inadequate_anc 1.9963043 1.3840720 2.8793525
preterm 1.3985420 0.9177535 2.1312037
แบบจำลองที่สอง: ระดับที่มีลำดับในทารกกลุ่มเดียวกัน
ชนิดของผลลัพธ์: หมวดที่มีลำดับ คือระดับภาวะแทรกซ้อนสี่ระดับ แบบแผนการสุ่มตัวอย่าง: ทารกแรกเกิดอิสระต่อกันกลุ่มเดิม แบบจำลอง: การถดถอยโลจิสติกสำหรับผลลัพธ์แบบลำดับ (ordinal logistic regression) ซึ่งใช้ทุกระดับพร้อมกันและรักษาลำดับของมันไว้ [2]
ผลลัพธ์สี่ระดับมี จุดตัด (cut point) สามจุด คือเส้นแบ่งระหว่างระดับที่อยู่ติดกัน จุดตัดแต่ละจุดแบ่งทารกเป็นสองกลุ่ม คือระดับ $j$ หรือแย่กว่า กับดีกว่าระดับ $j$ log odds ของการอยู่ที่ระดับ $j$ หรือแย่กว่าคือ cumulative logit ของจุดตัดนั้น
ordinal logistic regression ประมาณ cumulative logit หนึ่งค่าต่อจุดตัด แต่ละค่ามี intercept ของตัวเอง ทว่าใช้ชุดความชันร่วมกันชุดเดียว:
$$\operatorname{logit} P(Y \ge j) = \alpha_j + \beta x + \gamma z$$intercept $\alpha_j$ ลดลงเมื่อ $j$ เพิ่ม เพราะทารกที่ถึงระดับสูงขึ้นมีน้อยลง ความชัน $\beta$ และ $\gamma$ ไม่มี $j$ กำกับ odds ratio ของ ANC คือ $e^{\beta}$ จึงเท่ากันทุกจุดตัด นี่คือ ข้อสมมติ proportional odds (proportional odds assumption) ต่างจากแบบจำลองที่ใช้สร้างข้อมูล แบบจำลองนี้บังคับให้การคลอดก่อนกำหนดมีความชันค่าเดียวที่ทุกจุดตัด
คำสั่ง ologit ของ Stata และ polr() ของ R พิมพ์ intercept ด้วยเครื่องหมายตรงข้าม เป็นจุดตัด $\kappa_j = -\alpha_j$ ที่เพิ่มขึ้นตาม $j$ ความชันคงเครื่องหมายเดิม สัมประสิทธิ์บวกทำให้ระดับสูงมีโอกาสมากขึ้น
ใน cohort จำลอง odds ratio ของ ANC ที่ปรับด้วยการคลอดก่อนกำหนดคือ 2.01 (ช่วงเชื่อมั่น 95% 1.43 ถึง 2.82) ใกล้เคียง 2.00 ที่ใช้สร้างข้อมูล เป็น odds ratio ของการอยู่ในระดับที่สูงขึ้น ไม่ว่าจะลากเส้นแบ่งระหว่างระดับตรงไหน ช่วงของมันแคบกว่า 1.38 ถึง 2.88 ของแบบจำลองสองค่า เพราะแบบจำลองแบบลำดับใช้ระดับที่อยู่เหนือจุดตัดแรกด้วย
อ่านให้เป็นรูปธรรม odds ratio นี้ใช้กับทุกวิธีแบ่งระดับพร้อมกัน ในทารกที่มีสถานะคลอดก่อนกำหนดเหมือนกัน การฝากครรภ์ไม่เพียงพอเพิ่ม odds ของภาวะแทรกซ้อนใด ๆ ประมาณสองเท่า ของภาวะปานกลางหรือรุนแรง และของภาวะรุนแรง เช่นเดียวกัน การยึดสถานะการคลอดก่อนกำหนดไว้ให้เท่ากันเช่นนี้อาจละส่วนหนึ่งของผลของ ANC ไป ดังที่หัวข้อ สิ่งที่ odds ratio ทั้งสามมีร่วมกัน อธิบาย
เมื่อผลลัพธ์พบบ่อย การเพิ่ม odds เป็นสองเท่าไม่ได้เพิ่มความเสี่ยงเป็นสองเท่า ดังที่ความน่าจะเป็นด้านล่างแสดง จุดตัดคือ 0.458, 1.361 และ 2.267 พิจารณาทารกที่คลอดครบกำหนด ไม่ได้คลอดก่อนกำหนด และฝากครรภ์เพียงพอ แบบจำลองให้ทารกคนนี้มีความน่าจะเป็น 0.613 ที่จะไม่มีภาวะแทรกซ้อน และ 0.094 ที่จะมีภาวะรุนแรง
เมื่อฝากครรภ์ไม่เพียงพอ ค่าเหล่านี้กลายเป็น 0.441 และ 0.172 odds ratio ค่าเดียวจึงเลื่อนการกระจายของระดับทั้งหมดไปทางด้านรุนแรง การเลื่อนเพียงครั้งเดียวนี้คือสิ่งที่ข้อสมมติ proportional odds ให้ และเป็นสิ่งที่มันเสี่ยง
odds ratio ของการคลอดก่อนกำหนดคือ 1.84 (1.24 ถึง 2.74) ต่างจาก odds ratio ของ ANC ตรงที่มันสรุป odds ratio ซึ่งในความจริงจำลองเพิ่มจาก 1.35 ที่จุดตัดแรกเป็น 3.32 ที่จุดตัดที่สาม คำถามว่าตัวเลขเดียวสรุปได้เหมาะสมหรือไม่คือสิ่งที่การตรวจ proportional odds ถาม
Stata: ologit กับระดับทั้งสี่
ologit severity i.inadequate_anc i.preterm, nolog
. ologit severity i.inadequate_anc i.preterm, nolog
Ordered logistic regression Number of obs = 500
LR chi2(2) = 27.91
Prob > chi2 = 0.0000
Log likelihood = -597.08706 Pseudo R2 = 0.0228
----------------------------------------------------------------------------------
severity | Coefficient Std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | .6966207 .1732433 4.02 0.000 .3570702 1.036171
1.preterm | .6124173 .2012158 3.04 0.002 .2180416 1.006793
-----------------+----------------------------------------------------------------
/cut1 | .4581421 .1220414 .2189454 .6973388
/cut2 | 1.360579 .1350555 1.095875 1.625283
/cut3 | 2.266819 .1627571 1.947821 2.585817
----------------------------------------------------------------------------------
R: polr() กับระดับทั้งสี่
fo <- polr(factor(severity, levels = 0:3, ordered = TRUE) ~ inadequate_anc + preterm, data = o,
method = "logistic", Hess = TRUE)
print(summary(fo))
> fo <- polr(factor(severity, levels = 0:3, ordered = TRUE) ~
+ inadequate_anc + preterm, data = o, method = "logistic",
+ Hess = TRUE)
> print(summary(fo))
Call:
polr(formula = factor(severity, levels = 0:3, ordered = TRUE) ~
inadequate_anc + preterm, data = o, Hess = TRUE, method = "logistic")
Coefficients:
Value Std. Error t value
inadequate_anc 0.6966 0.1732 4.021
preterm 0.6124 0.2012 3.044
Intercepts:
Value Std. Error t value
0|1 0.4581 0.1220 3.7540
1|2 1.3606 0.1351 10.0742
2|3 2.2668 0.1628 13.9276
Residual Deviance: 1194.174
AIC: 1204.174
ตรวจ proportional odds ทีละจุดตัด
ตรวจข้อสมมติ proportional odds ได้โดยตรง ประมาณ logistic regression ทั่วไปที่จุดตัดแต่ละจุดทั้งสามจุด สามแบบจำลองเปรียบเทียบระดับ 1 หรือแย่กว่ากับไม่มี ระดับ 2 หรือแย่กว่ากับต่ำกว่านั้น และระดับ 3 กับต่ำกว่านั้น ถ้าข้อสมมติเป็นจริง odds ratio ของตัวแปรร่วมแต่ละตัวควรใกล้เคียงกันทั้งสามจุด
ตารางด้านล่างแสดงว่า ANC เป็นไปตามที่ข้อสมมติต้องการ โดยมี odds ratio 2.00, 2.01 และ 2.01 ที่สามจุดตัด การคลอดก่อนกำหนดไม่เป็นเช่นนั้น คือ 1.40, 2.05 และ 3.31 ใกล้ความจริงจำลองที่ 1.35, 2.01 และ 3.32 แบบจำลองที่จุดตัดแรกคือแบบจำลองสองค่าข้างต้น ตัวเลขจึงตรงกับแบบจำลองนั้น
| จุดตัด | odds ratio ของ ANC ค่าที่ประมาณได้ | odds ratio ของ ANC ความจริงจำลอง | odds ratio ของการคลอดก่อนกำหนด ค่าที่ประมาณได้ | odds ratio ของการคลอดก่อนกำหนด ความจริงจำลอง |
|---|---|---|---|---|
| ระดับ 1 หรือแย่กว่า เทียบกับไม่มี | 2.00 | 2.00 | 1.40 | 1.35 |
| ระดับ 2 หรือแย่กว่า เทียบกับระดับ 0 ถึง 1 | 2.01 | 2.00 | 2.05 | 2.01 |
| ระดับ 3 เทียบกับระดับ 0 ถึง 2 | 2.01 | 2.00 | 3.31 | 3.32 |
| แบบจำลองแบบลำดับ odds ratio ค่าเดียวสำหรับทุกจุดตัด | 2.01 | 2.00 | 1.84 | ไม่มีค่าเดียว |
Stata: logistic หนึ่งแบบจำลองต่อจุดตัด
* the proportional-odds check, step 1: a separate binary logit at each cut (severity >= 1, >= 2, >= 3);
* plain 0/1 covariates, so that e(b) holds inadequate_anc, preterm and the constant in that order
forvalues j = 1/3 {
generate byte ge`j' = severity >= `j'
logit ge`j' inadequate_anc preterm, nolog
matrix bcut`j' = e(b)
predict double pi`j', pr
canon cut`j'.b_anc 4 _b[inadequate_anc]
canon cut`j'.or_anc 4 exp(_b[inadequate_anc])
canon cut`j'.b_preterm 4 _b[preterm]
canon cut`j'.or_preterm 4 exp(_b[preterm])
}
. * the proportional-odds check, step 1: a separate binary logit at each cut (severity >= 1, >= 2, >= 3);
. * plain 0/1 covariates, so that e(b) holds inadequate_anc, preterm and the constant in that order
. forvalues j = 1/3 {
2. generate byte ge`j' = severity >= `j'
3. logit ge`j' inadequate_anc preterm, nolog
4. matrix bcut`j' = e(b)
5. predict double pi`j', pr
6. canon cut`j'.b_anc 4 _b[inadequate_anc]
7. canon cut`j'.or_anc 4 exp(_b[inadequate_anc])
8. canon cut`j'.b_preterm 4 _b[preterm]
9. canon cut`j'.or_preterm 4 exp(_b[preterm])
10. }
Logistic regression Number of obs = 500
LR chi2(2) = 17.79
Prob > chi2 = 0.0001
Log likelihood = -337.53237 Pseudo R2 = 0.0257
--------------------------------------------------------------------------------
ge1 | Coefficient Std. err. z P>|z| [95% conf. interval]
---------------+----------------------------------------------------------------
inadequate_anc | .6912976 .1868748 3.70 0.000 .3250298 1.057565
preterm | .3354302 .2149308 1.56 0.119 -.0858265 .756687
_cons | -.404775 .1254744 -3.23 0.001 -.6507003 -.1588497
--------------------------------------------------------------------------------
CANON c7.cut1.b_anc 0.6913
CANON c7.cut1.or_anc 1.9963
CANON c7.cut1.b_preterm 0.3354
CANON c7.cut1.or_preterm 1.3985
Logistic regression Number of obs = 500
LR chi2(2) = 24.29
Prob > chi2 = 0.0000
Log likelihood = -287.11966 Pseudo R2 = 0.0406
--------------------------------------------------------------------------------
ge2 | Coefficient Std. err. z P>|z| [95% conf. interval]
---------------+----------------------------------------------------------------
inadequate_anc | .6964978 .2038964 3.42 0.001 .2968681 1.096127
preterm | .7169957 .2238947 3.20 0.001 .2781702 1.155821
_cons | -1.411917 .1515253 -9.32 0.000 -1.708901 -1.114933
--------------------------------------------------------------------------------
CANON c7.cut2.b_anc 0.6965
CANON c7.cut2.or_anc 2.0067
CANON c7.cut2.b_preterm 0.7170
CANON c7.cut2.or_preterm 2.0483
Logistic regression Number of obs = 500
LR chi2(2) = 29.04
Prob > chi2 = 0.0000
Log likelihood = -191.55985 Pseudo R2 = 0.0705
--------------------------------------------------------------------------------
ge3 | Coefficient Std. err. z P>|z| [95% conf. interval]
---------------+----------------------------------------------------------------
inadequate_anc | .6983656 .2648017 2.64 0.008 .1793638 1.217367
preterm | 1.1974 .2683201 4.46 0.000 .6715027 1.723298
_cons | -2.501102 .2175031 -11.50 0.000 -2.9274 -2.074803
--------------------------------------------------------------------------------
CANON c7.cut3.b_anc 0.6984
CANON c7.cut3.or_anc 2.0105
CANON c7.cut3.b_preterm 1.1974
CANON c7.cut3.or_preterm 3.3115
Brant test: ความชันเท่ากันทุกจุดตัดหรือไม่
การเปรียบเทียบ odds ratio สามค่าด้วยตาไม่ใช่การทดสอบ เพราะค่าที่ประมาณแต่ละค่ามีความคลาดเคลื่อนจากการสุ่มตัวอย่าง Brant test ถามว่าความชันของแบบจำลองแยกตามจุดตัดต่างกันเกินกว่าความคลาดเคลื่อนนั้นจะอธิบายได้หรือไม่ [3] เป็น การทดสอบ Wald (Wald test) คือชั่งผลต่างของความชันเทียบกับความแปรปรวนจากการสุ่มตัวอย่างที่ประมาณได้ และให้สถิติไคสแควร์
แบบจำลองแยกใช้ทารกกลุ่มเดียวกัน ความชันของมันจึงสัมพันธ์กัน วิธีของ Brant ประมาณความสัมพันธ์นั้นจากแบบจำลองทั้งสามจุดตัด ขณะที่การเปรียบเทียบแบบตรงไปตรงมาจะมองข้ามมัน [3] ทดสอบได้ทั้งตัวแปรร่วมทั้งหมดพร้อมกันหรือทีละตัว องศาอิสระ (degrees of freedom) นับผลต่างของความชันที่ทดสอบ ในที่นี้คือสองค่าต่อตัวแปรร่วมหนึ่งตัว
สำหรับตัวแปรร่วมทั้งหมดรวมกัน สถิติเท่ากับ 11.72 ที่องศาอิสระ 4 (P = 0.020) การทดสอบพบหลักฐานต่อข้อสมมติสำหรับแบบจำลองโดยรวม สำหรับการคลอดก่อนกำหนดเพียงอย่างเดียว สถิติเท่ากับ 11.65 ที่องศาอิสระ 2 (P = 0.003)
สำหรับ ANC เพียงอย่างเดียว สถิติเท่ากับ 0.001 ที่องศาอิสระ 2 (P = 0.999) เพราะจำนวนถูกสร้างโดยไม่มีสัญญาณรบกวนจากการสุ่มตัวอย่าง สถิติของ ANC จึงอยู่เกือบที่ศูนย์ ในตัวอย่างจริงที่ข้อสมมติเป็นจริง สถิติแปรผันไปตามแต่ละตัวอย่าง และค่า P อาจตกที่ใดก็ได้ระหว่าง 0 ถึง 1
การทดสอบรวมบอกว่ามีบางอย่างผิดปกติ การทดสอบรายตัวแปรร่วมบอกว่าอยู่ที่ใด มีคำสั่งเสริมที่ทำการทดสอบนี้อัตโนมัติ คือ brant ใน Stata จากแพ็กเกจ SPost ที่ผู้ใช้เขียน และแพ็กเกจ brant ใน R โค้ดด้านล่างไม่ใช้ทั้งสอง แต่คำนวณการทดสอบโดยตรงจากแบบจำลองสามจุดตัด จึงเห็นทุกขั้นตอน
Stata: Brant test ที่คำนวณใน Mata
* step 2, the Brant (1990) Wald test that the slopes are equal across the three cuts. For cuts j and l the
* covariance of the two slope vectors is A_j^-1 [X' diag(pi_hi (1 - pi_lo)) X] A_l^-1, with
* A_j = X' diag(pi_j (1 - pi_j)) X, lo = min(j, l) and hi = max(j, l); the contrasts compare cut 1 with
* cuts 2 and 3 (4 df for both covariates, 2 df for each covariate)
mata:
X = st_data(., ("inadequate_anc", "preterm")), J(st_nobs(), 1, 1)
P = st_data(., ("pi1", "pi2", "pi3"))
B1 = st_matrix("bcut1")
B2 = st_matrix("bcut2")
B3 = st_matrix("bcut3")
b = (B1[|1, 1 \ 1, 2|], B2[|1, 1 \ 1, 2|], B3[|1, 1 \ 1, 2|])'
V = J(6, 6, 0)
for (j = 1; j <= 3; j++) {
for (l = 1; l <= 3; l++) {
lo = min((j, l))
hi = max((j, l))
Ajl = invsym(cross(X, P[., j]:*(1 :- P[., j]), X)) * cross(X, P[., hi]:*(1 :- P[., lo]), X) * invsym(cross(X, P[., l]:*(1 :- P[., l]), X))
V[|2*j - 1, 2*l - 1 \ 2*j, 2*l|] = Ajl[|1, 1 \ 2, 2|]
}
}
D = (1, 0, -1, 0, 0, 0 \ 0, 1, 0, -1, 0, 0 \ 1, 0, 0, 0, -1, 0 \ 0, 1, 0, 0, 0, -1)
Da = D[(1 \ 3), .]
Dp = D[(2 \ 4), .]
st_numscalar("sc_brant_all", (D*b)' * invsym(D*V*D') * (D*b))
st_numscalar("sc_brant_anc", (Da*b)' * invsym(Da*V*Da') * (Da*b))
st_numscalar("sc_brant_pre", (Dp*b)' * invsym(Dp*V*Dp') * (Dp*b))
end
canon brant.chi2_all 4 scalar(sc_brant_all)
canon brant.df_all 0 4
canon brant.p_all 4 chi2tail(4, scalar(sc_brant_all))
canon brant.chi2_anc 4 scalar(sc_brant_anc)
canon brant.p_anc 4 chi2tail(2, scalar(sc_brant_anc))
canon brant.chi2_preterm 4 scalar(sc_brant_pre)
canon brant.p_preterm 4 chi2tail(2, scalar(sc_brant_pre))
canon brant.df_each 0 2
. * step 2, the Brant (1990) Wald test that the slopes are equal across the three cuts. For cuts j and l the
. * covariance of the two slope vectors is A_j^-1 [X' diag(pi_hi (1 - pi_lo)) X] A_l^-1, with
. * A_j = X' diag(pi_j (1 - pi_j)) X, lo = min(j, l) and hi = max(j, l); the contrasts compare cut 1 with
. * cuts 2 and 3 (4 df for both covariates, 2 df for each covariate)
. mata:
------------------------------------------------- mata (type end to exit) -------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
: X = st_data(., ("inadequate_anc", "preterm")), J(st_nobs(), 1, 1)
: P = st_data(., ("pi1", "pi2", "pi3"))
: B1 = st_matrix("bcut1")
: B2 = st_matrix("bcut2")
: B3 = st_matrix("bcut3")
: b = (B1[|1, 1 \ 1, 2|], B2[|1, 1 \ 1, 2|], B3[|1, 1 \ 1, 2|])'
: V = J(6, 6, 0)
: for (j = 1; j <= 3; j++) {
> for (l = 1; l <= 3; l++) {
> lo = min((j, l))
> hi = max((j, l))
> Ajl = invsym(cross(X, P[., j]:*(1 :- P[., j]), X)) * cross(X, P[., hi]:*(1 :- P[., lo]), X) * invsym(cross(X, P[., l]:*(1 :- P[., l]), X))
> V[|2*j - 1, 2*l - 1 \ 2*j, 2*l|] = Ajl[|1, 1 \ 2, 2|]
> }
> }
: D = (1, 0, -1, 0, 0, 0 \ 0, 1, 0, -1, 0, 0 \ 1, 0, 0, 0, -1, 0 \ 0, 1, 0, 0, 0, -1)
: Da = D[(1 \ 3), .]
: Dp = D[(2 \ 4), .]
: st_numscalar("sc_brant_all", (D*b)' * invsym(D*V*D') * (D*b))
: st_numscalar("sc_brant_anc", (Da*b)' * invsym(Da*V*Da') * (Da*b))
: st_numscalar("sc_brant_pre", (Dp*b)' * invsym(Dp*V*Dp') * (Dp*b))
: end
---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
. canon brant.chi2_all 4 scalar(sc_brant_all)
CANON c7.brant.chi2_all 11.7179
. canon brant.df_all 0 4
CANON c7.brant.df_all 4
. canon brant.p_all 4 chi2tail(4, scalar(sc_brant_all))
CANON c7.brant.p_all 0.0196
. canon brant.chi2_anc 4 scalar(sc_brant_anc)
CANON c7.brant.chi2_anc 0.0011
. canon brant.p_anc 4 chi2tail(2, scalar(sc_brant_anc))
CANON c7.brant.p_anc 0.9994
. canon brant.chi2_preterm 4 scalar(sc_brant_pre)
CANON c7.brant.chi2_preterm 11.6476
. canon brant.p_preterm 4 chi2tail(2, scalar(sc_brant_pre))
CANON c7.brant.p_preterm 0.0030
. canon brant.df_each 0 2
CANON c7.brant.df_each 2
R: แบบจำลองแยกตามจุดตัดและ Brant test
# the proportional-odds check, step 1: a separate binary logit at each cut (severity >= 1, >= 2, >= 3)
fits <- lapply(1:3, function(j) {
o$ge <- as.integer(o$severity >= j)
glm(ge ~ inadequate_anc + preterm, family = binomial, data = o)
})
for (j in 1:3) {
f <- fits[[j]]
canon(sprintf("cut%d.b_anc", j), coef(f)[["inadequate_anc"]])
canon(sprintf("cut%d.or_anc", j), exp(coef(f)[["inadequate_anc"]]))
canon(sprintf("cut%d.b_preterm", j), coef(f)[["preterm"]])
canon(sprintf("cut%d.or_preterm", j), exp(coef(f)[["preterm"]]))
}
# step 2, the Brant (1990) Wald test that the slopes are equal across the three cuts. For cuts j and l the
# covariance of the two slope vectors is A_j^-1 [X' diag(pi_hi (1 - pi_lo)) X] A_l^-1, with
# A_j = X' diag(pi_j (1 - pi_j)) X, lo = min(j, l) and hi = max(j, l); the contrasts compare cut 1 with
# cuts 2 and 3 (4 df for both covariates, 2 df for each covariate)
X <- cbind(inadequate_anc = o$inadequate_anc, preterm = o$preterm, cons = 1)
P <- sapply(fits, fitted)
b <- unlist(lapply(fits, function(f) coef(f)[c("inadequate_anc", "preterm")]))
V <- matrix(0, 6, 6)
for (j in 1:3) {
for (l in 1:3) {
lo <- min(j, l)
hi <- max(j, l)
Ajl <- solve(crossprod(X, X * (P[, j] * (1 - P[, j])))) %*%
crossprod(X, X * (P[, hi] * (1 - P[, lo]))) %*%
solve(crossprod(X, X * (P[, l] * (1 - P[, l]))))
V[(2 * j - 1):(2 * j), (2 * l - 1):(2 * l)] <- Ajl[1:2, 1:2]
}
}
D <- rbind(c(1, 0, -1, 0, 0, 0), c(0, 1, 0, -1, 0, 0), c(1, 0, 0, 0, -1, 0), c(0, 1, 0, 0, 0, -1))
wald <- function(Dm) {
db <- Dm %*% b
as.numeric(t(db) %*% solve(Dm %*% V %*% t(Dm)) %*% db)
}
chi_all <- wald(D)
chi_anc <- wald(D[c(1, 3), ])
chi_pre <- wald(D[c(2, 4), ])
canon("brant.chi2_all", chi_all)
canon("brant.df_all", 4, 0)
canon("brant.p_all", pchisq(chi_all, 4, lower.tail = FALSE))
canon("brant.chi2_anc", chi_anc)
canon("brant.p_anc", pchisq(chi_anc, 2, lower.tail = FALSE))
canon("brant.chi2_preterm", chi_pre)
canon("brant.p_preterm", pchisq(chi_pre, 2, lower.tail = FALSE))
canon("brant.df_each", 2, 0)
> fits <- lapply(1:3, function(j) {
+ o$ge <- as.integer(o$severity >= j)
+ glm(ge ~ inadequate_anc + preterm, family = binomial, data = o)
+ })
> for (j in 1:3) {
+ f <- fits[[j]]
+ canon(sprintf("cut%d.b_anc", j), coef(f)[["inadequate_anc"]])
+ canon(sprintf("cut%d.or_anc", j), exp(coef(f)[["inadequate_anc"]]))
+ canon(sprintf("cut%d.b_preterm", j), coef(f)[["preterm"]])
+ canon(sprintf("cut%d.or_preterm", j), exp(coef(f)[["preterm"]]))
+ }
CANON c7.cut1.b_anc 0.6913
CANON c7.cut1.or_anc 1.9963
CANON c7.cut1.b_preterm 0.3354
CANON c7.cut1.or_preterm 1.3985
CANON c7.cut2.b_anc 0.6965
CANON c7.cut2.or_anc 2.0067
CANON c7.cut2.b_preterm 0.7170
CANON c7.cut2.or_preterm 2.0483
CANON c7.cut3.b_anc 0.6984
CANON c7.cut3.or_anc 2.0105
CANON c7.cut3.b_preterm 1.1974
CANON c7.cut3.or_preterm 3.3115
> X <- cbind(inadequate_anc = o$inadequate_anc, preterm = o$preterm,
+ cons = 1)
> P <- sapply(fits, fitted)
> b <- unlist(lapply(fits, function(f) coef(f)[c("inadequate_anc",
+ "preterm")]))
> V <- matrix(0, 6, 6)
> for (j in 1:3) {
+ for (l in 1:3) {
+ lo <- min(j, l)
+ hi <- max(j, l)
+ Ajl <- solve(crossprod(X, X * (P[, j] * (1 - P[, j])))) %*%
+ crossprod(X, X * (P[, hi] * (1 - P[, lo]))) %*% solve(crossprod(X,
+ X * (P[, l] * (1 - P[, l]))))
+ V[(2 * j - 1):(2 * j), (2 * l - 1):(2 * l)] <- Ajl[1:2,
+ 1:2]
+ }
+ }
> D <- rbind(c(1, 0, -1, 0, 0, 0), c(0, 1, 0, -1, 0,
+ 0), c(1, 0, 0, 0, -1, 0), c(0, 1, 0, 0, 0, -1))
> wald <- function(Dm) {
+ db <- Dm %*% b
+ as.numeric(t(db) %*% solve(Dm %*% V %*% t(Dm)) %*% db)
+ }
> chi_all <- wald(D)
> chi_anc <- wald(D[c(1, 3), ])
> chi_pre <- wald(D[c(2, 4), ])
> canon("brant.chi2_all", chi_all)
CANON c7.brant.chi2_all 11.7179
> canon("brant.df_all", 4, 0)
CANON c7.brant.df_all 4
> canon("brant.p_all", pchisq(chi_all, 4, lower.tail = FALSE))
CANON c7.brant.p_all 0.0196
> canon("brant.chi2_anc", chi_anc)
CANON c7.brant.chi2_anc 0.0011
> canon("brant.p_anc", pchisq(chi_anc, 2, lower.tail = FALSE))
CANON c7.brant.p_anc 0.9994
> canon("brant.chi2_preterm", chi_pre)
CANON c7.brant.chi2_preterm 11.6476
> canon("brant.p_preterm", pchisq(chi_pre, 2, lower.tail = FALSE))
CANON c7.brant.p_preterm 0.0030
> canon("brant.df_each", 2, 0)
CANON c7.brant.df_each 2
เมื่อ proportional odds ไม่เป็นจริง
ตัวแปรร่วมที่ถูกชี้ไม่ได้ทำให้แบบจำลองแบบลำดับใช้ไม่ได้ ในที่นี้การตรวจชี้เฉพาะการคลอดก่อนกำหนด odds ratio ของ ANC ที่ 2.01 ยังสรุป odds ratio ของจุดตัดสามค่าที่เกือบเท่ากัน ส่วน odds ratio ของการคลอดก่อนกำหนดที่ 1.84 ไม่เป็นเช่นนั้น เพราะอยู่ระหว่าง odds ratio ของจุดตัดที่ไล่จาก 1.40 ถึง 3.31
มีทางออกสามแบบที่ผ่อนข้อสมมตินี้ แต่ละแบบถูกระบุชื่อไว้ที่นี่และอธิบายในบทความที่ลิงก์ไว้ด้านล่าง partial proportional odds model ให้ตัวแปรร่วมที่ถูกชี้ ในที่นี้คือการคลอดก่อนกำหนด มี odds ratio ของตัวเองที่แต่ละจุดตัด generalized ordered logit model ปล่อยตัวแปรร่วมทุกตัวในแบบเดียวกัน multinomial logistic regression ทิ้งลำดับและเปรียบเทียบแต่ละระดับกับระดับอ้างอิง
ใน cohort นี้ รายงานที่สมเหตุสมผลแบบหนึ่งคือคง odds ratio แบบลำดับของ ANC ไว้ และให้การคลอดก่อนกำหนดมี odds ratio เฉพาะจุดตัดสามค่า partial proportional odds model ให้ได้ทั้งสองในการประมาณครั้งเดียว ไม่ว่าทางใด รายงานควรระบุว่าได้ทำการตรวจและพบอะไร
บทความแยกเรื่อง ordinal logistic regression และเรื่อง multinomial เทียบกับ ordinal regression อธิบายทางออกเหล่านี้ ไม่ว่าเลือกแบบใด การแสดง odds ratio เฉพาะจุดตัดควบคู่กับค่าสรุปอาจช่วยผู้อ่าน ค่า P ของ Brant ยังขึ้นกับขนาดตัวอย่างด้วย การศึกษาขนาดใหญ่อาจชี้ความเบี่ยงเบนเล็กน้อย และการศึกษาขนาดเล็กอาจมองข้ามความเบี่ยงเบนขนาดใหญ่
แบบจำลองที่สาม: คู่ case-control ที่จับคู่
ชนิดของผลลัพธ์: สองค่าอีกครั้ง คือ case หรือ control แบบแผนการสุ่มตัวอย่าง: ชุดที่จับคู่ แต่ละชุดคือ case หนึ่งคนกับ control ที่ถูกสุ่มมาให้ แบบจำลอง: conditional logistic regression ซึ่งเปรียบเทียบ case แต่ละคนกับ control ของมันเท่านั้น [4]
case ทั้ง 20 คนมีภาวะแทรกซ้อนภายใน 28 วันหลังเกิด control ของแต่ละคนคือทารกที่ยังไม่มีภาวะแทรกซ้อนในวันที่ case มี คู่จึงถูกจับคู่ตามอายุเป็นวัน และ control อาจกลายเป็น case ในภายหลังได้ การเลือก control จากผู้ที่ยังอยู่ในกลุ่มเสี่ยงเรียกว่า การสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยง (risk-set sampling) อธิบายไว้ในบทความแยกเรื่องการเลือกกลุ่มควบคุม
case กับ control ของมันก่อให้เกิด ชุดที่จับคู่ (matched set) ในที่นี้คือหนึ่งคู่ ซึ่งมีรหัสระบุของตัวเองในข้อมูล คู่ 20 คู่ให้ 40 แถว หนึ่งแถวต่อสมาชิกหนึ่งคน สิ่งสัมผัสคือการฝากครรภ์ไม่เพียงพอ ซึ่งบันทึกไว้สำหรับสมาชิกทั้งสองคนของทุกคู่
ประเภทของการจับคู่อธิบายไว้ในบทความเรื่อง conditional logistic regression เช่นเดียวกับ overmatching คือการจับคู่ตามปัจจัยที่สัมพันธ์กับสิ่งสัมผัสแต่ไม่ใช่ตัวกวน ซึ่งอาจเสียความแม่นยำ
| ANC ของ case | ANC ของ control | จำนวนคู่ | ชนิดของคู่ |
|---|---|---|---|
| ไม่เพียงพอ | ไม่เพียงพอ | 8 | สอดคล้องกัน |
| ไม่เพียงพอ | เพียงพอ | 4 | ไม่สอดคล้องกัน |
| เพียงพอ | ไม่เพียงพอ | 1 | ไม่สอดคล้องกัน |
| เพียงพอ | เพียงพอ | 7 | สอดคล้องกัน |
| ทุกคู่ | ทุกแบบ | 20 | สอดคล้องกัน 15 คู่ ไม่สอดคล้องกัน 5 คู่ |
คู่ที่สอดคล้องกันและคู่ที่ไม่สอดคล้องกัน
คู่ที่สอดคล้องกัน (concordant pair) คือคู่ที่ case และ control มีสถานะสิ่งสัมผัสเหมือนกัน มี 15 คู่ โดย 8 คู่สัมผัสทั้งสองคน และ 7 คู่ไม่สัมผัสทั้งสองคน คู่ที่ไม่สอดคล้องกัน (discordant pair) คือคู่ที่ทั้งสองต่างกัน มี 4 คู่ที่เฉพาะ case สัมผัส และ 1 คู่ที่เฉพาะ control สัมผัส การแบ่งนี้เป็นตัวขับเคลื่อนทุกอย่างที่ตามมา
| ANC | Case | Control |
|---|---|---|
| ไม่เพียงพอ | 12 | 9 |
| เพียงพอ | 8 | 11 |
| รวม | 20 | 20 |
ตัวอย่างคำนวณด้วยมือ: รวมคู่เข้าด้วยกัน
ตัวอย่างคำนวณด้วยมือ ชนิดของผลลัพธ์: สองค่า แบบแผนการสุ่มตัวอย่าง: คู่ที่จับคู่ ซึ่งแบบจำลองนี้ปฏิบัติต่อเสมือนทารกอิสระต่อกัน 40 คน
logistic regression ทั่วไปไม่สนใจว่า control คนไหนอยู่กับ case คนไหน มองเห็น case 20 คนและ control 20 คน และรวมเป็นตารางข้างบน คือสิ่งสัมผัสจำแนกตามสถานะ case
-
Case ที่สัมผัส
\[ 8 + 4 = 12 \]
คู่ที่สัมผัสทั้งสองคน บวกคู่ที่เฉพาะ case สัมผัส case อีก 7 + 1 = 8 คนฝากครรภ์เพียงพอ
-
Control ที่สัมผัส
\[ 8 + 1 = 9 \]
คู่ที่สัมผัสทั้งสองคน บวกคู่ที่เฉพาะ control สัมผัส control อีก 7 + 4 = 11 คนฝากครรภ์เพียงพอ
-
odds ของการสัมผัสใน case
\[ \frac{12}{8} = 1.50 \]
ฝากครรภ์ไม่เพียงพอเทียบกับเพียงพอใน case 20 คน
-
odds ของการสัมผัสใน control
\[ \frac{9}{11} = 0.818 \]
odds เดียวกันใน control 20 คน
-
odds ratio แบบไม่จับคู่
\[ \mathrm{OR} = \frac{1.50}{0.818} = 1.83 \]
อัตราส่วนของ odds การสัมผัสสองค่า เช่นเดียวกับตาราง case-control แบบไม่จับคู่ทั่วไป
ผลลัพธ์: เมื่อรวมคู่เข้าด้วยกัน odds ของการฝากครรภ์ไม่เพียงพอใน case สูงเป็น 1.83 เท่าของ control
Conditional logistic regression: เปรียบเทียบภายในคู่
conditional logistic regression ถามคำถามที่แคบกว่า เมื่อกำหนดว่าสมาชิกของคู่เป็น case เพียงคนเดียว สิ่งสัมผัสเปลี่ยนว่าสมาชิกคนใดเป็น case อย่างไร มีเพียงคู่ที่สมาชิกต่างกันในสิ่งสัมผัสเท่านั้นที่ตอบได้
ในคู่ที่สอดคล้องกัน สมาชิกทั้งสองมีสิ่งสัมผัสเหมือนกัน สิ่งสัมผัสจึงอธิบายไม่ได้ว่าทำไมคนหนึ่งเป็น case คู่เช่นนี้จึงไม่มีข้อมูลเกี่ยวกับสัมประสิทธิ์ของ ANC ในแบบจำลองที่มีตัวแปรร่วมมากขึ้น คู่ที่สอดคล้องกันยังให้ข้อมูลแก่สัมประสิทธิ์ของตัวแปรร่วมที่ต่างกันภายในคู่ได้
คู่ที่ไม่สอดคล้องกัน 5 คู่มีข้อมูลทั้งหมดเกี่ยวกับสิ่งสัมผัส ใน 4 คู่ในนั้น สมาชิกที่สัมผัสเป็น case ส่วนอีก 1 คู่ สมาชิกที่ไม่สัมผัสเป็น case พิจารณาคู่แบบ 1:1 ที่มีสิ่งสัมผัสแบบสองค่าหนึ่งตัวและไม่มีตัวแปรร่วมอื่น เมื่อนั้น odds ratio แบบมีเงื่อนไข (conditional OR) หรือ odds ratio จับคู่ $\mathrm{OR}_c$ ก็คืออัตราส่วนของสองจำนวนนี้:
$$\mathrm{OR}_c = \frac{4}{1} = 4.00$$ตัวเศษนับคู่ที่เฉพาะ case สัมผัส และตัวส่วนนับคู่ที่เฉพาะ control สัมผัส ตัวอย่างคำนวณด้วยมือถัดไปแสดงว่าอัตราส่วนนี้มาจากไหน [5]
ตัวอย่างคำนวณด้วยมือ: ทำไม intercept ของคู่จึงหักล้างกัน
ตัวอย่างคำนวณด้วยมือ ให้คู่ที่ $k$ มี intercept ของตัวเองคือ $\omega_k$ แทนทุกสิ่งที่สมาชิกของคู่มีร่วมกัน และมีความชันของสิ่งสัมผัสร่วมกันคือ $\beta$ เขียน $x_1$ แทนสิ่งสัมผัสของสมาชิกที่กลายเป็น case และ $x_2$ แทนของสมาชิกอีกคน ภาวะน่าจะเป็น (likelihood) คือความน่าจะเป็นของข้อมูลที่สังเกต อ่านเป็นฟังก์ชันของสัมประสิทธิ์ conditional likelihood ถามว่าสมาชิกคนที่ 1 ไม่ใช่คนที่ 2 เป็น case น่าจะเป็นเท่าใด เมื่อกำหนดว่าคู่มี case เพียงคนเดียว [5]
-
ความน่าจะเป็นของหนึ่งคู่
\[ P_k = \frac{e^{\omega_k + \beta x_1}}{e^{\omega_k + \beta x_1} + e^{\omega_k + \beta x_2}} \]
ภายใต้แบบจำลองโลจิสติก $e^{\omega_k + \beta x}$ คือ odds ที่สมาชิกจะเป็น case เมื่อกำหนดว่าคู่มี case หนึ่งคน ความน่าจะเป็นของ case คือสัดส่วนของ odds ของสมาชิกทั้งสอง
-
intercept ของคู่หักล้างกัน
\[ P_k = \frac{e^{\beta x_1}}{e^{\beta x_1} + e^{\beta x_2}} \]
การหารทั้งเศษและส่วนด้วย $e^{\omega_k}$ ทำให้ intercept หายไป จึงไม่ต้องประมาณสิ่งที่สมาชิกของคู่มีร่วมกัน
-
คู่ที่สอดคล้องกัน
\[ x_1 = x_2: \quad P_k = \frac{e^{\beta x_1}}{2\, e^{\beta x_1}} = 0.5 \]
ไม่ว่า $\beta$ จะเป็นเท่าใด ความน่าจะเป็นคือ 0.5 คู่ที่สอดคล้องกัน 15 คู่จึงไม่มีข้อมูลเกี่ยวกับ $\beta$
-
คู่ที่ไม่สอดคล้องกัน
\[ q = \frac{e^{\beta}}{1 + e^{\beta}} \]
เขียน $q$ แทนความน่าจะเป็นที่สมาชิกที่สัมผัสเป็น case สมาชิกที่ไม่สัมผัสเป็น case ด้วยความน่าจะเป็น $1 - q$
-
conditional likelihood
\[ L(\beta) \propto q^{4}\,(1 - q)^{1} \]
ใน 4 คู่ที่ไม่สอดคล้องกัน สมาชิกที่สัมผัสเป็น case และใน 1 คู่ สมาชิกที่ไม่สัมผัสเป็น case คู่ที่สอดคล้องกันเพียงคูณ $L(\beta)$ ด้วยค่าคงที่
-
ค่าสูงสุด
\[ q = \frac{4}{5} = 0.80 \]
ผลคูณ $q^{4}(1 - q)$ มีค่ามากที่สุดเมื่อ $q$ เท่ากับสัดส่วนของคู่ที่ไม่สอดคล้องกันซึ่งสมาชิกที่สัมผัสเป็น case เมื่อนั้น $1 - q = 0.20$
-
กลับสู่ odds ratio
\[ e^{\beta} = \frac{q}{1 - q} = \frac{0.80}{0.20} = 4.00 \]
แก้ขั้นของคู่ที่ไม่สอดคล้องกันหา $e^{\beta}$ ได้ odds $q/(1 - q)$ ซึ่งคือ odds ratio จับคู่ $\mathrm{OR}_c$
ผลลัพธ์: conditional likelihood ให้ odds ratio จับคู่ 4.00 คืออัตราส่วน 4/1 ของคู่ที่ไม่สอดคล้องกันสองชนิด ส่วนคู่ที่สอดคล้องกันไม่มีบทบาท
คู่ทั้ง 20 คู่เดียวกันใน Stata และ R
สคริปต์จำลองพิมพ์จำนวนคู่สี่แบบ ขยายเป็น 40 แถว และประมาณทั้งสองแบบจำลอง logistic regression ทั่วไปให้ 1.83 ช่วงเชื่อมั่น 95% 0.52 ถึง 6.43 conditional logistic regression ให้ 4.00 ช่วงเชื่อมั่น 95% 0.45 ถึง 35.79 โค้ด Stata ยังประมาณแบบจำลองทั่วไปซ้ำโดยจัดกลุ่ม standard error ตามคู่ ซึ่งจะกลับมากล่าวถึงในความเข้าใจผิดที่พบบ่อยใกล้ท้ายบทความ
ค่าประมาณจากการจับคู่ไม่แม่นยำ สัมประสิทธิ์ของมันคือ 1.386 ซึ่งเป็น log ของ 4.00 บนมาตรวัด log นั้น standard error ของมันคือ $\sqrt{1/4 + 1/1} = 1.118$ สูตรนี้ใช้ได้กับคู่ 1:1 ที่มีสิ่งสัมผัสแบบสองค่าหนึ่งตัวและไม่มีตัวแปรร่วมอื่น
standard error ยกกำลังสอง คือความแปรปรวน เท่ากับผลบวก $1/4 + 1/1$ ส่วนใหญ่มาจากคู่เดียวที่เฉพาะ control สัมผัส เมื่อมีคู่ที่ให้ข้อมูล 5 คู่ ช่วงจึงกว้างจากต่ำกว่า 1 ไปไกลเกิน 4.00 ข้อมูลจึงตัดความเป็นไปได้ที่ไม่มีความสัมพันธ์ไม่ได้
ช่วงแบบ Wald นี้ ซึ่งสร้างจากค่าประมาณและ standard error ของมัน เป็นเพียงค่าประมาณเมื่อมีคู่ที่ให้ข้อมูลน้อยเช่นนี้ ช่วงแบบแม่นตรงซึ่งอิงสัดส่วน 4 ต่อ 1 โดยตรงจะกว้างกว่านี้อีก
Stata เรียกผลลัพธ์ของ clogit ว่า conditional (fixed-effects) logistic regression ผลคงที่ (fixed effects) คือ intercept ของคู่ $\omega_k$ ซึ่ง conditional likelihood กำจัดออกแทนการประมาณ ผลลัพธ์จึงไม่มีพจน์คงที่
clogit() ของ R พิมพ์การเรียก coxph() ซึ่งคือการถดถอย Cox ของ R สำหรับข้อมูลเวลาถึงเหตุการณ์ เพราะประมาณ conditional likelihood ผ่านรูทีนนั้น บทความเรื่องการเลือกกลุ่มควบคุมตามเรื่องลิงก์นั้นและความหมายต่อ odds ratio จับคู่
Stata: ป้อนคู่ทั้ง 20 คู่
clear
input byte anc_case byte anc_ctrl int n_pairs
1 1 8
0 0 7
1 0 4
0 1 1
end
quietly summarize n_pairs
canon pairs.n 0 r(sum)
quietly summarize n_pairs if anc_case == anc_ctrl
canon pairs.concordant 0 r(sum)
quietly summarize n_pairs if anc_case != anc_ctrl
canon pairs.discordant 0 r(sum)
quietly summarize n_pairs if anc_case == 1 & anc_ctrl == 0
canon pairs.case_only_inadequate 0 r(sum)
quietly summarize n_pairs if anc_case == 0 & anc_ctrl == 1
canon pairs.control_only_inadequate 0 r(sum)
expand n_pairs
drop n_pairs
generate int set_id = _n
reshape long anc_, i(set_id) j(role) string
rename anc_ inadequate_anc
generate byte case = role == "case"
bysort set_id: assert _N == 2
bysort set_id (case): assert case[1] == 0 & case[2] == 1
quietly count
canon pairs.rows 0 r(N)
* the pooled 2 x 2 table that ordinary logistic regression sees
tabulate inadequate_anc case
. clear
. input byte anc_case byte anc_ctrl int n_pairs
anc_case anc_ctrl n_pairs
1. 1 1 8
2. 0 0 7
3. 1 0 4
4. 0 1 1
5. end
. quietly summarize n_pairs
. canon pairs.n 0 r(sum)
CANON c7.pairs.n 20
. quietly summarize n_pairs if anc_case == anc_ctrl
. canon pairs.concordant 0 r(sum)
CANON c7.pairs.concordant 15
. quietly summarize n_pairs if anc_case != anc_ctrl
. canon pairs.discordant 0 r(sum)
CANON c7.pairs.discordant 5
. quietly summarize n_pairs if anc_case == 1 & anc_ctrl == 0
. canon pairs.case_only_inadequate 0 r(sum)
CANON c7.pairs.case_only_inadequate 4
. quietly summarize n_pairs if anc_case == 0 & anc_ctrl == 1
. canon pairs.control_only_inadequate 0 r(sum)
CANON c7.pairs.control_only_inadequate 1
. expand n_pairs
(16 observations created)
. drop n_pairs
. generate int set_id = _n
. reshape long anc_, i(set_id) j(role) string
(j = case ctrl)
Data Wide -> Long
-----------------------------------------------------------------------------
Number of observations 20 -> 40
Number of variables 3 -> 3
j variable (2 values) -> role
xij variables:
anc_case anc_ctrl -> anc_
-----------------------------------------------------------------------------
. rename anc_ inadequate_anc
. generate byte case = role == "case"
. bysort set_id: assert _N == 2
. bysort set_id (case): assert case[1] == 0 & case[2] == 1
. quietly count
. canon pairs.rows 0 r(N)
CANON c7.pairs.rows 40
. * the pooled 2 x 2 table that ordinary logistic regression sees
. tabulate inadequate_anc case
inadequate | case
_anc | 0 1 | Total
-----------+----------------------+----------
0 | 11 8 | 19
1 | 9 12 | 21
-----------+----------------------+----------
Total | 20 20 | 40
Stata: logistic regression ทั่วไปกับคู่ที่จับคู่
* ordinary logistic regression ignores the matched sets
logistic case i.inadequate_anc, nolog
. * ordinary logistic regression ignores the matched sets
. logistic case i.inadequate_anc, nolog
Logistic regression Number of obs = 40
LR chi2(1) = 0.91
Prob > chi2 = 0.3413
Log likelihood = -27.27303 Pseudo R2 = 0.0163
----------------------------------------------------------------------------------
case | Odds ratio Std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | 1.833333 1.174412 0.95 0.344 .522362 6.434448
_cons | .7272727 .3379347 -0.69 0.493 .292535 1.808076
----------------------------------------------------------------------------------
Note: _cons estimates baseline odds.
Stata: การประมาณแบบรวมด้วย cluster-robust standard error
* clustering the standard error on the matched set changes the interval, not the odds ratio
logistic case i.inadequate_anc, vce(cluster set_id) nolog
. * clustering the standard error on the matched set changes the interval, not the odds ratio
. logistic case i.inadequate_anc, vce(cluster set_id) nolog
Logistic regression Number of obs = 40
Wald chi2(1) = 1.82
Prob > chi2 = 0.1770
Log pseudolikelihood = -27.27303 Pseudo R2 = 0.0163
(Std. err. adjusted for 20 clusters in set_id)
----------------------------------------------------------------------------------
| Robust
case | Odds ratio std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | 1.833333 .823118 1.35 0.177 .7604568 4.419858
_cons | .7272727 .1778604 -1.30 0.193 .4503282 1.174534
----------------------------------------------------------------------------------
Note: _cons estimates baseline odds.
R: คู่ที่จับคู่และการประมาณ logistic regression ทั่วไป
pairs <- data.frame(anc_case = c(1, 0, 1, 0), anc_ctrl = c(1, 0, 0, 1), n_pairs = c(8, 7, 4, 1))
canon("pairs.n", sum(pairs$n_pairs), 0)
canon("pairs.concordant", sum(pairs$n_pairs[pairs$anc_case == pairs$anc_ctrl]), 0)
canon("pairs.discordant", sum(pairs$n_pairs[pairs$anc_case != pairs$anc_ctrl]), 0)
canon("pairs.case_only_inadequate", sum(pairs$n_pairs[pairs$anc_case == 1 & pairs$anc_ctrl == 0]), 0)
canon("pairs.control_only_inadequate", sum(pairs$n_pairs[pairs$anc_case == 0 & pairs$anc_ctrl == 1]), 0)
pp <- pairs[rep(seq_len(nrow(pairs)), pairs$n_pairs), ]
pp$set_id <- seq_len(nrow(pp))
m <- rbind(data.frame(set_id = pp$set_id, case = 1L, inadequate_anc = pp$anc_case),
data.frame(set_id = pp$set_id, case = 0L, inadequate_anc = pp$anc_ctrl))
m <- m[order(m$set_id, m$case), ]
stopifnot(all(table(m$set_id) == 2), all(tapply(m$case, m$set_id, sum) == 1))
canon("pairs.rows", nrow(m), 0)
# the pooled 2 x 2 table that ordinary logistic regression sees
print(table(inadequate_anc = m$inadequate_anc, case = m$case))
a <- sum(m$case == 1 & m$inadequate_anc == 1)
bb <- sum(m$case == 1 & m$inadequate_anc == 0)
cc <- sum(m$case == 0 & m$inadequate_anc == 1)
dd <- sum(m$case == 0 & m$inadequate_anc == 0)
canon("pooled.cases_inadequate", a, 0)
canon("pooled.cases_adequate", bb, 0)
canon("pooled.controls_inadequate", cc, 0)
canon("pooled.controls_adequate", dd, 0)
canon("pooled.odds_cases", a / bb)
canon("pooled.odds_controls", cc / dd)
# ordinary logistic regression ignores the matched sets
fp <- glm(case ~ inadequate_anc, family = binomial, data = m)
print(summary(fp)$coefficients)
bp <- coef(fp)[["inadequate_anc"]]
sp <- sqrt(vcov(fp)["inadequate_anc", "inadequate_anc"])
canon("pooled.or", exp(bp))
canon("pooled.se", sp)
canon("pooled.or_lb", exp(bp - z * sp))
canon("pooled.or_ub", exp(bp + z * sp))
> pairs <- data.frame(anc_case = c(1, 0, 1, 0), anc_ctrl = c(1,
+ 0, 0, 1), n_pairs = c(8, 7, 4, 1))
> canon("pairs.n", sum(pairs$n_pairs), 0)
CANON c7.pairs.n 20
> canon("pairs.concordant", sum(pairs$n_pairs[pairs$anc_case ==
+ pairs$anc_ctrl]), 0)
CANON c7.pairs.concordant 15
> canon("pairs.discordant", sum(pairs$n_pairs[pairs$anc_case !=
+ pairs$anc_ctrl]), 0)
CANON c7.pairs.discordant 5
> canon("pairs.case_only_inadequate", sum(pairs$n_pairs[pairs$anc_case ==
+ 1 & pairs$anc_ctrl == 0]), 0)
CANON c7.pairs.case_only_inadequate 4
> canon("pairs.control_only_inadequate", sum(pairs$n_pairs[pairs$anc_case ==
+ 0 & pairs$anc_ctrl == 1]), 0)
CANON c7.pairs.control_only_inadequate 1
> pp <- pairs[rep(seq_len(nrow(pairs)), pairs$n_pairs),
+ ]
> pp$set_id <- seq_len(nrow(pp))
> m <- rbind(data.frame(set_id = pp$set_id, case = 1,
+ inadequate_anc = pp$anc_case), data.frame(set_id = pp$set_id,
+ case = 0, inadequate_anc = pp$anc_ctrl))
> m <- m[order(m$set_id, m$case), ]
> stopifnot(all(table(m$set_id) == 2), all(tapply(m$case,
+ m$set_id, sum) == 1))
> canon("pairs.rows", nrow(m), 0)
CANON c7.pairs.rows 40
> print(table(inadequate_anc = m$inadequate_anc, case = m$case))
case
inadequate_anc 0 1
0 11 8
1 9 12
> a <- sum(m$case == 1 & m$inadequate_anc == 1)
> bb <- sum(m$case == 1 & m$inadequate_anc == 0)
> cc <- sum(m$case == 0 & m$inadequate_anc == 1)
> dd <- sum(m$case == 0 & m$inadequate_anc == 0)
> canon("pooled.cases_inadequate", a, 0)
CANON c7.pooled.cases_inadequate 12
> canon("pooled.cases_adequate", bb, 0)
CANON c7.pooled.cases_adequate 8
> canon("pooled.controls_inadequate", cc, 0)
CANON c7.pooled.controls_inadequate 9
> canon("pooled.controls_adequate", dd, 0)
CANON c7.pooled.controls_adequate 11
> canon("pooled.odds_cases", a/bb)
CANON c7.pooled.odds_cases 1.5000
> canon("pooled.odds_controls", cc/dd)
CANON c7.pooled.odds_controls 0.8182
> fp <- glm(case ~ inadequate_anc, family = binomial,
+ data = m)
> print(summary(fp)$coefficients)
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.3184537 0.4646602 -0.6853476 0.4931246
inadequate_anc 0.6061358 0.6405884 0.9462172 0.3440378
> bp <- coef(fp)[["inadequate_anc"]]
> sp <- sqrt(vcov(fp)["inadequate_anc", "inadequate_anc"])
> canon("pooled.or", exp(bp))
CANON c7.pooled.or 1.8333
> canon("pooled.se", sp)
CANON c7.pooled.se 0.6406
> canon("pooled.or_lb", exp(bp - z * sp))
CANON c7.pooled.or_lb 0.5224
> canon("pooled.or_ub", exp(bp + z * sp))
CANON c7.pooled.or_ub 6.4344
Stata: conditional logistic regression กับคู่ที่จับคู่
* conditional logistic regression compares the case with its own control inside each matched set
clogit case i.inadequate_anc, group(set_id) or nolog
. * conditional logistic regression compares the case with its own control inside each matched set
. clogit case i.inadequate_anc, group(set_id) or nolog
Conditional (fixed-effects) logistic regression Number of obs = 40
LR chi2(1) = 1.93
Prob > chi2 = 0.1650
Log likelihood = -12.89922 Pseudo R2 = 0.0695
----------------------------------------------------------------------------------
case | Odds ratio Std. err. z P>|z| [95% conf. interval]
-----------------+----------------------------------------------------------------
1.inadequate_anc | 4 4.472136 1.24 0.215 .4470826 35.78757
----------------------------------------------------------------------------------
R: clogit() กับคู่ที่จับคู่
# conditional logistic regression compares the case with its own control inside each matched set
fc <- clogit(case ~ inadequate_anc + strata(set_id), data = m)
print(summary(fc))
> fc <- clogit(case ~ inadequate_anc + strata(set_id),
+ data = m)
> print(summary(fc))
Call:
coxph(formula = Surv(rep(1, 40L), case) ~ inadequate_anc + strata(set_id),
data = m, method = "exact")
n= 40, number of events= 20
coef exp(coef) se(coef) z Pr(>|z|)
inadequate_anc 1.386 4.000 1.118 1.24 0.215
exp(coef) exp(-coef) lower .95 upper .95
inadequate_anc 4 0.25 0.4471 35.79
Concordance= 0.575 (se = 0.075 )
Likelihood ratio test= 1.93 on 1 df, p=0.2
Wald test = 1.54 on 1 df, p=0.2
Score (logrank) test = 1.8 on 1 df, p=0.2
ทำไม 1.83 และ 4.00 จึงต่างกัน
ตัวเลขสองค่านี้ตอบคำถามต่างกัน odds ratio แบบไม่จับคู่เปรียบเทียบ odds ของการสัมผัสของ case ทั้งหมดกับ control ทั้งหมด odds ratio แบบจับคู่เปรียบเทียบการสัมผัสภายในคู่ ระหว่างทารกที่มีสิ่งที่สมาชิกของคู่มีร่วมกันเหมือนกัน
การรวมทำให้การเปรียบเทียบภายในคู่เลือนไป คู่ที่สอดคล้องกันทุกคู่เพิ่มสิ่งสัมผัสเดียวกันเข้าทั้งสองคอลัมน์ของตารางรวม คู่ที่สัมผัสทั้งสองคนเพิ่มสมาชิกที่สัมผัส 8 คนเข้าแต่ละคอลัมน์ คู่ที่ไม่สัมผัสทั้งสองคนเพิ่มสมาชิกที่ไม่สัมผัส 7 คนเข้าแต่ละคอลัมน์ ในที่นี้ 15 จาก 20 คู่สอดคล้องกัน
เมื่อ case 12 จาก 20 คน และ control 9 จาก 20 คนฝากครรภ์ไม่เพียงพอ ถ้าสิ่งสัมผัสของสมาชิกสองคนในคู่ไม่สัมพันธ์กัน จะมีคู่ที่สอดคล้องกันราว 10 คู่ (0.60 x 0.45 + 0.40 x 0.55 = 0.49 ดังนั้น 0.49 x 20 = 9.8) เทียบกับ 15 คู่ที่พบจริง ผลนี้เข้ากันได้กับการที่ตัวแปรที่ใช้จับคู่สัมพันธ์กับสิ่งสัมผัส แม้คู่ 20 คู่จะแสดงเรื่องนี้ไม่ได้
เมื่อตัวแปรที่ใช้จับคู่ ในที่นี้คืออายุเป็นวัน สัมพันธ์กับสิ่งสัมผัส การจับคู่จะดึงสิ่งสัมผัสของ control เข้าหาของ case การวิเคราะห์ที่ไม่สนใจชุดจึงมักทำให้ odds ratio ลำเอียงเข้าหา 1 [5] ส่วนที่เพิ่มร่วมกันข้างต้นคือวิธีที่ความลำเอียงนั้นปรากฏในตารางรวม
เมื่อตัวแปรที่ใช้จับคู่ไม่สัมพันธ์กับสิ่งสัมผัส การวิเคราะห์ทั้งสองแบบประมาณปริมาณเดียวกัน แล้วต่างกันหลัก ๆ ที่ความแม่นยำ เมื่อมีคู่ที่ให้ข้อมูลเพียง 5 คู่ โอกาสอาจเปิดช่องว่างระหว่างสองค่าได้เช่นกัน
non-collapsibility ซึ่งนิยามไว้ด้านล่างสำหรับแบบจำลอง cohort คือช่องว่างระหว่าง odds ratio ภายในกลุ่มย่อยกับ odds ratio ที่เฉลี่ยข้ามกลุ่มย่อย ซึ่งเกิดขึ้นเมื่อกลุ่มย่อยทำนายผลลัพธ์ได้แม้ไม่มีตัวกวน ในที่นี้ทุกคู่มี case หนึ่งคนและ control หนึ่งคนที่อายุเป็นวันเท่ากัน ตัวคู่เองและอายุจึงทำนายไม่ได้ว่าใครเป็น case ช่องว่างระหว่าง 1.83 กับ 4.00 จึงไม่ได้มาจาก non-collapsibility ช่องว่างนี้สะท้อนแบบแผนการศึกษา คือการเปรียบเทียบภายในคู่กับการเปรียบเทียบแบบรวม ซึ่งการจับคู่อาจดึงสิ่งสัมผัสของ control เข้าหาของ case ร่วมกับความบังเอิญเมื่อคู่ที่ให้ข้อมูลมีน้อยเช่นนี้
แล้วตัวเลขไหนควรอยู่ในรายงาน การศึกษาออกแบบมาเพื่อเปรียบเทียบภายในคู่ odds ratio จับคู่จึงมักเป็นค่าที่ควรรายงาน ควรมาพร้อมช่วงเชื่อมั่นของมัน 0.45 ถึง 35.79 ซึ่งเข้ากันได้ทั้งกับการไม่มีความสัมพันธ์และกับความสัมพันธ์ที่แรง
การจับคู่จัดการได้เฉพาะตัวแปรที่ใช้จับคู่ ตัวกวนอื่นยังต้องอยู่ในแบบจำลองมีเงื่อนไขเป็นตัวแปรร่วม และสัมประสิทธิ์ของมันประมาณภายในชุดที่จับคู่ odds ratio จับคู่ที่ไม่มีตัวแปรเหล่านั้นปรับเฉพาะอายุเป็นวัน และไม่ปรับอย่างอื่น
link เดียว แบบจำลองสามแบบ: ตารางสำหรับตัดสินใจ
ตารางรวมแบบจำลองทั้งสาม อ่านแต่ละแถวจากซ้าย ชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่างมาก่อน ตามด้วยแบบจำลองและคำสั่ง คอลัมน์สุดท้ายบอกความหมายของสัมประสิทธิ์ที่ยกกำลัง ซึ่งเป็นจุดที่แบบจำลองทั้งสามต่างกันมากที่สุด
แถวสุดท้ายใช้ อัตราส่วนฮาซาร์ด (hazard ratio) ซึ่งเปรียบเทียบอัตราการเกิดเหตุการณ์ ณ ขณะใดขณะหนึ่งในผู้ที่ยังไม่เกิดเหตุการณ์ ภายใต้ ข้อสมมติฮาซาร์ดเป็นสัดส่วน (proportional hazards) อัตราส่วนนั้นคงที่ตลอดช่วงติดตาม
| ชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่าง | แบบจำลอง | Stata | R | สัมประสิทธิ์ที่ยกกำลัง |
|---|---|---|---|---|
| ผลลัพธ์สองค่า คนอิสระต่อกัน เช่น cohort หรือการศึกษา case-control ที่ไม่จับคู่ | Logistic regression | logistic y i.x | glm(y ~ x, family = binomial) | odds ratio ของผู้สัมผัสเทียบกับผู้ไม่สัมผัส โดยคงตัวแปรร่วมอื่นไว้ ในการศึกษา case-control สิ่งที่ค่านี้ประมาณขึ้นกับวิธีสุ่มกลุ่มควบคุม (เป็น odds ratio เฉพาะเมื่อกลุ่มควบคุมมาจากผู้ที่ยังไม่เกิดผลลัพธ์เมื่อสิ้นสุดการติดตาม) |
| หมวดที่มีลำดับ คนอิสระต่อกัน | Ordinal logistic regression (proportional odds) | ologit y i.x | MASS::polr(factor(y, ordered = TRUE) ~ x) | odds ratio ค่าเดียวสำหรับหมวดที่สูงขึ้น สมมติว่าเท่ากันทุกจุดตัด ตรวจด้วย Brant test |
| ผลลัพธ์สองค่า (case หรือ control) ชุดที่จับคู่ รวมถึง control ที่สุ่มจากกลุ่มเสี่ยง | Conditional logistic regression | clogit case i.x, group(set_id) or | survival::clogit(case ~ x + strata(set_id)) | odds ratio ภายในชุดที่จับคู่ เมื่อสุ่มจากกลุ่มเสี่ยงจะประมาณ hazard ratio ภายใต้ proportional hazards ตัวแปรที่ใช้จับคู่ไม่มีสัมประสิทธิ์ |
สิ่งที่ odds ratio ทั้งสามมีร่วมกัน
ในแบบจำลองทั้งสาม สัมประสิทธิ์คือผลต่างของ log odds และเลขชี้กำลังของมันคือ odds ratio สิ่งที่ odds ratio นั้นมีเงื่อนไขอยู่ต่างกัน logistic regression และ ordinal logistic regression มีเงื่อนไขตามตัวแปรร่วมอื่น ส่วน conditional logistic regression มีเงื่อนไขตามชุดที่จับคู่ด้วย ชุดเงื่อนไขแต่ละชุดอาจนิยาม ปริมาณเป้าหมายของการประมาณ (estimand) ของตัวเอง คือปริมาณที่แน่ชัดซึ่งการศึกษาตั้งใจประมาณ
ใน birth cohort odds ratio ของ ANC ที่ 2.00 และ 2.01 มีเงื่อนไข (conditional) ตามการคลอดก่อนกำหนด การคลอดก่อนกำหนดเป็น ปัจจัยพยากรณ์โรค (prognostic factor) คือลักษณะที่ทำนายผลลัพธ์ แบบจำลองที่ไม่มีมันจะให้ odds ratio ระดับประชากร (marginal) ซึ่งเฉลี่ยตามสถานะการคลอดก่อนกำหนด odds ratio ไม่ยุบรวมได้ (non-collapsible) ดังนั้นทั้งสองค่าไม่จำเป็นต้องเท่ากัน แม้การคลอดก่อนกำหนดจะไม่สัมพันธ์กับ ANC [6]
odds ratio แบบมีเงื่อนไขกับแบบระดับประชากร รวมทั้ง hazard ratio ทั้งสองแบบ อาจต่างกันได้แม้ไม่มีตัวกวนเลย ในการทดลองแบบสุ่ม การปรับด้วยปัจจัยพยากรณ์โรคที่มีผลแรงมักทำให้ odds ratio ขยับห่างจาก 1 มากขึ้น นั่นคือการเปลี่ยน estimand ไม่ใช่การกำจัดความลำเอียง กรณีการทดลองแบบสุ่ม [7] คำนวณไว้ในตอนของชุดบทความนี้ว่าด้วย odds ratio แบบมีเงื่อนไขและระดับประชากร
รายงานจึงควรระบุว่า odds ratio แต่ละค่ามีเงื่อนไขตามอะไร ควบคู่กับแบบจำลอง การฝากครรภ์ไม่เพียงพออาจนำไปสู่การคลอดก่อนกำหนด ทำให้การคลอดก่อนกำหนดเป็น ตัวแปรสื่อกลาง (mediator) คือขั้นบนเส้นทางไปสู่ภาวะแทรกซ้อน ถ้าเป็นเช่นนั้น odds ratio 2.00 และ 2.01 ละเส้นทางนั้นไป ต่างจาก ผลรวม (total effect) ของ ANC ที่รวมมันไว้ ตัวแปรร่วมใดควรอยู่ในเงื่อนไขเป็นคำถามเชิงเหตุผล ซึ่งควรตัดสินด้วยแผนภาพเชิงเหตุผลก่อนการวิเคราะห์
ความเข้าใจผิดที่พบบ่อยและวิธีแก้
-
"clogit ให้ odds ratio ที่ใหญ่กว่า จึงแม่นยำกว่า"
ขนาดของค่าประมาณไม่ได้บอกอะไรเกี่ยวกับความแม่นยำ ในที่นี้ 4.00 จากการจับคู่อยู่ห่างจาก 1 มากกว่า 1.83 ที่ไม่จับคู่ แต่ช่วงของมัน 0.45 ถึง 35.79 กว้างกว่า 0.52 ถึง 6.43 มาก ในข้อมูลอื่น odds ratio จับคู่อาจเป็นค่าที่เล็กกว่า
วิธีแก้: ตัวเลขสองค่านี้ตอบคำถามต่างกัน: odds ratio แบบจับคู่เปรียบเทียบการได้รับสิ่งสัมผัสภายในคู่ ส่วน odds ratio แบบไม่คำนึงถึงการจับคู่รวมทุกคู่เข้าด้วยกัน ควรเลือกแบบจำลองจากรูปแบบการศึกษาและ estimand แล้วรายงานช่วงความเชื่อมั่น
-
"ologit ใช้ได้กับผลลัพธ์แบบลำดับทุกชนิด"
ใน cohort จำลอง แบบจำลองแบบลำดับให้การคลอดก่อนกำหนดมี odds ratio ค่าเดียวคือ 1.84 ขณะที่ odds ratio ตามจุดตัดอยู่ที่ 1.40 ถึง 3.31 Brant test ชี้การคลอดก่อนกำหนด (P = 0.003) และไม่ชี้ ANC (P = 0.999)
วิธีแก้: ologit สมมติ odds ratio ค่าเดียวที่ทุกจุดตัด ให้ตรวจข้อสมมตินั้น และบอกว่าทำอะไรเมื่อข้อสมมตินั้นไม่เป็นจริง
-
"clogit ทิ้งคู่ที่สอดคล้องกัน"
ทุกคู่เข้าสู่ conditional likelihood เมื่อมีสิ่งสัมผัสแบบสองค่าหนึ่งตัวและไม่มีตัวแปรร่วมอื่น คู่ที่สอดคล้องกันให้ความน่าจะเป็น 0.5 ไม่ว่าสัมประสิทธิ์เป็นเท่าใด จึงไม่มีข้อมูลเกี่ยวกับสัมประสิทธิ์นั้น เมื่อมีตัวแปรร่วมมากขึ้น มันยังให้ข้อมูลแก่ตัวแปรที่ต่างกันภายในคู่ได้
วิธีแก้: เก็บทุกชุดที่จับคู่ไว้ในข้อมูล และรายงานว่ามีกี่ชุดที่ไม่สอดคล้องกันในสิ่งสัมผัส เพราะชุดเหล่านั้นกำหนดความแม่นยำของ odds ratio ของสิ่งสัมผัส
-
"cluster-robust standard error แก้การวิเคราะห์ที่ไม่สนใจการจับคู่ได้"
cluster-robust standard error ยอมให้หน่วยสังเกตภายในกลุ่ม ในที่นี้คือคู่ มีความสัมพันธ์กันได้ การจัดกลุ่มตามคู่ใน logistic regression แบบรวม ดังผลลัพธ์ Stata ที่แสดงพร้อมโค้ดของคู่ที่จับคู่ ทำให้ช่วงเชื่อมั่น 95% แคบลงเป็น 0.76 ถึง 4.42 จาก 0.52 ถึง 6.43 ช่วงแคบลงในที่นี้เพราะแต่ละคู่มี case เพียงหนึ่งคนโดยการออกแบบ odds ratio ยังเป็น 1.83 ซึ่งยังเป็นการเปรียบเทียบแบบรวม ช่วงทั้งสองจึงไม่ตอบคำถามภายในคู่
วิธีแก้: วิเคราะห์ชุดที่จับคู่ด้วย conditional logistic regression standard error ที่ต่างไปเปลี่ยนช่วง ไม่เปลี่ยนค่าประมาณหรือคำถามที่มันตอบ
-
"logistic regression ทั่วไปที่ใส่ตัวบ่งชี้หนึ่งตัวต่อคู่ให้ odds ratio จับคู่"
ตัวบ่งชี้ของแต่ละคู่ถูกประมาณจากทารกเพียงสองคน วิธี maximum likelihood ทั่วไป ซึ่งเป็นวิธีประมาณค่าที่ใช้กันบ่อย จึงผลัก odds ratio ให้ห่างจาก 1 ยิ่งขึ้น สำหรับคู่ 1:1 ที่มีสิ่งสัมผัสแบบสองค่าหนึ่งตัว odds ratio ของแบบจำลองนั้นคือกำลังสองของค่าจับคู่ [5]
วิธีแก้: ใช้ clogit วางเงื่อนไขตามชุดที่จับคู่ conditional likelihood กำจัด intercept ของคู่แทนการประมาณมัน
-
"odds ratio 4.00 หมายถึงความเสี่ยงของภาวะแทรกซ้อนเป็นสี่เท่า"
odds ratio จับคู่คืออัตราส่วนของ odds ภายในคู่ ไม่ใช่อัตราส่วนของความเสี่ยง ในที่นี้ control แต่ละคนมาจากทารกที่ยังอยู่ในกลุ่มเสี่ยง ณ เวลาที่ case เกิดเหตุการณ์ ภายใต้ proportional hazards แบบแผนนี้ทำให้ค่านั้นประมาณ hazard ratio ไม่ใช่ risk ratio ที่ 28 วัน บทความเรื่องการเลือกกลุ่มควบคุมคำนวณเรื่องนี้ไว้
วิธีแก้: ระบุมาตรวัดที่แบบแผนการศึกษาประมาณ และรายงาน 4.00 พร้อมช่วงเชื่อมั่น 95% ของมัน 0.45 ถึง 35.79
ทำอย่างไรในการวิเคราะห์ของคุณเอง
- เขียนชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่างลงในแผนการวิเคราะห์ก่อนเลือกคำสั่ง
- สำหรับผลลัพธ์แบบลำดับ ควรพิจารณา ordinal logistic regression แทนการยุบระดับ และรายงาน odds ratio เฉพาะจุดตัดควบคู่กับค่าร่วม
- ตรวจ proportional odds ด้วย Brant test ทั้งโดยรวมและรายตัวแปรร่วม และระบุว่าทำอะไรกับตัวแปรร่วมที่ถูกชี้
- สำหรับ control ที่จับคู่หรือสุ่มจากกลุ่มเสี่ยง ให้สร้างรหัสของชุดตั้งแต่ขั้นสุ่มตัวอย่าง ไม่ใช่ภายหลัง และวิเคราะห์ชุดเหล่านั้นด้วย conditional logistic regression
- รายงานจำนวนชุดที่ไม่สอดคล้องกันและช่วงเชื่อมั่น 95% ควบคู่กับ odds ratio จับคู่
- ถ้าการวิเคราะห์แบบจับคู่กับแบบไม่จับคู่ไม่ตรงกัน ควรพิจารณารายงานการวิเคราะห์แบบจับคู่เป็นค่าประมาณหลัก odds ratio แบบไม่จับคู่อธิบายตารางรวมและมองข้ามการจับคู่ จึงอาจลำเอียงเข้าหา 1 เมื่อตัวแปรที่ใช้จับคู่สัมพันธ์กับสิ่งสัมผัส
อภิธานศัพท์
- logit link
- link ที่เป็น ln{p/(1 - p)} หรือ log odds ซึ่งผูกตัวทำนายเชิงเส้นเข้ากับความน่าจะเป็น สัมประสิทธิ์ของมันคือ log odds ratio
- cut point (จุดตัดระหว่างหมวดของผลลัพธ์แบบลำดับ)
- เส้นแบ่งระหว่างหมวดสองหมวดที่อยู่ติดกันของผลลัพธ์แบบลำดับ จุดตัดแต่ละจุดมี intercept ของตัวเองใน ordinal logistic regression
- cumulative logit
- log odds ของการอยู่ที่หมวดหนึ่งหรือสูงกว่าของผลลัพธ์แบบลำดับ มีหนึ่งค่าต่อจุดตัดแต่ละจุด
- ordinal logistic regression
- แบบจำลองสำหรับผลลัพธ์แบบลำดับที่ประมาณ cumulative logit ที่จุดตัดทุกจุดโดยใช้ชุดความชันร่วมกันชุดเดียว
- proportional odds assumption (ข้อสมมติ proportional odds)
- ข้อสมมติที่ว่า odds ratio ของตัวแปรร่วมเท่ากันที่ทุกจุดตัดของผลลัพธ์แบบลำดับ
- Brant test
- Wald test ของข้อสมมติ proportional odds ที่เปรียบเทียบความชันของแบบจำลองแยกตามจุดตัด ทั้งโดยรวมและทีละตัวแปรร่วม
- matched set (ชุดที่จับคู่)
- case หนึ่งคนพร้อม control ที่สุ่มหรือจับคู่ให้มัน คือหนึ่งคู่ในการจับคู่ 1:1
- concordant pair (คู่ที่สอดคล้องกัน)
- คู่ที่จับคู่ซึ่ง case และ control มีสิ่งสัมผัสเหมือนกัน เมื่อมีสิ่งสัมผัสแบบสองค่าหนึ่งตัว คู่นี้ไม่ให้ข้อมูลเกี่ยวกับสิ่งสัมผัสใน conditional logistic regression
- discordant pair (คู่ที่ไม่สอดคล้องกัน)
- คู่ที่จับคู่ซึ่ง case และ control ต่างกันในสิ่งสัมผัส เมื่อจับคู่ 1:1 มีสิ่งสัมผัสแบบสองค่าหนึ่งตัวและไม่มีตัวแปรร่วมอื่น odds ratio จับคู่คืออัตราส่วนของคู่สองชนิด
- conditional likelihood
- ความน่าจะเป็นที่สมาชิกคนใดของชุดที่จับคู่เป็น case เมื่อกำหนดว่าชุดมี case เพียงหนึ่งคน intercept ของชุดหักล้างกันหมดไป
- conditional logistic regression
- แบบจำลองโลจิสติกที่ให้ intercept แยกกันสำหรับชุดที่จับคู่แต่ละชุด ประมาณผ่าน conditional likelihood เพื่อให้ case แต่ละคนถูกเปรียบเทียบกับ control ของมันเท่านั้น
- risk-set sampling (การสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยง)
- การเลือก control ของ case แต่ละคนจากผู้ที่ยังอยู่ในกลุ่มเสี่ยง ณ เวลาที่ case เกิดเหตุการณ์
- non-collapsibility (ความไม่ยุบรวมได้)
- คุณสมบัติของ odds ratio และ hazard ratio ที่ค่าระดับประชากรอาจต่างจากค่ามีเงื่อนไขร่วมโดยไม่มีตัวกวนเลย
- mediator (ตัวแปรสื่อกลาง)
- ตัวแปรบนเส้นทางเหตุและผลจากสิ่งสัมผัสไปสู่ผลลัพธ์ ซึ่งผลส่วนหนึ่งอาจผ่านไป
- cluster-robust standard error
- standard error ที่ยอมให้หน่วยสังเกตภายในกลุ่ม เช่นคู่ที่จับคู่ สัมพันธ์กันได้ เปลี่ยนช่วงแต่ไม่เปลี่ยนค่าประมาณ
เอกสารอ้างอิง
- StataCorp. Stata Base Reference Manual: [R] logit. Stata documentation (PDF) at stata.com. Accessed 11 October 2026. https://www.stata.com/manuals/rlogit.pdf
- StataCorp. Stata Base Reference Manual: [R] ologit. Stata documentation (PDF) at stata.com. Accessed 11 October 2026. https://www.stata.com/manuals/rologit.pdf
- Brant R. Assessing proportionality in the proportional odds model for ordinal logistic regression. Biometrics. 1990;46(4):1171 ff. doi:10.2307/2532457 https://doi.org/10.2307/2532457
- StataCorp. Stata Base Reference Manual: [R] clogit. Stata documentation (PDF) at stata.com. Accessed 11 October 2026. https://www.stata.com/manuals/rclogit.pdf
- Breslow NE, Day NE. Statistical Methods in Cancer Research, Volume I: The Analysis of Case-Control Studies. IARC Scientific Publications No. 32. Lyon: International Agency for Research on Cancer; 1980. https://publications.iarc.who.int/Book-And-Report-Series/Iarc-Scientific-Publications/Statistical-Methods-In-Cancer-Research-Volume-I-The-Analysis-Of-Case-Control-Studies-1980
- Greenland S, Pearl J, Robins JM. Confounding and collapsibility in causal inference. Stat Sci. 1999;14(1):29-46. doi:10.1214/ss/1009211805 https://doi.org/10.1214/ss/1009211805
- Gail MH, Wieand S, Piantadosi S. Biased estimates of treatment effect in randomized experiments with nonlinear regressions and omitted covariates. Biometrika. 1984;71(3):431-444. doi:10.1093/biomet/71.3.431 https://doi.org/10.1093/biomet/71.3.431
ประเด็นสำคัญ
- logistic regression, ordinal logistic regression และ conditional logistic regression ใช้ logit link ร่วมกัน สัมประสิทธิ์ที่ยกกำลังแต่ละค่าจึงเป็น odds ratio
- ชนิดของผลลัพธ์และแบบแผนการสุ่มตัวอย่าง ที่ระบุก่อนพิมพ์คำสั่งใด ๆ เป็นตัวกำหนดว่าจะใช้แบบจำลองใดในสามแบบ
- ordinal logistic regression สมมติ odds ratio ค่าเดียวที่ทุกจุดตัด ให้ตรวจทั้งโดยรวมและรายตัวแปรร่วมด้วย Brant test
- ในคู่ที่จับคู่ 1:1 ที่มีสิ่งสัมผัสแบบสองค่าหนึ่งตัว มีเพียงคู่ที่ไม่สอดคล้องกันที่ให้ข้อมูลแก่ odds ratio ของสิ่งสัมผัส 4 ต่อ 1 ให้ 4.00 ขณะที่การรวมคู่ให้ 1.83
- ค่าประมาณจากการจับคู่ตอบคำถามของแบบแผนการศึกษา แต่เมื่อมีคู่ที่ไม่สอดคล้องกัน 5 คู่ ช่วง 0.45 ถึง 35.79 แสดงว่าไม่แม่นยำ
อ่านต่อในวิกิ: [[conditional-vs-marginal-or-th]] [[case-control-sampling-schemes-stata-and-r-th]] [[ordinal-logistic-regression]] [[conditional-logistic-regression]] [[multinomial-vs-ordinal-regression]]