Survival – Phân tích sự tồn tại _ Phân tích sống còn_ Phân tích sống sót_ Phân tích sự kiện
Phân tích sống sót (PTSS) hoặc phân tích sự kiện khi nhà nghiên cứu muốn tìm hiểu ảnh hưởng đến các biến kết cục (biến phụ thuộc) mang tính thời gian.
Khi so sánh 2 phương pháp điều trị cho các bệnh có tần số tử vong cao như bệnh AIDS, các bệnh ung thư… Nếu mô hình phân tích như phân tích hồi qui logistic, chỉ để ý đến biến kết cục (sống/chết hoặc khỏi bệnh/không khỏi bệnh) mà không quan tâm đến yếu tố thời gian thì đôi không tìm thấy sự khác biệt giữa 2 phương pháp điều trị vì tỉ lệ tử vong gần như nhau, nhưng thời gian dẫn đến tử vong ở 2 nhóm có thể khác nhau. Vì vậy chúng ta phải sử dụng mô hình PTSS thì mới thấy sự khác biệt này. Như vậy mô hình nghiên cứu mô tả kết cục là biến nhị phân (sống/chết) tuy quan trọng nhưng không chính xác.
Mô hình chính để thể hiện mối liên hệ giữa thời gian dẫn đến bệnh (hay không bệnh) và các yếu tố nguy cơ (risk factors) là mô hình có tên là “survival analysis – Phân tích sự sông sót”. Cụm từ “survival analysis” xuất phát từ nghiên cứu trong bảo hiểm, và giới nghiên cứu y khoa từ đó dùng cụm từ cho bộ môn của mình. Nhưng như nói trên, sống/chết không phải là biến duy nhất, vì trong thực tế chúng ta cũng có những biến cố có bệnh hay không bệnh, xảy ra hay không xảy ra. Ngoài ra, trong các bộ môn kĩ thuật, người ta dùng một cụm từ khác “reliability analysis – Phân tích độ tin cậy”, dể chỉ cho khái niệm survival analysis. Tuy nhiên trong đề tài này chúng tôi sẽ dùng cụm từ Phân tích sống còn.
Mô hình hóa dữ liệu thời gian đến sự kiện là một chủ đề quan trọng với nhiều ứng dụng trong các lĩnh vực khác nhau. Tập hợp các phương pháp để phân tích dữ liệu đó được gọi là phân tích sự tồn tại, phân tích lịch sử sự kiện hoặc phân tích thời lượng. Phân tích sống còn được áp dụng rộng rãi vì định nghĩa về một ‘sự kiện’ có thể rất đa dạng và các ví dụ bao gồm chết và sống, thành công và phá sản, tiếp tục hoặc hủy bỏ. Do đó, các lĩnh vực ứng dụng bao gồm từ y học và xã hội học đến tiếp thị và kinh tế. Trong bài tiểu luận này, chúng ta xem xét các vấn đề cơ bản về lý thuyết của phân tích sống còn bao gồm các công cụ ước tính cho các chức năng sinh tồn và rủi ro và cách thực hành trên R dựa trên package survival.
Trong thời đại ngày nay - Công nghiệp 4.0 là thời đại cho sự phát triển kinh tế xã hội và công nghệ kĩ thuật, việc ứng dụng của phân tích sống sót trong tài chính doanh nghiệp trở thành một vấn đề có ý nghĩa thời sự, đặc biệt trong bối cảnh sự phát triển đa ngành nghề một cách sâu rộng và nhanh chóng như hiện nay.
Thấy được tầm quan trọng cũng như tiềm năng khi phân tích sự tồn tại của package survival trong thời điểm hiện tại, chúng tôi mong muốn có thể kết hợp các ngôn ngữ lập trình để giải quyết vấn đề trong kinh tế xã hội hiện và mô hình này sẽ giúp ích cho việc phân tích các nguy cơ, sự sống sót có thể xảy ra từ đó có thể dễ dàng giải quyết được các vấn đề quan trọng.
Một ưu điểm cùa Phân tích sống còn là xử lý được các trường hợp đối tượng nghiên cứu bỏ cuộc giữa chừng (như mất dấu theo dõi, ngưng điều trị do tác dụng phụ của thuốc hoặc tử vong do bệnh lý khác…).
Trong mô hình phân tích này các đối tượng còn sống kể cả đối tượng bỏ cuộc được gọi là censored hoặc sự kiện chưa xảy ra. Các đối tượng tử vong hoặc hết sống được gọi là events hoặc sự kiện đã kết thúc. Mục tiêu suy luận đối với phân tích sống còn là khoảng thời gian giữa thời điểm bắt đầu và thời điểm sự kiện xảy ra. Trong nghiên cứu y học hiện nay, phân tích sống còn được sử dụng rộng rãi trong các nghiên cứu lâm sàng để đánh giá hiệu quả của một phương phương điều trị hoặc để đánh giá tình trạng sống còn của một số các biện pháp điều trị ung thư. Trong kinh tế xã hội, phân tích sống còn được dùng để phân tích rủi ro/ nguy cơ để đánh giá sự tồn tại, của một doanh nghiệp nào đó.
Đối tượng nghiên cứu
Phân tích sống còn tập trung mô tả cho một cá thể hay một nhóm cá thể nhất định. Một điểm xác định của một sự kiện được gọi là failure (như là xuất hiện bệnh, chữa khỏi bệnh, tử vong, tái phát sau khi đáp ứng với điều trị…) mà xảy ra sau một khoảng thời gian được gọi là failure time (thời gian dẫn đến sự kiện) (hoặc follow-up time (thời gian theo dõi) trong nghiên cứu thuần tập/nghiên cứu dựa vào dân số) trong suốt thời gian các cá thể được quan sát. Để xác định thời gian dẫn đến sự kiện, chúng ta cần xác định thời điểm bắt đầu (có thể là ngày nhận vào, ngày chẩn đoán…).
Phạm vi nghiên cứu
Trong bài báo này, chúng tôi xem xét các vấn đề cơ bản về lý thuyết của phân tích tỷ lệ sống sót bao gồm các công cụ ước tính cho các chức năng sinh tồn và rủi ro. Chúng tôi thảo luận chi tiết về Mô hình mối nguy theo tỷ lệ Cox và cả các phương pháp để kiểm tra giả định về mối nguy theo tỷ lệ (PH). Hơn nữa, chúng tôi thảo luận về các mô hình Cox phân tầng cho các trường hợp khi giả định PH không đúng. Phần thảo luận của chúng tôi được bổ sung bằng một ví dụ hoạt động sử dụng ngôn ngữ lập trình thống kê R để cho phép ứng dụng thực tế của phương pháp này. Phép ước tính thường được dùng để Phân tích sống còn được gọi là ước tính KaplanMeier. Phép ước tính này giúp ta tính được xác suất sống sót tích lũy tại các mốc thời gian khác nhau (sẽ minh họa trong các ví dụ sau). Nếu muốn so sánh sự khác biệt giữa 2 nhóm điều trị, dùng kiểm định log-rank bằng cách so sánh 2 hàm xác suất tích lũy của 2 nhóm.
Package “survival” chứa các quy trình phân tích sự sống sót cốt lõi, bao gồm định nghĩa về các đối tượng Surv, đường cong Kaplan-Meier và Aalen-Johansen (đa trạng thái), mô hình Cox và mô hình thời gian thất bại tăng tốc theo tham số.
Phiên bản mới nhất: 3.5-5
Phiên bản R hỗ trợ: R(>3.5.0)
Nhóm tác giá: Terry M Therneau, Thomas Lumley, Atkinson Elizabeth, Crowson Cynthia
Người phụ trách chính: Terry M Therneau
Tài liệu hướng dẫn: https://cran.r-project.org/web/packages/survival/survival.pdf
Phân tích sống còn được sử dụng phổ biến trong sinh học, y học, kỹ thuật, tiếp thị, khoa học xã hội hoặc khoa học. Các tên tương đồng với “Phân tích sống còn” như là phân tích lịch sử sự kiện, lý thuyết độ tin cậy hoặc phân tích thời gian.
Có hai phương pháp đóng góp quan trọng cho sự phát triển của lĩnh vực này. Đầu tiên là phương pháp của Kaplan và Meier, người đã giới thiệu một công cụ “Ước tính xác suất sống sót”. Thứ hai là từ Cox, người đã giới thiệu một mô hình được gọi là “Mô hình rủi ro theo tỷ lệ Cox” (CPHM), là một mô hình hồi quy. Thật đặc biệt khi cả hai mô hình này đều được sử dụng nhiều cho đến nay và luôn là lực chọn của mọi nhà khoa học khi phân tích về dữ liệu sống còn (survival data).
Trong những năm qua, nhiều báo cáo đã được viết về khảo sát phân tích tỷ lệ sống còn. Tuy nhiên, do sự phức tạp của các phương pháp, đặc biệt là đối với dữ liệu đa biến, các cuộc thảo luận dễ dẫn đến nhầm lẫn. Ngoài ra, ở cấp độ giới thiệu có sẵn các đánh giá, đây hoàn toàn là lý thuyết hoặc sử dụng các ngôn ngữ lập trình (như stata hoặc SAS) được sử dụng chủ yếu trong các thống kê sinh học.
Ngược lại, đánh giá của một số người đã kết hợp trình bày lý thuyết với thực tế sử dụng ngôn ngữ lập trình thống kê R. R là ngôn ngữ được sử dụng rộng rãi cho các vấn đề chung trong khoa học dữ liệu vì nó kết hợp các tính năng từ các mô hình lập trình khác nhau. Chúng ta hướng tới một cách trình bày toàn diện và đủ rộng để bao gồm tất cả các chủ đề cần thiết cho phân tích sống còn đa biến nhưng đồng thời cũng dễ hiểu.
Để thực hiện tất cả điều này, cần bổ sung cho việc trình bày các phương pháp với đầy đủ thông tin cơ bản. Giải thích chi tiết việc kiểm chứng và xử lý dữ liệu về các sự kiện thời gian, bởi vì về cơ bản tất cả các phương pháp đều sử dụng kiểm chứng để rút ra các ước tính hiệu quả. Cuối cùng, chúng ta thêm một ví dụ đã hoạt động cho thấy cách thực hiện phân tích khả năng sống sót thực tế với R.
Trong phần tiếp theo, chúng ta cung cấp một dữ liệu để phân tích sự sống còn và ý nghĩa của việc phân tích chúng. Tiếp theo, chúng ta mô tả các đặc điểm chung của các hàm sinh tồn, ước lượng cho các hàm sinh tồn và so sánh hai đường cong sinh tồn. Sau đó, chúng ta giới thiệu Mô hình nguy cơ theo tỷ lệ Cox. Nếu giả định rủi ro theo tỷ lệ không đúng, người ta cần sử dụng mô hình Cox phân tầng. Cuối cùng, chúng ta trình bày một phần về phân tích tỷ lệ sống thực tế bằng R. Bài báo kết thúc với phần tóm tắt ngắn gọn và kết luận.
Dưới đây là một số hàm cơ bản trong package survival:
Surv(time, event): được sử dụng để tạo ra đối tượng “Surv” trong R. Đối tượng này được sử dụng để định nghĩa và lưu trữ dữ liệu sống sót và thời gian theo dõi tương ứng.
Survfit(formula, data): được sử dụng để tính toán tỷ lệ sống sót hoặc phân tích thời gian đến sự kiện cho các nhóm dữ liệu khác nhau.
Coxph(Surv (time, event) ~ I, data) : mô hình cox tiêu chuẩn, mô hình Cox được sử dụng để xác định yếu tố ảnh hưởng đến thời gian sống sót hoặc thời gian đến sự kiện.
Survreg(formula, data): Hàm này được sử dụng để phân tích các yếu tố ảnh hưởng đến thời gian sống sót hoặc thời gian đến sự kiện trong mô hình hồi quy.
Survexp(): Hàm này được sử dụng để tính toán tỷ lệ sống sót hoặc phân tích thời gian đến sự kiện cho các nhóm dữ liệu khác nhau trong mô hình phân phối mũ.
Để sử dụng package survival trên RStudio, bạn có thể thực hiện các bước sau:
Bước 1: Cài đặt package survival
Bạn có thể cài đặt package survival bằng cách chạy lệnh
sau trong RStudio:
install.packages("survival")
Nếu package đã được cài đặt trước đó, bạn có thể bỏ qua bước này.
Bước 2: Tải package survival
Sau khi cài đặt, bạn cần tải package survival bằng lệnh
sau:
library(survival)
Bước 3: Sử dụng các hàm của package survival
Sau khi đã tải package survival, bạn có thể sử dụng các
hàm của package để phân tích dữ liệu thời gian đến sự kiện. Các hàm này
bao gồm:
survfit(): tính toán các ước lượng sống sót bằng phương
pháp Kaplan-Meier hoặc phương pháp Nelson-Aalen.coxph(): phân tích mô hình hồi quy Cox để tìm các yếu
tố ảnh hưởng đến thời gian đến sự kiện.survreg(): phân tích mô hình hồi quy Cox để tìm các yếu
tố ảnh hưởng đến thang đo thời gian đến sự kiện không phải là dạng số
nguyên (ví dụ: thời gian sống sót được đo bằng giờ, ngày, tháng,…).cox.zph(): kiểm định giả thuyết về sự đồng nhất của các
hệ số hồi quy Cox theo thời gian.Trong suốt phần này, chúng tôi sẽ sử dụng dataset lung từ package survival làm dữ liệu phân tích.
Dataset lung là một trong những tập dữ liệu tiêu biểu của package survival được sử dụng trong nhiều ví dụ và hướng dẫn về phân tích dữ liệu. Dataset này được tạo ra từ nghiên cứu lâm sàng về ung thư phổi và bao gồm thông tin về 228 bệnh nhân ung thư phổi với các biến như tuổi, giới tính, số lượng góc nhìn của các khối u, phân loại tổn thương và tỷ lệ sống sót sau 1 năm.
Dataset lung được sử dụng để minh họa cho nhiều khía cạnh khác nhau của phân tích dữ liệu, bao gồm mô hình hồi quy tuyến tính, mô hình hồi quy logistic, phân tích phân nhóm và đặc biệt là phân tích sống sót.
Để truy cập dataset lung từ package survival trong R, ta sử dụng lệnh sau:
data(cancer)
Dữ liệu ung thư phổi có cấu trúc như sau:
str(lung)
## 'data.frame': 228 obs. of 10 variables:
## $ inst : num 3 3 3 5 1 12 7 11 1 7 ...
## $ time : num 306 455 1010 210 883 ...
## $ status : num 2 2 1 2 2 1 2 2 2 2 ...
## $ age : num 74 68 56 57 60 74 68 71 53 61 ...
## $ sex : num 1 1 1 1 1 1 2 2 1 1 ...
## $ ph.ecog : num 1 0 0 1 0 1 2 2 1 2 ...
## $ ph.karno : num 90 90 90 90 100 50 70 60 70 70 ...
## $ pat.karno: num 100 90 90 60 90 80 60 80 80 70 ...
## $ meal.cal : num 1175 1225 NA 1150 NA ...
## $ wt.loss : num NA 15 15 11 0 0 10 1 16 34 ...
Dataset lung cũng có thể được tải từ gói “survival” trong R và chứa các biến sau:
inst: Mã số của bệnh viện
time: Thời gian sống sót hoặc thời gian theo dõi (trong ngày)
status: Tình trạng sống sau thời gian theo dõi (1 = Sống sót, 2 = Qua đời)
age: Tuổi của bệnh nhân
sex: Giới tính của bệnh nhân (1 = Nam, 2 = Nữ)
ph.ecog: Mức độ tàn phế của bệnh nhân (0 = Tốt, 1 = Tàn phế nhẹ, 2 = Tàn phế trung bình, 3 = Tàn phế nặng)
ph.karno: Điểm Karnofsky của bệnh nhân, đo lường sức khỏe tổng thể (100 điểm là sức khỏe tốt nhất, 0 điểm là bệnh nhân đã qua đời)
pat.karno: Điểm Karnofsky của bệnh nhân dựa trên người chăm sóc (nếu có)
meal.cal: Lượng calo tiêu thụ hàng ngày của bệnh nhân
wt.loss: Mức độ giảm cân của bệnh nhân trong 6 tháng trước đó
Các biến inst, time và status được sử dụng trong phân tích sống còn, trong khi các biến còn lại được sử dụng để đánh giá các yếu tố ảnh hưởng đến sống sót của bệnh nhân. Các biến ph.ecog và ph.karno cung cấp thông tin về sức khỏe tổng thể của bệnh nhân, trong khi các biến meal.cal và wt.loss cung cấp thông tin về chế độ ăn uống và mức độ giảm cân của bệnh nhân.
Lưu ý rằng trạng thái được mã hóa theo cách không chuẩn trong bộ dữ liệu này. Thông thường, bạn sẽ thấy 1=sự kiện xảy ra, 0=bị kiểm duyệt. Tôi sẽ thay đổi lại biến “status” để tránh nhầm lẫn
library(dplyr)
lung <- lung %>% mutate(status = recode(status, `1` = 0, `2` = 1))
Bây giờ chúng ta có: 2. time: Thời gian sống sót hoặc thời gian theo dõi (trong ngày) 3. status: Tình trạng sống sau thời gian theo dõi (0 = Sống sót, 1 = Đã chết) 4. age: Tuổi của bệnh nhân
head(lung[, c("time", "status", "sex")])
## time status sex
## 1 306 1 1
## 2 455 1 1
## 3 1010 0 1
## 4 210 1 1
## 5 883 1 1
## 6 1022 0 1
Lưu ý: hàm Surv() trong gói {survival} mặc định chấp nhận TRUE/FALSE, trong đó TRUE là sự kiện và FALSE là kiểm duyệt; 1/0 trong đó 1 là sự kiện và 0 là kiểm duyệt; hoặc 2/1 trong đó 2 là sự kiện và 1 là kiểm duyệt. Hãy cẩn thận để đảm bảo dữ liệu được định dạng đúng.
Trong phân tích sống còn, để thực hiện các phương pháp phân tích dữ liệu sống sót, ta cần tạo đối tượng sống còn (survival object) dựa trên các thông tin về thời gian sống sót và trạng thái sự kiện của các đối tượng trong một nghiên cứu. Sau đó, ta có thể sử dụng các hàm để vẽ đường sống sót và ước lượng các mô hình sống sót.
Trong R, ta có thể tạo đối tượng sống còn bằng hàm Surv(). Cú pháp của hàm này như sau:
Surv(time, event, type = "right")
Trong đó:
time: là vector chứa thông tin về thời gian sống sót hoặc thời gian theo dõi của các đối tượng. Đây có thể là thời gian từ lần khám đầu tiên đến lần khám cuối cùng, thời gian từ điều trị đến sự kiện xảy ra, hoặc thời gian từ lần phát hiện căn bệnh đến sự kiện xảy ra.
event: là vector chứa thông tin về trạng thái sự kiện của các đối tượng. Giá trị của event có thể là 0 (còn sống) hoặc 1 (chết hoặc xảy ra sự kiện quan tâm).
type: là tham số xác định cách tính toán thời gian sống sót. Giá trị mặc định của type là “right”, có nghĩa là thời gian sống sót được tính từ thời điểm khám cuối cùng trừ đi thời điểm xảy ra sự kiện. Nếu giá trị của type là “left”, thời gian sống sót sẽ được tính từ thời điểm khám đầu tiên trừ đi thời điểm xảy ra sự kiện.
Ví dụ, ta có dữ liệu về thời gian sống sót và trạng thái sự kiện của 10 bệnh nhân đầu tiên như sau:
Surv(lung$time, lung$status)[1:10]
## [1] 306 455 1010+ 210 883 1022+ 310 361 218 166
Ta thấy đối tượng 1 có sự kiện tại thời điểm 306 ngày, đối tượng 2 có sự kiện tại thời điểm 455 ngày, đối tượng 3 bị kiểm duyệt tại thời điểm 1010 ngày, v.v.
Hàm survfit() là một hàm dùng để tính toán và vẽ đường sống sót ước tính (survival curves) dựa trên dữ liệu về thời gian sống sót và trạng thái sự kiện của các đối tượng trong một nghiên cứu. Hàm survfit() là một trong những hàm cơ bản của package survival được dùng để phân tích sống còn trên R.
Cú pháp của hàm survfit() như sau:
survfit(formula, data, weights, ...)
Trong đó:
formula: là công thức mô tả mối quan hệ giữa đối tượng sống còn và các biến giải thích. Công thức này có dạng Surv(time, event) ~ x1 + x2 + …
data: là data frame chứa các biến trong công thức.
weights: là vector chứa các trọng số cho mỗi quan sát trong mô hình.
Tôi sẽ dùng lại đối tượng sinh tồn Surv(lung\(time, lung\)status) để tạo các đường cong sinh tồn bằng phương pháp Kaplan-Meier. Tôi tạo đường cong sinh tồn tổng thể cho toàn bộ nhóm, gán nó cho đối tượng surfitlung và xem xét cấu trúc bằng str():
survfitlung <- survfit(Surv(time, status) ~ 1, data = lung)
str(survfitlung)
## List of 16
## $ n : int 228
## $ time : num [1:186] 5 11 12 13 15 26 30 31 53 54 ...
## $ n.risk : num [1:186] 228 227 224 223 221 220 219 218 217 215 ...
## $ n.event : num [1:186] 1 3 1 2 1 1 1 1 2 1 ...
## $ n.censor : num [1:186] 0 0 0 0 0 0 0 0 0 0 ...
## $ surv : num [1:186] 0.996 0.982 0.978 0.969 0.965 ...
## $ std.err : num [1:186] 0.0044 0.00885 0.00992 0.01179 0.01263 ...
## $ cumhaz : num [1:186] 0.00439 0.0176 0.02207 0.03103 0.03556 ...
## $ std.chaz : num [1:186] 0.00439 0.0088 0.00987 0.01173 0.01257 ...
## $ type : chr "right"
## $ logse : logi TRUE
## $ conf.int : num 0.95
## $ conf.type: chr "log"
## $ lower : num [1:186] 0.987 0.966 0.959 0.947 0.941 ...
## $ upper : num [1:186] 1 1 0.997 0.992 0.989 ...
## $ call : language survfit(formula = Surv(time, status) ~ 1, data = lung)
## - attr(*, "class")= chr "survfit"
Một số thành phần chính của đối tượng survfit này sẽ được sử dụng để tạo các đường cong sinh tồn bao gồm:
thời gian: các mốc thời gian mà tại đó đường cong có một bước, tức là ít nhất một sự kiện đã xảy ra
surv: ước tính tỷ lệ sống sót tại thời điểm tương ứng