Ứng dụng mô hình hồi qui logistic

Tuan V. Nguyen Senior Principal Research Fellow, Garvan Institute of Medical Research Professor, UNSW School of Public Health and Community Medicine Professor of Predictive Medicine, University of Technology Sydney Adj. Professor of Epidemiology and Biostatistics, School of Medicine Sydney, University of Notre Dame Australia

Phân tích dữ liệu và ứng dụng | Đại học Dược Hà Nội | 12/6 to 17/6/2019

© Tuan V. Nguyen

Ba ứng dụng hình hồi qui logistic

• Đánh giá mối liên quan, ảnh hưởng (association)

• Hiệu chỉnh (adjustment) cho yếu tố nhiễu

• Tiên lượng (prediction)

Ứng dụng 1: Đánh giá mối liên quan, ảnh hưởng

Cân nặng của trẻ sơ sinh và mẹ hút thuốc lá

• Câu hỏi nghiên cứu: có mối liên quan giữa mẹ hút thuốc lá và trọng lượng

sơ sinh?

• Nghiên cứu cắt ngang trên 189 bà mẹ và trẻ sơ sinh

Biến outcome: low, biến tiên lượng: smoke

bw = read.csv("~/Dropbox/_Conferences and Workshops/Dai hoc Duoc 6- 2019/Datasets/birthwt.csv")

head(bw)

id low age lwt race smoke ptl ht ui ftv bwt

1 85 0 19 182 2 0 0 0 1 0 2523

2 86 0 33 155 3 0 0 0 0 3 2551

3 87 0 20 105 1 1 0 0 0 1 2557

Hiển thị mối liên quan bằng biểu đồ

Summary:

bw$smoke = as.factor(bw$smoke)

bs$low = as.factor(bw$low)

library(DescTools)

Desc(bw$smoke ~ bw$low)

estimate lwr.ci upr.ci' odds ratio 2.022 1.081 3.783 rel. risk (col1) 1.258 1.013 1.561 rel. risk (col2) 0.622 0.409 0.945 bw$low 0 1 Sum bw$smoke

0 freq 86 29

115

p.row 74.8% 25.2%

.

p.col 66.2% 49.2%

.

1 freq 44 30

74

p.col 33.8% 50.8%

.

Sum freq 130 59

189

p.row . .

.

p.col . .

.

Mô hình hồi qui logistic

Mô hình liên quan giữa smoke và low như sau: gọi P là xác suất trọng lượng thấp

(cid:230) (cid:246)

(cid:247) = a + b smoke

log

(cid:231)

P 1- P

Ł ł

# Tri n khai b ng R ể ằ

m = glm(low ~ smoke, family=binomial, data=bw) summary(m)

# tính odds ratio library(epiDisplay) logistic.display(m)

Kết quả phân tích

Mô hình liên quan giữa smoke và low như sau: gọi P là xác suất trọng lượng thấp

(cid:230) (cid:246)

(cid:247) = a + b smoke

log

(cid:231)

P 1- P

Ł ł

Coefficients:

Estimate Std. Error z value Pr(>|z|)

(Intercept) -1.0871 0.2147 -5.062 4.14e-07 ***

smoke1 0.7041 0.3196 2.203 0.0276 *

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(cid:230) (cid:246)

log

(cid:247) = - 1.087 + 0.704 · smoke

(cid:231)

P 1- P

Ł ł

Kết quả phân tích

> logistic.display(m)

Logistic regression predicting low : 1 vs 0 OR(95%CI) P(Wald's test) P(LR-test) smoke: 1 vs 0 2.02 (1.08,3.78) 0.028 0.027 Log-likelihood = -114.9023 No. of observations = 189 AIC value = 233.8046

Diễn giải: có mối liên quan giữa bà mẹ hút thuốc lá và sinh con thiếu cân. Tính trung bình, những bà mẹ hút thuốc lá có odds sinh con thiếu cân cao 2 lần (KTC 95%: 1.1 đến 3.8) so với những bà mẹ không hút thuốc lá.

Ứng dụng 2: Hiệu chỉnh cho yếu tố nhiễu

Yếu tố nhiễu (confounding factors)

• Câu hỏi nghiên cứu: có mối liên quan giữa bà mẹ hút thuốc lá và

sinh con thiếu cân ?

• Yếu tố nhiễu: có liên quan đến yếu tố nguy cơ (risk factor) và

outcome

• Mối liên quan giữa hút thuốc lá và thiếu cân có thể do yếu tố

chủng tộc và độ tuổi?

• Câu hỏi mới: mối liên quan giữa hút thuốc lá và thiếu cân có độc

lập với yếu tố chủng tộc và độ tuổi?

Yếu tố nhiễu: chủng tộc

> Desc(factor(bw$race) ~ bw$low)

> Desc(factor(bw$race) ~ bw$smoke)

bw$low 0 1 Sum factor(bw$race)

bw$smoke 0 1 Sum factor(bw$race)

1 freq 73 23 96 p.row 76.0% 24.0% . p.col 56.2% 39.0% .

1 freq 44 52 96 p.row 45.8% 54.2% . p.col 38.3% 70.3% .

2 freq 15 11 26 p.row 57.7% 42.3% . p.col 11.5% 18.6% .

2 freq 16 10 26 p.row 61.5% 38.5% . p.col 13.9% 13.5% .

3 freq 42 25 67 p.row 62.7% 37.3% . p.col 32.3% 42.4% .

3 freq 55 12 67 p.row 82.1% 17.9% . p.col 47.8% 16.2% .

Mô hình hiệu chỉnh cho yếu tố nhiễu

• Câu hỏi mới: mối liên quan giữa hút thuốc lá và thiếu cân có độc

lập với yếu tố chủng tộc và độ tuổi?

• Giải pháp: phân tích mô hình hồi qui logistic

• Mô hình mới:

b

(cid:230) (cid:246)

(cid:247) = a

log

1smoke + b

2race + b

3age

· (cid:231)

P 1- P

Ł ł

m = glm(low ~ smoke + factor(race) + age, family = binomial, data=bw)

logistic.display(m)

Kết quả phân tích

> m = glm(low ~ smoke + factor(race) + age, family = binomial, data=bw)

> logistic.display(m)

OR lower95ci upper95ci Pr(>|Z|)

smoke1 3.0058203 1.4499820 6.231081 0.003087348

factor(race)2 2.7494834 1.0453178 7.231924 0.040385185

factor(race)3 2.8769483 1.2983141 6.375061 0.009239817

age 0.9657186 0.9045206 1.031057 0.296337134

Diễn giải: Sau khi hiệu chỉnh cho yếu tố chủng tộc và tuổi, những bà mẹ hút thuốc lá có odds sinh con thiếu cân tăng gấp 3 lần (KTC 95%: 1.45 đến 6.23) so với những bà mẹ không hút thuốc lá. Những bà mẹ da đen (OR 2.7; KTC95: 1.05 – 7.23) và Hispanics (OR 2.88; 95% CI: 1.30 – 6.37) có odds sinh con thiếu cân cao hơn những bà mẹ da trắng. Độ tuổi có mẹ không có liên quan có ý nghĩa thống kê với sinh con thiếu cân.

Yếu tố nào có liên quan đến sinh con thiếu cân?

• Có 8 yếu tố nguy cơ trong nghiên cứu: age, race, lwt,

smoke, ptl, ht, ui, ftv

• Câu hỏi mới: yếu tố nào có liên quan đến sinh con thiếu cân?

• Có ít nhất 2^8 = 256 mô hình khả dĩ !

• Giải pháp: tìm mô hình tối ưu (model selection, feature selection)

• Phương pháp: Bayesian model averaging (BMA)

Tiêu chí cho một mô hình đa biến tối ưu

• Đơn giản (simplicity): ít biến số

• Đầy đủ (adequacy): mô tả dữ liệu một cách thỏa đáng

• Có ý nghĩa thực tế (practicality): phải được yểm trợ bằng lí

thuyết hay có ý nghĩa sinh học

Thách thức: Tìm mô hình với ít biến số nhất nhưng giải

thích nhiều dữ liệu nhất (nguyên tắc parsimony)

Tiêu chuẩn thống kê: AIC

• AIC = Akaike Information Criterion.

AIC = Deviance + 2x(thông số trong mô hình)

• Mô hình tối ưu nhất có AIC thấp nhất. AIC “phạt” những mô hình có

nhiều thông số.

• Nhiều thuật toán thống kê đã được phát triển đi tìm một mô hình với

k thông số sao cho AIC thấp nhất.

• Tìm một mô hình tối ưu là một thách thức rất lớn!

Thuật toán tìm mô hình tối ưu

• Stepwise (forward) algorithm

– Bắt đầu với mô hình 1 thông số

– Từng bước thêm vào 1 thông số cho đến khi deviance không còn giảm

nữa.

• Backward algorithm

– Bắt đầu với mô hình gồm k thông số.

– Loại bỏ từng bước những thông số không có ý nghĩa thống kê

cho đến khi deviance thấp nhất.

• AIC based algorithm: Tìm tất cả các mô hình với AIC thấp nhất

• Bayesian model average (BMA). Tìm tổ các thông số sao cho BIC

thấp nhất

Áp dụng BMA vào nghiên cứu thiếu cân

# Gọi package "BMA". Nếu chưa có, cần phải cài đặt trước library(BMA)

# Chuẩn bị dữ liệu, biến y và biến x, loại bỏ biến id (cột số 1) và y (cột 2) và bwt (cột xvars = bw[, -c(1, 2, 11)]

yvar = bw[, 2]

# Tìm mô hình tối ưu bằng BMA

m = bic.glm(xvars, yvar, strict=F, OR=20, glm.family="binomial")

summary(m)

imageplot.bma(m)

Kết quả phân tích BMA

> summary(m)

84 models were selected

Best 5 models (cumulative posterior probability = 0.2531 ):

p!=0 EV SD model 1 model 2 model 3 model 4 model 5

Intercept 100 -0.390128 1.575728 1.45068 1.09291 -2.32488 -0.35754 1.06795

age 10.4 -0.004815 0.018070 . . . . .

lwt 54.8 -0.008473 0.009253 -0.01865 -0.01707 . -0.01535 -0.01692

race 44.3 0.212462 0.280368 . . 0.55898 0.48955 .

smoke 52.1 0.484523 0.552668 . . 1.11668 1.08002 .

ptl 41.2 0.291512 0.410590 . 0.72560 . . .

ht 59.7 1.011382 0.999519 1.85551 1.85604 . 1.74427 1.96157

ui 30.0 0.263111 0.470489 . . . . 0.93000

ftv 2.0 -0.001015 0.024588 . . . . .

nVar 2 3 2 4 3

BIC -753.82285 -753.75940 -753.62525 -753.44086 -753.11035

post prob 0.058 0.056 0.052 0.048 0.040

BMA phân tích tất cả 84 mô hình, nhưng chỉ báo cáo 5 mô hình tốt nhất

Diễn giải kết quả BMA

p!=0 EV SD

Intercept 100 -0.390128 1.575728

age 10.4 -0.004815 0.018070

lwt 54.8 -0.008473 0.009253

race 44.3 0.212462 0.280368

smoke 52.1 0.484523 0.552668

ptl 41.2 0.291512 0.410590

ht 59.7 1.011382 0.999519

ui 30.0 0.263111 0.470489

ftv 2.0 -0.001015 0.024588

xác suất biến số có liên quan với y

Cột 1 các biến số tiên lượng p!=0 EV expected value: giá trị trung bình của hệ số hồi qui logistic (𝛃) SD standard deviation: độ lệch chuẩn của 𝛃

p!=0 EV SD model 1 model 2 model 3 model 4 model 5 Intercept 100 -0.390128 1.575728 1.45068 1.09291 -2.32488 -0.35754 1.06795 age 10.4 -0.004815 0.018070 . . . . . lwt 54.8 -0.008473 0.009253 -0.01865 -0.01707 . -0.01535 -0.01692 race 44.3 0.212462 0.280368 . . 0.55898 0.48955 . smoke 52.1 0.484523 0.552668 . . 1.11668 1.08002 . ptl 41.2 0.291512 0.410590 . 0.72560 . . . ht 59.7 1.011382 0.999519 1.85551 1.85604 . 1.74427 1.96157 ui 30.0 0.263111 0.470489 . . . . 0.93000 ftv 2.0 -0.001015 0.024588 . . . . . nVar 2 3 2 4 3 BIC -753.82285 -753.75940 -753.62525 -753.44086 -753.11035 post prob 0.058 0.056 0.052 0.048 0.040

số biến tiên lượng của một mô hình

Bayesian information criterion (càng thấp càng tốt)

nVar BIC post prob Model 1 Xác suất hậu định (posterior probability). bao gồm biến lwt (cân nặng của mẹ) và ht (cao huyết áp). Mô hình này có xác

suất hậu định ~6% và được xem là tốt nhất trong các mô hình

Model 2 Model 3 BMA 'đề nghị' biến lwt, ptl, và ht được nhận dạng với các biến race và smoke

Kết quả BMA qua biểu đồ

Models selected by BMA

age

lwt

race

smoke

ptl

ht

ui

ftv

1

2

3

4

5

6

7 8

10

12 14

17

20 23

27 31

36

42 48 56 65 76

Model #

Phương pháp BMA

• Không đưa ra một mô hình duy nhất, mà đề nghị hơn 1 mô hình tốt

nhất

• Mô hình nào tốt nhất còn tuỳ thuộc vào kiến thức chuyên ngành

• BMA được 'chứng minh' là phương pháp tốt nhất khi chọn mô hình

tối ưu

• Dựa vào kết quả, chúng ta chọn mô hình nào để đánh giá mối liên quan giữa các yếu tố nguy cơ và nguy cơ sinh con thiếu cân?

Chọn mô hình" tối ưu"

• Xây dựng mô hình đa biến vừa là một nghệ thuật, vừa là một khoa

học

• Một mô hình đẹp là mô hình mô tả sát với thực tế.

• Một mô hình phản ánh 100% thực tế thì đó không còn là “mô hình”

nữa, hay quá phức tạp không thể ứng dụng được.

• Ngược lại một mô hình chỉ mô tả thực tế khoảng 1% thì cũng không

thể sử dụng được.

“Mô hình nào cũng sai so với thực tế, nhưng trong số các mô hình sai đó, có một vài mô hình có ích (George Box)”.