การเลือกกลุ่มควบคุมใน Stata และ R: Cohort เดียว กลุ่มควบคุมสามชุด ค่าประมาณสามค่า
On this page
Read the English version
บทคัดย่อ
การศึกษาแบบ case-control วัดสิ่งสัมผัสในผู้ป่วยและในตัวอย่างของกลุ่มควบคุม และสิ่งที่ผลคูณไขว้ (cross-product) ของมัน ซึ่งเป็นอัตราส่วนที่โปรแกรมพิมพ์เป็น odds ratio ประมาณได้นั้นขึ้นกับว่ากลุ่มควบคุมมาจากที่ใด บทความนี้สุ่มกลุ่มควบคุมสามชุดจาก cohort จำลองหนึ่งชุดที่มี 10,000 คนและผู้ป่วย 1,560 ราย กลุ่มควบคุมจากผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตามประมาณ odds ratio ของ cohort คือ 3.14 กลุ่มควบคุมจากทุกคนตอนเริ่มต้น รวมผู้ที่ภายหลังกลายเป็นผู้ป่วย ประมาณอัตราส่วนความเสี่ยง (risk ratio) คือ 2.50 กลุ่มควบคุมจากผู้ที่ยังมีความเสี่ยงอยู่ ณ เวลาเกิดเหตุการณ์ของผู้ป่วยแต่ละราย เมื่อวิเคราะห์เป็นชุดจับคู่ ประมาณอัตราส่วนฮาซาร์ด (hazard ratio) ซึ่งเป็นอัตราส่วนของอัตราการเกิดเหตุการณ์ ณ ขณะใดขณะหนึ่ง คือ 2.78 ตัวอย่างที่บันทึกไว้อย่างละหนึ่งชุดให้ odds ratio 2.99 อัตราส่วนความเสี่ยง (ผลคูณไขว้) 2.53 และอัตราส่วนฮาซาร์ด 2.61 โดยช่วงเชื่อมั่นแต่ละช่วงครอบคลุมเป้าหมายของมัน บทความนี้สรุปว่าจะอ่านค่าประมาณแบบ case-control ได้ก็ต่อเมื่อระบุชุดกลุ่มควบคุมและวิธีวิเคราะห์ของมันแล้ว
งบคลังตัวอย่างชีวภาพ กับสามวิธีเลือกกลุ่มควบคุม
ทีมวิจัยทีมหนึ่งมี cohort ผู้ใหญ่ 10,000 คนที่ติดตามเป็นเวลา 5 ปีโดยไม่มีผู้ใดหลุดจากการติดตาม ในช่วงนั้นมี 1,560 คนเป็นโรคที่ศึกษา เลือดที่เก็บไว้ตอนรับเข้าร่วมการศึกษาบอกได้ว่าใครมีสิ่งสัมผัสที่สงสัย แต่การตรวจวัด (assay) ทุกครั้งมีค่าใช้จ่าย งบประมาณรองรับการตรวจวัด 4,680 ครั้ง คือผู้ป่วยทั้ง 1,560 รายบวกกลุ่มควบคุมสองคนต่อผู้ป่วยหนึ่งราย รวมกลุ่มควบคุม 3,120 คน
การศึกษาแบบ case-control (case-control study) ทำงานดังนี้ คือวัดสิ่งสัมผัสในผู้ป่วย และในตัวอย่างของสมาชิก cohort คนอื่น ซึ่งเรียกว่ากลุ่มควบคุม (control) แทนการวัดในทุกคน กฎที่ใช้เลือกกลุ่มควบคุม รวมกับกลุ่มควบคุมที่ได้ เรียกว่าชุดกลุ่มควบคุม (control series) เพื่อนร่วมงานสามคนเสนอกฎสามแบบ:
- สุ่มกลุ่มควบคุมจากผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม เรียกว่าการเลือกกลุ่มควบคุมแบบสะสม (cumulative sampling หรือ exclusive sampling);
- สุ่มกลุ่มควบคุมจากทั้ง cohort ตามที่เป็นอยู่ตอนเริ่มต้น รวมผู้ที่ภายหลังกลายเป็นผู้ป่วย เรียกว่าการเลือกกลุ่มควบคุมจากฐานประชากรทั้งหมด (case-base sampling หรือ inclusive sampling);
- ทุกครั้งที่เกิดผู้ป่วยหนึ่งราย สุ่มกลุ่มควบคุมจากผู้ที่ยังมีความเสี่ยงอยู่ ณ ขณะนั้น เรียกว่าการสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยง (risk-set sampling หรือ concurrent sampling)
ทั้งสามแบบจะตรวจวัดผู้ป่วยชุดเดียวกัน 1,560 ราย แต่จะไม่ได้ตัวเลขเดียวกัน และแต่ละแบบจะถูกต้องสำหรับปริมาณคนละอย่าง งานวิจัยแบบ case-control ที่ตีพิมพ์จำนวนมากไม่ได้ระบุชัดว่าสุ่มกลุ่มควบคุมอย่างไร จึงไม่ชัดว่า odds ratio ของงานเหล่านั้นประมาณอะไร [1] บทความนี้ทดลองใช้ทั้งสามแบบกับ cohort จำลองหนึ่งชุดใน Stata และ R แล้วตรวจค่าประมาณแต่ละค่ากับเป้าหมายของมัน
Cohort จำลองหนึ่งชุด
cohort ในบทความนี้เป็นข้อมูลจำลอง ในผู้ที่ได้รับสิ่งสัมผัส 2,000 คน มี 600 คนเป็นผู้ป่วยตลอด 8,500 คนปี ในผู้ที่ไม่ได้รับสิ่งสัมผัส 8,000 คน มี 960 คนเป็นผู้ป่วยตลอด 37,600 คนปี คนปี (person-years) คือผลรวมของเวลาที่แต่ละคนอยู่ในการติดตามโดยยังไม่เป็นโรค คนที่ไม่เคยเป็นผู้ป่วยจึงเพิ่มให้ 5 ปี
ชุดกลุ่มควบคุมทุกชุดในบทความนี้สุ่มจาก cohort ที่กำหนดไว้ชัดเจนชุดเดียว ทุกแบบจึงซ้อนอยู่ใน cohort นั้น ตามการใช้ทั่วไป การศึกษาแบบ nested case-control (nested case-control design) หมายถึงกฎแบบที่สาม กลุ่มควบคุมมาจากกลุ่มเสี่ยง (risk set) ของผู้ป่วยแต่ละราย คือผู้ที่ยังมีความเสี่ยงอยู่ ณ ขณะที่ผู้ป่วยรายนั้นเกิดโรค ยอดรวมของ cohort ตรงกับตัวอย่างในบทความเรื่องการศึกษาแบบ nested case-control ซึ่งอธิบายเหตุผลของการออกแบบไว้โดยละเอียด
บทความนี้เน้นที่คำสั่งและการอ่านผลลัพธ์ ไฟล์ข้อมูลแต่ละไฟล์มีหนึ่งแถวต่อคนหรือต่อคนที่ถูกสุ่ม โดยมีตัวระบุ id สิ่งสัมผัส exposed (1 หรือ 0) และสถานะผู้ป่วย สำหรับแต่ละแบบ Stata และ R วิเคราะห์กลุ่มควบคุมที่สุ่มไว้ชุดเดียวกัน ซึ่งในที่นี้เรียกว่าตัวอย่างที่บันทึกไว้ (frozen sample) ไฟล์เหล่านี้ไม่ได้เผยแพร่ ช่องโค้ดแสดงคำสั่งและผลลัพธ์จริง
กลุ่มควบคุมเป็นตัวแทนของอะไร
ทุกการวิเคราะห์ด้านล่างเริ่มจากตาราง 2 คูณ 2 ของสิ่งสัมผัสในผู้ป่วยและในกลุ่มควบคุม ให้ $a$ และ $b$ แทนผู้ป่วยที่ได้รับและไม่ได้รับสิ่งสัมผัส และ $c$ และ $d$ แทนกลุ่มควบคุมที่ได้รับและไม่ได้รับสิ่งสัมผัส ผลคูณไขว้ (cross-product) ของตาราง หรือเรียกว่าอัตราส่วนผลคูณไขว้ (cross-product ratio) คือ
$$\text{cross-product} = \frac{a \times d}{b \times c} = \frac{a / b}{c / d}$$รูปแบบที่สองบอกว่ามันคืออะไร คือ odds ของการได้รับสิ่งสัมผัสในผู้ป่วย หารด้วย odds ของการได้รับสิ่งสัมผัสในกลุ่มควบคุม ฝั่งผู้ป่วยเป็น 600 ต่อ 960 ในทุกแบบ คำตอบจึงขึ้นกับว่ากลุ่มควบคุมเป็นตัวแทนของอะไร คำสั่ง cc และ logistic ของ Stata เมื่อมีสิ่งสัมผัสที่มีสองค่า (binary) หนึ่งตัว ให้ตัวเลขนี้เหมือนกัน และพิมพ์ออกมาเป็น odds ratio ส่วนใน R ค่าเอกซ์โพเนนเชียลของสัมประสิทธิ์จากการถดถอยโลจิสติกชุดเดียวกันให้ตัวเลขเดียวกัน
ผลคูณไขว้ประมาณอะไรขึ้นกับสิ่งที่กลุ่มควบคุมสะท้อน [2, 3, 4]:
- ผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม: odds ratio;
- ทุกคนตอนเริ่มต้น: อัตราส่วนความเสี่ยง (risk ratio);
- เวลาของผู้ที่อยู่ในภาวะเสี่ยง (person-time at risk): อัตราส่วนอัตราอุบัติการณ์ (rate ratio) และเมื่อวิเคราะห์เป็นชุดจับคู่ (matched set) คือผู้ป่วยแต่ละรายพร้อมกลุ่มควบคุมที่สุ่มให้รายนั้น จะได้อัตราส่วนฮาซาร์ด (hazard ratio) ซึ่งเป็นอัตราส่วนของอัตราการเกิดเหตุการณ์ ณ ขณะใดขณะหนึ่ง
สองความหมายแรกสมมติให้เป็นcohort ปิด (closed cohort) คือทุกคนเข้ามาตอนเริ่มต้นและถูกติดตามเป็นระยะเวลาคงที่เท่ากันโดยไม่มีผู้หลุดจากการติดตาม cohort นี้เป็น closed cohort เป้าหมายสองค่านั้นจึงเป็น odds ratio 5 ปี และอัตราส่วนความเสี่ยง 5 ปีของ cohort นี้ การสุ่มจากกลุ่มเสี่ยงไม่ต้องการเงื่อนไขนั้น เพราะกลุ่มเสี่ยงแต่ละชุดมีเฉพาะผู้ที่ยังอยู่ในการติดตาม ณ เวลาที่ผู้ป่วยรายนั้นเกิดโรค เป้าหมายของมัน คืออัตราส่วนอัตราอุบัติการณ์หรืออัตราส่วนฮาซาร์ด จึงยังชัดเจนเมื่อผู้คนเข้าร่วมในเวลาต่างกันหรือหลุดจากการติดตาม
เป้าหมายเหล่านี้ไม่ต้องการให้โรคพบน้อย ความพบน้อยสำคัญเฉพาะเมื่ออ่าน odds ratio เป็นอัตราส่วนความเสี่ยง
ตัวอย่างคำนวณด้วยมือ: cohort เดียว ตัวหารสามแบบ
ตัวอย่างคำนวณด้วยมือ จากยอดรวมของ cohort ด้านบน (ข้อมูลจำลอง) ตัวห้อย 1 แทนกลุ่มที่ได้รับสิ่งสัมผัส และตัวห้อย 0 แทนกลุ่มที่ไม่ได้รับ ขั้นที่ 1 ถึง 3 ใช้ทั้ง cohort ขั้นที่ 4 ถึง 6 ใช้ตารางประกอบที่ปัดเศษแล้ว ตารางละกลุ่มควบคุม 1,000 คน ซึ่งเป็นตารางเดียวกับในบทความเรื่องการศึกษาแบบ nested case-control
-
ความเสี่ยงและอัตราส่วนความเสี่ยง
\[ \text{risk}_1 = \frac{600}{2000} = 0.300, \;\; \text{risk}_0 = \frac{960}{8000} = 0.120 \]
อัตราส่วนความเสี่ยงคือ 0.300 / 0.120 = 2.50 ตัวหารของมันคือคนตอนเริ่มต้น
-
Odds และ odds ratio
\[ \text{odds}_1 = \frac{600}{1400} = 0.429, \;\; \text{odds}_0 = \frac{960}{7040} = 0.136 \]
Odds หารจำนวนผู้ป่วยด้วย 1,400 และ 7,040 คนที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม odds ratio คือ (600 × 7,040) / (1,400 × 960) = 3.14
-
อัตราและอัตราส่วนอัตราอุบัติการณ์
\[ \text{rate}_1 = \frac{600}{8500} = 0.0706, \;\; \text{rate}_0 = \frac{960}{37600} = 0.0255 \]
อัตราคือจำนวนผู้ป่วยต่อคนปี จากอัตราที่ยังไม่ปัดเศษ อัตราส่วนอัตราอุบัติการณ์คือ 2.76
-
กลุ่มควบคุมแบบ exclusive
\[ \frac{600 \times 834}{960 \times 166} = 3.14 \]
สุ่มกลุ่มควบคุม 1,000 คนจาก 8,440 คนที่ยังไม่เป็นโรค แบ่งสัดส่วนเหมือน 1,400 และ 7,040 คือได้รับสิ่งสัมผัส 166 คน ไม่ได้รับ 834 คน ผลคูณไขว้ให้ odds ratio
-
กลุ่มควบคุมแบบ inclusive
\[ \frac{600 \times 800}{960 \times 200} = 2.50 \]
สุ่มกลุ่มควบคุม 1,000 คนจากทั้ง 10,000 คนตอนเริ่มต้น แบ่งสัดส่วนเหมือน 2,000 และ 8,000 คือได้รับสิ่งสัมผัส 200 คน ไม่ได้รับ 800 คน ผลคูณไขว้ให้อัตราส่วนความเสี่ยง
-
กลุ่มควบคุมแบบ concurrent
\[ \frac{600 \times 816}{960 \times 184} = 2.77 \]
สุ่มกลุ่มควบคุม 1,000 คนตามสัดส่วนของ 8,500 และ 37,600 คนปี คือได้รับสิ่งสัมผัส 184 คน ไม่ได้รับ 816 คน หลังปัดเป็นจำนวนคนเต็ม ผลคูณไขว้ 2.77 ต่างจากอัตราส่วนอัตราอุบัติการณ์ 2.76 เพียงเพราะการปัดเศษนั้น
ผลลัพธ์: ผู้ป่วยชุดเดียวกัน 1,560 รายให้ 3.14, 2.50 หรือ 2.77 ขึ้นกับเพียงว่ากลุ่มควบคุมเป็นตัวแทนของตัวหารแบบใด
ตารางรวมเหล่านี้เป็นตัวอย่างประกอบ การสุ่มจากกลุ่มเสี่ยงจริงเก็บชุดจับคู่ไว้และวิเคราะห์ทีละชุด ดังที่หัวข้อกลุ่มควบคุมแบบ concurrent แสดง
ทั้ง cohort: อัตราส่วนสี่ค่าเป็นเกณฑ์เทียบ
ตารางด้านล่างแสดงอัตราส่วนสี่ค่าของทั้ง cohort พร้อมช่วงเชื่อมั่น 95% (CI) จากแบบจำลองในช่องโค้ด ค่าเหล่านี้เป็นเกณฑ์เทียบของตัวอย่างทั้งสามชุด กลุ่มควบคุมแบบ exclusive มุ่งไปที่ odds ratio 5 ปี กลุ่มควบคุมแบบ inclusive มุ่งไปที่อัตราส่วนความเสี่ยง 5 ปี และกลุ่มควบคุมจากกลุ่มเสี่ยงมุ่งไปที่อัตราส่วนฮาซาร์ด
ฮาซาร์ด (hazard) คืออัตราการเกิดเหตุการณ์ ณ ขณะใดขณะหนึ่งในผู้ที่ยังไม่เป็นโรค แบบจำลอง Cox (Cox model) ซึ่งเป็นการถดถอยมาตรฐานสำหรับข้อมูลเวลาถึงเหตุการณ์ สรุปอัตราส่วนของฮาซาร์ดของสองกลุ่มตลอดการติดตามเป็นอัตราส่วนฮาซาร์ดค่าเดียว เมื่ออัตราของแต่ละกลุ่มคงที่ตลอดเวลา อัตราส่วนฮาซาร์ดกับอัตราส่วนอัตราอุบัติการณ์เป็นปริมาณเดียวกัน
ใน cohort จำลองนี้ ฮาซาร์ดของกลุ่มที่ได้รับสิ่งสัมผัสสูงขึ้นตลอดการติดตาม ขณะที่ของกลุ่มที่ไม่ได้รับเกือบคงที่ อัตราส่วนของสองฮาซาร์ดจึงเลื่อนสูงขึ้นเอง และ 2.78 คือค่าสรุปของมันจากแบบจำลอง Cox ตลอดการติดตาม ส่วนอัตราส่วนอัตราอุบัติการณ์ 2.76 รวมผู้ป่วยและคนปีทั้งหมดเข้าด้วยกัน สองค่าสรุปจึงต่างกันเล็กน้อย ค่าเวลาซ้ำกัน (ties) ไม่เกี่ยวกับความต่างนี้ เพราะไม่มีผู้ป่วยสองรายที่มีเวลาเกิดเหตุการณ์ตรงกัน ดังที่หัวผลลัพธ์ของ Stata ว่า 'Cox regression with no ties' ยืนยัน
| มาตรวัด (หน่วย) | ได้รับสิ่งสัมผัส | ไม่ได้รับสิ่งสัมผัส | อัตราส่วน (ช่วงเชื่อมั่น 95%) |
|---|---|---|---|
| ความเสี่ยง (ผู้ป่วยต่อคนตอนเริ่มต้น) | 0.300 | 0.120 | อัตราส่วนความเสี่ยง 5 ปี 2.50 (2.29 ถึง 2.73) |
| Odds (ผู้ป่วยต่อคนที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม) | 0.429 | 0.136 | odds ratio 5 ปี 3.14 (2.80 ถึง 3.53) |
| อัตรา (ผู้ป่วยต่อคนปี) | 0.0706 | 0.0255 | อัตราส่วนอัตราอุบัติการณ์ 2.76 (2.50 ถึง 3.06) |
| ฮาซาร์ด (แบบจำลอง Cox) | สูงขึ้นตลอดการติดตาม | เกือบคงที่ตลอดการติดตาม | อัตราส่วนฮาซาร์ด 2.78 (2.51 ถึง 3.07) |
Stata: อัตราส่วนฮาซาร์ดของทั้ง cohort
* 5. Full-cohort hazard ratio: the target of risk-set sampling
stset time, failure(case) id(id)
stcox exposed, nolog. * 5. Full-cohort hazard ratio: the target of risk-set sampling
. stset time, failure(case) id(id)
Survival-time data settings
ID variable: id
Failure event: case!=0 & case<.
Observed time interval: (time[_n-1], time]
Exit on or before: failure
--------------------------------------------------------------------------
10,000 total observations
0 exclusions
--------------------------------------------------------------------------
10,000 observations remaining, representing
10,000 subjects
1,560 failures in single-failure-per-subject data
46,100 total analysis time at risk and under observation
At risk from t = 0
Earliest observed entry t = 0
Last observed exit t = 5
. stcox exposed, nolog
Failure _d: case
Analysis time _t: time
ID variable: id
Cox regression with no ties
No. of subjects = 10,000 Number of obs = 10,000
No. of failures = 1,560
Time at risk = 46,100
LR chi2(1) = 343.79
Log likelihood = -14067.77 Prob > chi2 = 0.0000
------------------------------------------------------------------------------
_t | Haz. ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.776434 .1445672 19.61 0.000 2.507066 3.074743
------------------------------------------------------------------------------R: อัตราส่วนสี่ค่าของทั้ง cohort
# 3. Full-cohort risk ratio (log-binomial) and odds ratio (logistic)
fit_rr <- glm(case ~ exposed, family = binomial(link = "log"), data = cohort)
canon_ratio("cohort.rr", coef(fit_rr), vcov(fit_rr))
fit_or <- glm(case ~ exposed, family = binomial, data = cohort)
canon_ratio("cohort.or", coef(fit_or), vcov(fit_or))
# 4. Full-cohort rate ratio: Poisson model with log person-years as the offset
fit_irr <- glm(case ~ exposed + offset(log(time)), family = poisson, data = cohort)
canon_ratio("cohort.irr", coef(fit_irr), vcov(fit_irr))
# 5. Full-cohort hazard ratio: the target of risk-set sampling (Breslow ties, as in Stata)
fit_hr <- coxph(Surv(time, case) ~ exposed, data = cohort, ties = "breslow")
canon_ratio("cohort.hr", coef(fit_hr), vcov(fit_hr))> fit_rr <- glm(case ~ exposed, family = binomial(link = "log"),
+ data = cohort)
> canon_ratio("cohort.rr", coef(fit_rr), vcov(fit_rr))
CANON w7.cohort.rr 2.5000
CANON w7.cohort.rr_lb 2.2861
CANON w7.cohort.rr_ub 2.7340
> fit_or <- glm(case ~ exposed, family = binomial, data = cohort)
> canon_ratio("cohort.or", coef(fit_or), vcov(fit_or))
CANON w7.cohort.or 3.1429
CANON w7.cohort.or_lb 2.7958
CANON w7.cohort.or_ub 3.5330
> fit_irr <- glm(case ~ exposed + offset(log(time)),
+ family = poisson, data = cohort)
> canon_ratio("cohort.irr", coef(fit_irr), vcov(fit_irr))
CANON w7.cohort.irr 2.7647
CANON w7.cohort.irr_lb 2.4966
CANON w7.cohort.irr_ub 3.0616
> fit_hr <- coxph(Surv(time, case) ~ exposed, data = cohort,
+ ties = "breslow")
> canon_ratio("cohort.hr", coef(fit_hr), vcov(fit_hr))
CANON w7.cohort.hr 2.7764
CANON w7.cohort.hr_lb 2.5071
CANON w7.cohort.hr_ub 3.0747กลุ่มควบคุมแบบ exclusive: ผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม
การเลือกแบบสะสม (cumulative หรือ exclusive sampling) สุ่มกลุ่มควบคุมเฉพาะจากผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม ใน cohort นี้กลุ่มดังกล่าวมี 8,440 คน เป็นผู้ได้รับสิ่งสัมผัส 1,400 คน และไม่ได้รับ 7,040 คน odds ของการได้รับสิ่งสัมผัสในกลุ่มควบคุมประมาณ odds ในกลุ่มนั้น ผลคูณไขว้จึงประมาณ odds ratio 5 ปีของ cohort คือ 3.14 [2, 4] odds ratio นั้นเข้าใกล้อัตราส่วนความเสี่ยงเฉพาะเมื่อโรคพบน้อย และความเสี่ยง 0.300 ในกลุ่มที่ได้รับสิ่งสัมผัสนั้นห่างจากคำว่าน้อยมาก
ตัวอย่างที่บันทึกไว้สุ่มกลุ่มควบคุม 3,120 คนจากกลุ่มนั้น เป็นผู้ได้รับสิ่งสัมผัส 539 คน และไม่ได้รับ 2,581 คน เมื่อรวมกับผู้ป่วย 1,560 ราย ผลคูณไขว้คือ (600 × 2,581) / (960 × 539) = 2.99 ในไฟล์ตัวอย่าง y เป็น 1 สำหรับผู้ป่วยและ 0 สำหรับกลุ่มควบคุม และ logistic ประมาณ
ในสมการนี้ $\operatorname{logit}(p) = \ln\{p/(1 - p)\}$ คือ log odds และ $i$ ระบุแถว $\beta_0$ คือค่าจุดตัดแกน (intercept) และ $\beta_1$ คือสัมประสิทธิ์ของสิ่งสัมผัส เมื่อมีสิ่งสัมผัสที่มีสองค่าหนึ่งตัว $\exp(\beta_1)$ เท่ากับผลคูณไขว้ คือ 2.99 ช่วงเชื่อมั่น 95% ของมัน 2.61 ถึง 3.44 ครอบคลุมเป้าหมาย 3.14
Stata: ตัวอย่างแบบ exclusive
cc y exposed
canon4 exclusive.or_lb_exact.stata r(lb_or)
canon4 exclusive.or_ub_exact.stata r(ub_or)
logistic y exposed, nolog. cc y exposed
Proportion
| Exposed Unexposed | Total exposed
-----------------+------------------------+------------------------
Cases | 600 960 | 1560 0.3846
Controls | 539 2581 | 3120 0.1728
-----------------+------------------------+------------------------
Total | 1139 3541 | 4680 0.2434
| |
| Point estimate | [95% conf. interval]
|------------------------+------------------------
Odds ratio | 2.992811 | 2.600961 3.443775 (exact)
Attr. frac. ex. | .6658659 | .6155268 .709621 (exact)
Attr. frac. pop | .2561023 |
+-------------------------------------------------
chi2(1) = 253.49 Pr>chi2 = 0.0000
. canon4 exclusive.or_lb_exact.stata r(lb_or)
CANON w7.exclusive.or_lb_exact.stata 2.6010
. canon4 exclusive.or_ub_exact.stata r(ub_or)
CANON w7.exclusive.or_ub_exact.stata 3.4438
. logistic y exposed, nolog
Logistic regression Number of obs = 4,680
LR chi2(1) = 243.62
Prob > chi2 = 0.0000
Log likelihood = -2857.0778 Pseudo R2 = 0.0409
------------------------------------------------------------------------------
y | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.992811 .2105856 15.58 0.000 2.607267 3.435366
_cons | .3719489 .014061 -26.16 0.000 .345386 .4005546
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.กลุ่มควบคุมแบบ inclusive: ทุกคนตอนเริ่มต้น
การเลือกจากฐานประชากรทั้งหมด (case-base หรือ inclusive sampling) สุ่มกลุ่มควบคุมจากทั้ง cohort ตามที่เป็นอยู่ตอนเริ่มต้น ไม่ว่าภายหลังคนเหล่านั้นจะเป็นอย่างไร กลุ่มควบคุมจึงสะท้อนผู้ได้รับสิ่งสัมผัส 2,000 คนและผู้ไม่ได้รับ 8,000 คนตอนเริ่มต้น ซึ่งเป็นตัวหารของความเสี่ยง 5 ปีสองค่า ผลคูณไขว้จึงประมาณอัตราส่วนความเสี่ยง 5 ปี คือ 2.50 [2, 4]
กลุ่มควบคุมบางคนกลายเป็นผู้ป่วยภายหลัง และยังคงอยู่ในกลุ่ม ในตัวอย่างที่บันทึกไว้ 493 จาก 3,120 คนที่สุ่มจากทั้ง 10,000 คนกลายเป็นผู้ป่วยระหว่างการติดตาม คนเหล่านี้ปรากฏสองครั้ง ครั้งละหนึ่งบทบาท กลุ่มควบคุมมีผู้ได้รับสิ่งสัมผัส 619 คน และไม่ได้รับ 2,501 คน และอัตราส่วนความเสี่ยง (ผลคูณไขว้) คือ 2.53
คนที่ปรากฏสองครั้งทำให้แถวไม่เป็นอิสระต่อกัน ซึ่งกระทบค่าคลาดเคลื่อนมาตรฐาน (standard error) แต่ไม่กระทบค่าประมาณ สคริปต์จึงจัดกลุ่มค่าคลาดเคลื่อนมาตรฐานตามแต่ละคน (cluster) ด้วย vce(cluster id) ใน Stata และ vcovCL() จากแพ็กเกจ sandwich ใน R ช่วงเชื่อมั่น 95% คือ 2.25 ถึง 2.83 ครอบคลุมเป้าหมาย 2.50 Stata ยังคงตั้งหัวคอลัมน์ว่า 'Odds ratio' เพราะคำสั่งคำนวณผลคูณไขว้ และสิ่งที่ทำให้มันเป็นอัตราส่วนความเสี่ยงคือการออกแบบ
การเก็บตัวอย่างสุ่มตอนเริ่มต้นไว้เป็น subcohort แล้วติดตามตามเวลา ให้การออกแบบแบบ case-cohort ซึ่งมีวิธีวิเคราะห์แบบถ่วงน้ำหนักของตัวเอง [5, 6] บทความนี้วิเคราะห์ตัวอย่างแบบ inclusive เพียงเป็นตาราง 2 คูณ 2
Stata: ตัวอย่างแบบ inclusive
cc y exposed
* A person can appear as a case and as a control, so the standard error is clustered on the person
logistic y exposed, vce(cluster id) nolog. cc y exposed
Proportion
| Exposed Unexposed | Total exposed
-----------------+------------------------+------------------------
Cases | 600 960 | 1560 0.3846
Controls | 619 2501 | 3120 0.1984
-----------------+------------------------+------------------------
Total | 1219 3461 | 4680 0.2605
| |
| Point estimate | [95% conf. interval]
|------------------------+------------------------
Odds ratio | 2.525242 | 2.201884 2.896123 (exact)
Attr. frac. ex. | .6039984 | .5458434 .6547108 (exact)
Attr. frac. pop | .2323071 |
+-------------------------------------------------
chi2(1) = 187.22 Pr>chi2 = 0.0000
. * A person can appear as a case and as a control, so the standard error is cl
> ustered on the person
. logistic y exposed, vce(cluster id) nolog
Logistic regression Number of obs = 4,680
Wald chi2(1) = 250.45
Prob > chi2 = 0.0000
Log pseudolikelihood = -2888.3749 Pseudo R2 = 0.0304
(Std. err. adjusted for 4,187 clusters in id)
------------------------------------------------------------------------------
| Robust
y | Odds ratio std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.525242 .1478116 15.83 0.000 2.251537 2.83222
_cons | .3838465 .0132611 -27.72 0.000 .3587156 .4107379
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.กลุ่มควบคุมแบบ concurrent: กลุ่มเสี่ยง ณ เวลาเกิดเหตุการณ์ของผู้ป่วยแต่ละราย
การสุ่มจากกลุ่มเสี่ยง (risk-set หรือ concurrent sampling) ทำงานผ่านผู้ป่วยตามลำดับเวลา เมื่อเกิดผู้ป่วยหนึ่งราย จะสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยงของรายนั้น คือทุกคนที่ยังอยู่ในการติดตามและยังไม่เป็นโรค ณ ขณะนั้น ผู้ป่วยหนึ่งรายกับกลุ่มควบคุมของเขาเป็นชุดจับคู่ (matched set) กลุ่มควบคุมอาจกลายเป็นผู้ป่วยภายหลัง และคนหนึ่งอาจถูกสุ่มให้หลายชุด
คนที่ไม่เป็นโรคนานกว่าจะอยู่ในกลุ่มเสี่ยงมากกว่า กลุ่มควบคุมโดยรวมจึงสะท้อนเวลาของผู้ที่อยู่ในภาวะเสี่ยงได้โดยประมาณ การวิเคราะห์แบบจับคู่ ไม่ใช่ตารางรวม เป็นตัวที่ให้ค่าประมาณ เมื่อวิเคราะห์ทีละชุด การออกแบบนี้ประมาณอัตราส่วนฮาซาร์ดโดยไม่ต้องสมมติว่าโรคพบน้อย [3, 7]
ใน cohort นี้ เป้าหมายคืออัตราส่วนฮาซาร์ดจากแบบจำลอง Cox ของทั้ง cohort คือ 2.78 การวิเคราะห์แบบจับคู่มุ่งไปที่ค่านี้พอดีเมื่ออัตราส่วนของฮาซาร์ดคงที่ตลอดเวลา ซึ่งเรียกว่าฮาซาร์ดเป็นสัดส่วน (proportional hazards) เมื่ออัตราส่วนนี้เลื่อน เช่นในที่นี้ การวิเคราะห์แบบจับคู่มุ่งไปที่ค่าสรุปที่เกือบเป็นค่าเดียวกับแบบจำลอง Cox ของทั้ง cohort
การถดถอยโลจิสติกแบบมีเงื่อนไข (conditional logistic regression) คือการวิเคราะห์ทีละชุดนั้น เป็นแบบจำลองโลจิสติกที่มีจุดตัดแกนของตัวเองให้ทุกชุดจับคู่:
$$\operatorname{logit} \Pr(y_{ij} = 1) = \alpha_j + \beta\,\mathrm{exposed}_{ij}$$ในสมการนี้ $j$ ระบุชุดจับคู่ $i$ ระบุคนภายในชุด และ $y_{ij}$ เป็น 1 สำหรับผู้ป่วย พจน์ $\alpha_j$ ดูดซับทุกอย่างที่สมาชิกของชุด $j$ มีร่วมกัน รวมถึงเวลาที่ผู้ป่วยเกิดโรค การกำหนดเงื่อนไขว่าแต่ละชุดมีผู้ป่วยเพียงหนึ่งรายกำจัด $\alpha_j$ ทุกตัว ผู้ป่วยแต่ละรายจึงถูกเปรียบเทียบกับกลุ่มควบคุมของตัวเองเท่านั้น ส่วนที่เหลือมีรูปแบบเดียวกับส่วนที่กลุ่มเสี่ยงนั้นมีต่อแบบจำลอง Cox ดังนั้น $\exp(\beta)$ ประมาณอัตราส่วนฮาซาร์ด [7]
ตัวอย่างกลุ่มเสี่ยงที่บันทึกไว้มี 1,560 ชุดจับคู่ ชุดละผู้ป่วยหนึ่งรายกับกลุ่มควบคุมสองคน รวมกลุ่มควบคุมที่ได้รับสิ่งสัมผัส 609 คน และไม่ได้รับ 2,511 คน ในกลุ่มควบคุมเหล่านี้ 244 คนกลายเป็นผู้ป่วยภายหลัง คำสั่ง clogit y exposed, group(set) or nolog โดย set คือเลขของชุดจับคู่ ให้อัตราส่วนฮาซาร์ด 2.61 (ช่วงเชื่อมั่น 95% คือ 2.27 ถึง 3.01) ซึ่งครอบคลุมเป้าหมาย 2.78 ฟังก์ชัน clogit() ของ R กับชุดเดียวกันให้ค่าเดียวกัน
Stata: ตัวอย่างกลุ่มเสี่ยงที่บันทึกไว้
clogit y exposed, group(set) or nolog
canon4 riskset.hr exp(_b[exposed])
canon4 riskset.hr_lb exp(_b[exposed]-invnormal(0.975)*_se[exposed])
canon4 riskset.hr_ub exp(_b[exposed]+invnormal(0.975)*_se[exposed])
* For contrast only: pooling the risk sets into one unmatched 2 by 2 table ignores the sampling design
logistic y exposed, nolog. clogit y exposed, group(set) or nolog
Conditional (fixed-effects) logistic regression Number of obs = 4,680
LR chi2(1) = 188.91
Prob > chi2 = 0.0000
Log likelihood = -1619.3805 Pseudo R2 = 0.0551
------------------------------------------------------------------------------
y | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.613864 .186032 13.50 0.000 2.273536 3.005136
------------------------------------------------------------------------------
. canon4 riskset.hr exp(_b[exposed])
CANON w7.riskset.hr 2.6139
. canon4 riskset.hr_lb exp(_b[exposed]-invnormal(0.975)*_se[exposed])
CANON w7.riskset.hr_lb 2.2735
. canon4 riskset.hr_ub exp(_b[exposed]+invnormal(0.975)*_se[exposed])
CANON w7.riskset.hr_ub 3.0051
. * For contrast only: pooling the risk sets into one unmatched 2 by 2 table ig
> nores the sampling design
. logistic y exposed, nolog
Logistic regression Number of obs = 4,680
LR chi2(1) = 188.17
Prob > chi2 = 0.0000
Log likelihood = -2884.8011 Pseudo R2 = 0.0316
------------------------------------------------------------------------------
y | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.57697 .1775796 13.74 0.000 2.251402 2.949619
_cons | .3823178 .0145075 -25.34 0.000 .3549152 .4118361
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.การสุ่มของ Stata เอง: stset, sttocc และ clogit
Stata สุ่มกลุ่มเสี่ยงเองได้ในสามขั้น ขั้นแรก stset time, failure(case) id(id) ประกาศเวลาติดตาม เหตุการณ์ และตัวระบุแต่ละคน จากนั้น sttocc exposed, number(2) nodots สุ่มกลุ่มควบคุมสองคนจากกลุ่มเสี่ยงของผู้ป่วยแต่ละราย และเก็บตัวแปร exposed ไว้ในข้อมูลชุดใหม่ สุดท้าย clogit _case exposed, group(_set) or nolog วิเคราะห์ชุดจับคู่
sttocc เพิ่มตัวแปรสามตัว คือ _case (1 สำหรับผู้ป่วย 0 สำหรับกลุ่มควบคุม) _set (ชุดจับคู่) และ _time (เวลาที่ผู้ป่วยเกิดเหตุการณ์) การวิเคราะห์ต้องใช้ _case และ _set เพราะตัวบ่งชี้ผู้ป่วยเดิมของ cohort จะนับกลุ่มควบคุมที่ภายหลังเป็นโรคว่าเป็นผู้ป่วย การสุ่มของ Stata ให้อัตราส่วนฮาซาร์ด 2.64 (ช่วงเชื่อมั่น 95% คือ 2.30 ถึง 3.03) ซึ่งต่างจาก 2.61 ของตัวอย่างที่บันทึกไว้ เพราะเป็นการสุ่มคนละรอบจากกลุ่มเสี่ยงชุดเดียวกัน
Stata: sttocc แล้วตามด้วย clogit
* 6. Stata's own risk-set draw: sttocc keeps each case and samples 2 controls from its risk set
sttocc exposed, number(2) nodots
count
display "CANON w7.sttocc.n.stata " r(N)
tabulate _case exposed
* Conditional logistic regression keeps each matched set intact and estimates the hazard ratio
clogit _case exposed, group(_set) or nolog. tabulate _case exposed
0 for |
controls; |
1 for | exposed
cases | 0 1 | Total
-----------+----------------------+----------
0 | 2,534 586 | 3,120
1 | 960 600 | 1,560
-----------+----------------------+----------
Total | 3,494 1,186 | 4,680
. * Conditional logistic regression keeps each matched set intact and estimates
> the hazard ratio
. clogit _case exposed, group(_set) or nolog
Conditional (fixed-effects) logistic regression Number of obs = 4,680
LR chi2(1) = 198.56
Prob > chi2 = 0.0000
Log likelihood = -1614.5561 Pseudo R2 = 0.0579
------------------------------------------------------------------------------
_case | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
exposed | 2.642912 .1856319 13.84 0.000 2.303012 3.032977
------------------------------------------------------------------------------ชุดกลุ่มควบคุมสามชุดเทียบกัน
| ชุดกลุ่มควบคุม | สุ่มกลุ่มควบคุมจาก | สิ่งที่การวิเคราะห์ประมาณ | เป้าหมาย (ค่าจริงในข้อมูลจำลอง) | ตัวอย่างที่บันทึกไว้ (ช่วงเชื่อมั่น 95%) | Stata | R |
|---|---|---|---|---|---|---|
| แบบสะสม (exclusive) | 8,440 คนที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม | Odds ratio | 3.14 | odds ratio 2.99 (2.61 ถึง 3.44) | cc y exposed; logistic y exposed, nolog | glm(y ~ exposed, family = binomial, data = ex) |
| จากฐานประชากรทั้งหมด (inclusive) | ทั้ง 10,000 คนตอนเริ่มต้น | อัตราส่วนความเสี่ยง | 2.50 | อัตราส่วนความเสี่ยง (ผลคูณไขว้) 2.53 (2.25 ถึง 2.83) | logistic y exposed, vce(cluster id) nolog | glm(y ~ exposed, family = binomial, data = inc) แล้วใช้ vcovCL() จัดกลุ่มตาม id |
| จากกลุ่มเสี่ยง (concurrent) | ผู้ที่ยังมีความเสี่ยงอยู่ ณ เวลาเกิดเหตุการณ์ของผู้ป่วยแต่ละราย | อัตราส่วนฮาซาร์ด (คืออัตราส่วนอัตราอุบัติการณ์เมื่ออัตราคงที่) | 2.78 | อัตราส่วนฮาซาร์ด 2.61 (2.27 ถึง 3.01) | clogit y exposed, group(set) or nolog | clogit(y ~ exposed + strata(set), data = rs) |
R: ชุดกลุ่มควบคุมสามชุด
# 6. Exclusive controls (non-cases at the end of follow-up): the cross-product estimates the odds ratio
ex <- dat("W7_exclusive.csv")
canon_n("exclusive.n", nrow(ex))
canon_n("exclusive.controls_exposed", sum(ex$y == 0 & ex$exposed == 1))
canon_n("exclusive.controls_unexposed", sum(ex$y == 0 & ex$exposed == 0))
tab_ex <- table(exposed = factor(ex$exposed, 1:0), y = factor(ex$y, 1:0))
print(tab_ex)
fit_ex <- glm(y ~ exposed, family = binomial, data = ex)
canon_ratio("exclusive.or", coef(fit_ex), vcov(fit_ex))
# Exact limits from the noncentral hypergeometric distribution (point value here is the conditional estimate)
fx <- fisher.test(tab_ex)
canon("exclusive.or_cmle.r", unname(fx$estimate))
canon("exclusive.or_lb_exact.r", fx$conf.int[1])
canon("exclusive.or_ub_exact.r", fx$conf.int[2])
# 7. Inclusive controls (everyone at baseline): the cross-product estimates the risk ratio
inc <- dat("W7_inclusive.csv")
canon_n("inclusive.n", nrow(inc))
canon_n("inclusive.controls_exposed", sum(inc$y == 0 & inc$exposed == 1))
canon_n("inclusive.controls_unexposed", sum(inc$y == 0 & inc$exposed == 0))
canon_n("inclusive.controls_later_cases", sum(inc$y == 0 & inc$case == 1))
fit_inc <- glm(y ~ exposed, family = binomial, data = inc)
# A person can appear as a case and as a control, so the variance is clustered on the person
# (HC0 with the G/(G-1) cluster factor, the same small-sample factor Stata's logistic uses)
v_inc <- vcovCL(fit_inc, cluster = ~id, type = "HC0", cadjust = TRUE)
canon_ratio("inclusive.or", coef(fit_inc), v_inc)
# 8. Concurrent (risk-set) controls, frozen draw shared with Stata: matched analysis by set
rs <- dat("W7_riskset.csv")
canon_n("riskset.n", nrow(rs))
canon_n("riskset.sets", length(unique(rs$set)))
canon_n("riskset.controls_exposed", sum(rs$y == 0 & rs$exposed == 1))
canon_n("riskset.controls_unexposed", sum(rs$y == 0 & rs$exposed == 0))
canon_n("riskset.controls_later_cases", sum(rs$y == 0 & rs$case == 1))
fit_rs <- clogit(y ~ exposed + strata(set), data = rs)
canon_ratio("riskset.hr", coef(fit_rs), vcov(fit_rs))
# For contrast only: pooling the risk sets into one unmatched 2 by 2 table ignores the sampling design
fit_crude <- glm(y ~ exposed, family = binomial, data = rs)
canon("riskset.crude_or", exp(coef(fit_crude)[["exposed"]])) > ex <- dat("W7_exclusive.csv")
> canon_n("exclusive.n", nrow(ex))
CANON w7.exclusive.n 4680
> canon_n("exclusive.controls_exposed", sum(ex$y ==
+ 0 & ex$exposed == 1))
CANON w7.exclusive.controls_exposed 539
> canon_n("exclusive.controls_unexposed", sum(ex$y ==
+ 0 & ex$exposed == 0))
CANON w7.exclusive.controls_unexposed 2581
> tab_ex <- table(exposed = factor(ex$exposed, 1:0),
+ y = factor(ex$y, 1:0))
> print(tab_ex)
y
exposed 1 0
1 600 539
0 960 2581
> fit_ex <- glm(y ~ exposed, family = binomial, data = ex)
> canon_ratio("exclusive.or", coef(fit_ex), vcov(fit_ex))
CANON w7.exclusive.or 2.9928
CANON w7.exclusive.or_lb 2.6073
CANON w7.exclusive.or_ub 3.4354
> fx <- fisher.test(tab_ex)
> canon("exclusive.or_cmle.r", unname(fx$estimate))
CANON w7.exclusive.or_cmle.r 2.9920
> canon("exclusive.or_lb_exact.r", fx$conf.int[1])
CANON w7.exclusive.or_lb_exact.r 2.6007
> canon("exclusive.or_ub_exact.r", fx$conf.int[2])
CANON w7.exclusive.or_ub_exact.r 3.4436
> inc <- dat("W7_inclusive.csv")
> canon_n("inclusive.n", nrow(inc))
CANON w7.inclusive.n 4680
> canon_n("inclusive.controls_exposed", sum(inc$y ==
+ 0 & inc$exposed == 1))
CANON w7.inclusive.controls_exposed 619
> canon_n("inclusive.controls_unexposed", sum(inc$y ==
+ 0 & inc$exposed == 0))
CANON w7.inclusive.controls_unexposed 2501
> canon_n("inclusive.controls_later_cases", sum(inc$y ==
+ 0 & inc$case == 1))
CANON w7.inclusive.controls_later_cases 493
> fit_inc <- glm(y ~ exposed, family = binomial, data = inc)
> v_inc <- vcovCL(fit_inc, cluster = ~id, type = "HC0",
+ cadjust = TRUE)
> canon_ratio("inclusive.or", coef(fit_inc), v_inc)
CANON w7.inclusive.or 2.5252
CANON w7.inclusive.or_lb 2.2515
CANON w7.inclusive.or_ub 2.8322
> rs <- dat("W7_riskset.csv")
> canon_n("riskset.n", nrow(rs))
CANON w7.riskset.n 4680
> canon_n("riskset.sets", length(unique(rs$set)))
CANON w7.riskset.sets 1560
> canon_n("riskset.controls_exposed", sum(rs$y == 0 &
+ rs$exposed == 1))
CANON w7.riskset.controls_exposed 609
> canon_n("riskset.controls_unexposed", sum(rs$y ==
+ 0 & rs$exposed == 0))
CANON w7.riskset.controls_unexposed 2511
> canon_n("riskset.controls_later_cases", sum(rs$y ==
+ 0 & rs$case == 1))
CANON w7.riskset.controls_later_cases 244
> fit_rs <- clogit(y ~ exposed + strata(set), data = rs)
> canon_ratio("riskset.hr", coef(fit_rs), vcov(fit_rs))
CANON w7.riskset.hr 2.6139
CANON w7.riskset.hr_lb 2.2735
CANON w7.riskset.hr_ub 3.0051
> fit_crude <- glm(y ~ exposed, family = binomial, data = rs)
> canon("riskset.crude_or", exp(coef(fit_crude)[["exposed"]]))
CANON w7.riskset.crude_or 2.5770 ความเข้าใจผิดที่พบบ่อยและวิธีแก้
-
ตัดกลุ่มควบคุมแบบ inclusive ที่ภายหลังกลายเป็นผู้ป่วยทิ้ง
ในตัวอย่างแบบ inclusive ที่บันทึกไว้ 493 จาก 3,120 คนของกลุ่มควบคุมกลายเป็นผู้ป่วยภายหลัง การตัดคนเหล่านี้ออกเหลือเฉพาะกลุ่มควบคุมที่ยังไม่เป็นโรค ซึ่งก็คือกลุ่มแบบ exclusive และดันผลคูณไขว้จากอัตราส่วนความเสี่ยงไปทาง odds ratio
วิธีแก้: เก็บกลุ่มควบคุมตอนเริ่มต้นไว้ทุกคน แม้คนเดียวกันจะปรากฏเป็นผู้ป่วยด้วย และจัดกลุ่มค่าคลาดเคลื่อนมาตรฐานตามแต่ละคน
-
วิเคราะห์ตัวอย่างกลุ่มเสี่ยงด้วยแบบจำลองโลจิสติกรวมที่ไม่จับคู่แบบเดียว
การรวมทุกชุดเป็นตาราง 2 คูณ 2 ตารางเดียวทิ้งข้อมูลว่ากลุ่มควบคุมแต่ละคนถูกสุ่มเมื่อใด ในการสุ่มรอบนี้ ค่ารวม 2.58 บังเอิญตกใกล้ค่าแบบจับคู่ 2.61 นั่นเป็นลักษณะของ cohort นี้และการสุ่มรอบนี้ ไม่ใช่สิ่งที่รับประกันได้ odds ratio รวมที่ไม่จับคู่ไม่ใช่ค่าที่ถูกต้องโดยทั่วไปสำหรับตัวอย่างกลุ่มเสี่ยง
วิธีแก้: วิเคราะห์ชุดจับคู่ด้วยการถดถอยโลจิสติกแบบมีเงื่อนไข คือ
clogitทั้งใน Stata และ R -
ใช้
ccกับข้อมูลจับคู่ccสร้างตาราง 2 คูณ 2 ที่ไม่จับคู่ตารางเดียว เมื่อใช้กับตัวอย่างกลุ่มเสี่ยง มันให้ผลคูณไขว้รวมค่าเดียวกับแบบจำลองโลจิสติกที่ไม่จับคู่ และมีปัญหาเดียวกันวิธีแก้: ใช้
ccหรือlogisticสำหรับตัวอย่างแบบ exclusive ที่ไม่จับคู่ ใช้logisticพร้อมvce(cluster id)สำหรับตัวอย่างแบบ inclusive (ccจัดกลุ่มไม่ได้) และใช้clogitสำหรับชุดจับคู่ -
วิเคราะห์ผลลัพธ์ของ
sttoccด้วยตัวแปรเดิมของ cohortหลัง
sttoccผู้ป่วยของแต่ละชุดจับคู่ถูกทำเครื่องหมายด้วย_caseและชุดด้วย_setตัวบ่งชี้ผู้ป่วยเดิมของ cohort จะทำเครื่องหมายกลุ่มควบคุมที่ภายหลังเป็นผู้ป่วย และแบบจำลองที่ไม่มีgroup(_set)จะมองข้ามการจับคู่วิธีแก้: ต่อจาก
sttocc exposed, number(2) nodotsด้วยclogit _case exposed, group(_set) or nolog -
"odds ratio แบบ case-control ทุกค่าประมาณ odds ratio"
ผู้ป่วยชุดเดียวกัน 1,560 รายให้ odds ratio 2.99 อัตราส่วนความเสี่ยง (ผลคูณไขว้) 2.53 และอัตราส่วนฮาซาร์ด 2.61 เป้าหมายของมันคือ 3.14, 2.50 และ 2.78 Stata ติดป้ายค่าประมาณเหล่านี้ทุกค่าว่า odds ratio และใน R ค่าเอกซ์โพเนนเชียลของสัมประสิทธิ์จาก logistic หรือ clogit ก็อ่านกันตามธรรมเนียมว่าเป็น odds ratio ไม่ว่าการออกแบบจะเป็นแบบใด งานวิจัยที่ตีพิมพ์จำนวนมากไม่ระบุกฎการสุ่ม [1]
วิธีแก้: สิ่งที่ผลคูณไขว้ประมาณขึ้นกับวิธีสุ่มกลุ่มควบคุม คือ odds ratio สำหรับกลุ่มควบคุมจากผู้ที่ไม่เป็นโรคเมื่อสิ้นสุดการติดตาม อัตราส่วนความเสี่ยงสำหรับกลุ่มควบคุมจาก cohort ตอนเริ่มต้น และอัตราส่วนอัตราอุบัติการณ์หรืออัตราส่วนฮาซาร์ดสำหรับกลุ่มควบคุมจากกลุ่มเสี่ยง ณ เวลาเกิดเหตุการณ์ของผู้ป่วยแต่ละราย
สิ่งที่ควรทำในการวิเคราะห์ของคุณเอง
- เขียนกฎการสุ่มกลุ่มควบคุมไว้ในโปรโตคอล คือกลุ่มที่สุ่มจาก เวลาที่สุ่ม จำนวนกลุ่มควบคุมต่อผู้ป่วยหนึ่งราย กลุ่มควบคุมจะกลายเป็นผู้ป่วยภายหลังได้หรือไม่ และทุกคนถูกติดตามเป็นระยะเวลาเท่ากันโดยไม่มีผู้หลุดจากการติดตามหรือไม่
- ระบุมาตรวัดที่กฎนั้นและการติดตามนั้นมุ่งไปที่ (odds ratio อัตราส่วนความเสี่ยง หรืออัตราส่วนฮาซาร์ด) ก่อนประมาณแบบจำลองใด ๆ
- สำหรับตัวอย่างกลุ่มเสี่ยง เก็บชุดจับคู่ไว้และวิเคราะห์ด้วยการถดถอยโลจิสติกแบบมีเงื่อนไข ตารางรวมใช้บรรยายตัวอย่างได้ แต่มักไม่ใช่ค่าประมาณที่ควรรายงาน
- สำหรับตัวอย่างแบบ inclusive เก็บกลุ่มควบคุมที่ภายหลังกลายเป็นผู้ป่วยไว้ และมักจัดกลุ่มค่าคลาดเคลื่อนมาตรฐานตามแต่ละคน
- รายงานค่าประมาณแต่ละค่าภายใต้ชื่อของมาตรวัดที่มันประมาณ ไม่ว่าโปรแกรมจะพิมพ์ป้ายว่าอะไร และระบุกฎการสุ่มไว้ข้างค่านั้น
อภิธานศัพท์
- case-control study
- การศึกษาที่วัดสิ่งสัมผัสในผู้ป่วยและในตัวอย่างของกลุ่มควบคุมที่สุ่มมาเป็นตัวแทนประชากรต้นทางของผู้ป่วย แทนการวัดในทุกคน
- control series
- กลุ่มควบคุมที่การศึกษาสุ่มมา รวมกับกฎที่ใช้สุ่ม
- closed cohort
- cohort ที่ทุกคนเข้ามาตอนเริ่มต้นและถูกติดตามเป็นระยะเวลาคงที่เท่ากันโดยไม่มีผู้หลุดจากการติดตาม
- cumulative (exclusive) sampling (การเลือกกลุ่มควบคุมแบบสะสม)
- การสุ่มกลุ่มควบคุมจากผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตาม ใน closed cohort ผลคูณไขว้ประมาณ odds ratio
- case-base (inclusive) sampling (การเลือกกลุ่มควบคุมจากฐานประชากรทั้งหมด)
- การสุ่มกลุ่มควบคุมจากทั้ง cohort ตอนเริ่มต้น รวมผู้ที่ภายหลังกลายเป็นผู้ป่วย ใน closed cohort ผลคูณไขว้ประมาณอัตราส่วนความเสี่ยง
- risk-set (concurrent) sampling (การสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยง)
- การสุ่มกลุ่มควบคุม ณ เวลาเกิดเหตุการณ์ของผู้ป่วยแต่ละรายจากผู้ที่ยังมีความเสี่ยงอยู่ เมื่อวิเคราะห์เป็นชุดจับคู่ จะประมาณอัตราส่วนฮาซาร์ด
- risk set
- ทุกคนที่ยังอยู่ในการติดตามและยังไม่เป็นโรค ณ เวลาที่ผู้ป่วยรายหนึ่งเกิดโรค
- nested case-control design (การศึกษาแบบ nested case-control)
- การศึกษาแบบ case-control ที่สุ่มจาก cohort ที่กำหนดไว้ ตามการใช้ทั่วไปคือสุ่มกลุ่มควบคุมจากกลุ่มเสี่ยงของผู้ป่วยแต่ละราย
- cross-product
- ผลคูณของผู้ป่วยที่ได้รับสิ่งสัมผัสกับกลุ่มควบคุมที่ไม่ได้รับ หารด้วยผู้ป่วยที่ไม่ได้รับคูณกลุ่มควบคุมที่ได้รับ เป็นตัวเลขที่ cc และ logistic พิมพ์เป็น odds ratio
- conditional logistic regression
- แบบจำลองโลจิสติกที่มีจุดตัดแกนของตัวเองให้แต่ละชุดจับคู่ ประมาณเพื่อให้ผู้ป่วยแต่ละรายถูกเปรียบเทียบกับกลุ่มควบคุมของตัวเองเท่านั้น
- matched set
- ผู้ป่วยหนึ่งรายพร้อมกลุ่มควบคุมที่สุ่มให้รายนั้น
- hazard ratio (อัตราส่วนฮาซาร์ด)
- อัตราส่วนของอัตราการเกิดเหตุการณ์ ณ ขณะใดขณะหนึ่งระหว่างสองกลุ่ม ในผู้ที่ยังไม่เกิดเหตุการณ์
เอกสารอ้างอิง
- Knol MJ, Vandenbroucke JP, Scott P, Egger M. What do case-control studies estimate? Survey of methods and assumptions in published case-control research. Am J Epidemiol. 2008;168(9):1073-1081. doi:10.1093/aje/kwn217 https://doi.org/10.1093/aje/kwn217
- Rodrigues L, Kirkwood BR. Case-control designs in the study of common diseases: updates on the demise of the rare disease assumption and the choice of sampling scheme for controls. Int J Epidemiol. 1990;19(1):205-213. doi:10.1093/ije/19.1.205 https://doi.org/10.1093/ije/19.1.205
- Greenland S, Thomas DC. On the need for the rare disease assumption in case-control studies. Am J Epidemiol. 1982;116(3):547-553. doi:10.1093/oxfordjournals.aje.a113439 https://doi.org/10.1093/oxfordjournals.aje.a113439
- Pearce N. What does the odds ratio estimate in a case-control study? Int J Epidemiol. 1993;22(6):1189-1192. doi:10.1093/ije/22.6.1189 https://doi.org/10.1093/ije/22.6.1189
- Prentice RL. A case-cohort design for epidemiologic cohort studies and disease prevention trials. Biometrika. 1986;73(1):1-11. doi:10.1093/biomet/73.1.1 https://doi.org/10.1093/biomet/73.1.1
- Barlow WE, Ichikawa L, Rosner D, Izumi S. Analysis of case-cohort designs. J Clin Epidemiol. 1999;52(12):1165-1172. doi:10.1016/S0895-4356(99)00102-X https://doi.org/10.1016/S0895-4356(99)00102-X
- Langholz B, Goldstein L. Risk set sampling in epidemiologic cohort studies. Stat Sci. 1996;11(1):35-53. doi:10.1214/ss/1032209663 https://doi.org/10.1214/ss/1032209663
ประเด็นสำคัญ
- ผู้ป่วยเหมือนกันในทุกการออกแบบ กลุ่มควบคุมเป็นตัวกำหนดว่าผลคูณไขว้ประมาณอะไร เพราะกลุ่มควบคุมเป็นตัวแทนของตัวหาร
- ใน closed cohort กลุ่มควบคุมจากผู้ที่ยังไม่เป็นโรคเมื่อสิ้นสุดการติดตามประมาณ odds ratio 5 ปี ซึ่งในที่นี้คือ 3.14 และเข้าใกล้อัตราส่วนความเสี่ยงเฉพาะเมื่อโรคพบน้อย
- ใน closed cohort กลุ่มควบคุมจากทั้ง cohort ตอนเริ่มต้น รวมผู้ที่ภายหลังกลายเป็นผู้ป่วย ประมาณอัตราส่วนความเสี่ยง 5 ปี ซึ่งในที่นี้คือ 2.50 และมักจัดกลุ่มค่าคลาดเคลื่อนมาตรฐานตามแต่ละคน
- กลุ่มควบคุมจากกลุ่มเสี่ยงของผู้ป่วยแต่ละราย เมื่อวิเคราะห์เป็นชุดจับคู่ด้วยการถดถอยโลจิสติกแบบมีเงื่อนไข ประมาณอัตราส่วนฮาซาร์ด ซึ่งในที่นี้คือ 2.78 โดยไม่ต้องเป็น closed cohort และ odds ratio รวมที่ไม่จับคู่ไม่ใช่ค่าที่ถูกต้องโดยทั่วไป
- รายงานค่าประมาณแต่ละค่าภายใต้ชื่อของมาตรวัดที่มันประมาณ ไม่ว่า Stata หรือ R จะพิมพ์ป้ายว่าอะไร
อ่านต่อในวิกิ: [[nested-case-control-design-th]]