При анализе (сопровождении) испытуемых в эксперименте на выживание (например, перечисленных выше) обычно некоторые испытуемые имеют цензуру время выживания. То есть время выживания таких субъектов фактически не соблюдается. Так обстоит дело, например, когда речь идет о смерти как о событии. Другими словами, цензура возникает, когда мы не можем собрать полный набор данных об субъекте (например, потому, что он выбыл из исследования).

  • Время до цензурирования C — это продолжительность между начальным событием и цензурированием. Цензура обычно происходит, когда рассматриваемое исследование заканчивается.
  • Правильная цензура: цензура, которая происходит до события, т. е. никакая информация не может быть собрана о событии, представляющем интерес для этого субъекта.
  • левая цензура: происходит, когда истинное время выживания человека меньше или равно времени, в течение которого информация наблюдалась, например, времени до воздействия COVID.

Чтобы визуализировать это, мы можем построить искусственный цензурированный набор данных следующим образом:

Итак, мы создали набор данных выживания с цензурой справа. Давайте визуализируем этот синтетический набор данных:

На изображении выше субъекты 2, 5 и 6 действительно пострадали от события, в то время как след субъекта 7 был потерян. То есть мы не знаем, пострадал ли субъект 7 от какого-либо события или просто его/ее записи исчезли по какой-то причине (например, из-за того, что он удалился от места, где проводится исследование). Обратите внимание, что в случае изучения события смерти цензура учитывает незарегистрированные случаи смерти, такие как проблемы с дополнительным уведомлением о COVID-19 по всему миру.

Описание времени до события

Начнем с определения основных величин в анализе выживания:

Определение [время выживания]: пусть T ∈ ℝ⁺ — неотрицательная случайная величина (время выживания) и tконкретное произвольное значение для T. Эту величину часто называют время отказа или время до события, и она является центральной величиной в анализе выживаемости.

При работе с временем выживания (отказа) особое внимание следует уделить четкому и однозначному определению начала времени, от которого фактически измеряется время выживания. Например, мы можем рассматривать возраст как время выживания, а начало времени — это, естественно, время рождения субъектов. Другой пример можно найти в медицинских исследованиях, в таких случаях время происхождения может определяться конкретным событием, таким как участие в исследовании или диагностирование какого-либо заболевания. Излишне говорить, что определение фактического события (или сбоя) должно быть четким.

Вернемся к рисунку 1, где мы построили очень небольшой синтетический набор данных о выживании, есть субъекты без зарегистрированного события, это могло произойти либо потому, что человек не испытал событие, либо потому, что мы потеряли след субъекта. Такие наблюдения называются цензурированными, точнее:

Определение [переменная цензора]: пусть d ∈ {0,1} будет переменной ответа, такой что d=1, если событие произошло во время d=0 в случае цензуры; То есть событие не наблюдалось, потому что оно не произошло, или потому, что мы потеряли след субъекта в эксперименте.

Механизм правой цензуры называется независимым, если частота событий (неудачи) субъектов для любого произвольного времени выживания t›0 такая же, как и в случае, когда цензурирование никогда не происходило. .

Плотность смерти

При рассмотрении времени до наступления события (например, смерти) мы можем построить гистограмму количества событий как функции времени. Проще говоря, мы могли бы подобрать кривую к этой гистограмме (например, с помощью оценки плотности ядра) и получить функцию плотности вероятности смерти f(t). После нормализации в любой при заданном времени t доля наблюдаемых событий в популяции определяется по формуле:

Кумулятивная функция распределения смертей. По сути, мы предполагаем, что T (время до события) является случайной величиной с кумулятивной функцией распределения:

Чтобы проиллюстрировать это, мы будем использовать набор данных, который содержит информацию о странах с 1946 по 2008 год. Таким образом, это набор данных с правой цензурой. Среди информации о странах есть два специальных столбца с именами bornyear и endyear, которые обозначают годы, когда конкретная страна (с таким названием) родилась и закончилась. Мы можем использовать время существования стран в качестве случайной переменной времени выживания.

Функция выживания

Определение. Функция выживания S(t) дает вероятность того, что субъект проживет дольше определенного времени t, другими словами:

Дана функция плотности смерти f(t). Площадь под кривой справа от времени t представляет собой долю субъектов в популяции, переживших время t. Мы называем эту величину S(t), и она определяется как:

Функция выживания — это мгновенная норма риска, которую мы вводим далее.

Мгновенная функция опасности

Определение. Мгновенная скорость, с которой случайно выбранный человек, о котором известно, что он жив в момент времени t, умрет в момент времени t+Δt, называется условная частота отказов или мгновенная опасность, λ(t). Дано:

где Δt — бесконечно малый интервал времени.

Обратите внимание, что, в отличие от функции выживания, мгновенная степень опасности фокусируется на событии, а не на выживании до него. Приведенное выше выражение можно записать как:

Где во второй строке мы использовали следующий факт:

а затем мы использовали тот факт, что для очень малых Δt мы можем написать:

В последней строке: f(t) — это плотность смерти и S(t) — функция выживания.

Кроме того, отметив, что:

мы можем написать для функции опасности:

Пример: это соотношение позволяет нам проанализировать очень простой сценарий: сценарий, в котором риск опасности не меняется со временем t. Это: λ=c, где c — постоянное значение. Что оставляет для функции выживания:

Сказанное выше легко проиллюстрировать несколькими строками кода:

Определение. Также естественно определить функцию кумулятивной опасности следующим образом:

Мы можем переписать функции плотности выживания и смерти в терминах этой кумулятивной функции риска.

Обратите внимание, что функция опасности может представлять больший интерес, чем p.d.f. пациенту, который пережил определенный период времени и хотел узнать что-нибудь об их прогнозе. Есть и другие причины для введения функций опасности λ(t) и Λ(t):

  • Интерпретируемость. Предположим, что T обозначает время от операции по поводу рака молочной железы до рецидива. Затем, когда пациент, перенесший операцию, посещает своего врача, его больше интересуют условные вероятности, такие как «Учитывая, что у меня еще не было рецидива, каковы мои шансы на его появление в следующем году», чем безусловные вероятности. (как описано в pdf).
  • Аналитические упрощения. Когда данные подвергаются правильной цензуре, представления функции опасности часто упрощают анализ. Например, представьте, что вы собираете когорту из N пациентов, которым только что исполнилось 50 лет, а затем наблюдаете за ними в течение 1 года. Затем, если d мужчин умирают в течение года наблюдения, соотношение d/N оценивает (дискретная) функция риска T = возраст на момент смерти. Мы увидим, что Λ(·) обладает хорошими аналитическими свойствами.
  • Упрощение моделирования. Для многих биомедицинских явлений T таково, что λ(t) довольно медленно меняется в t. Таким образом, λ(·) хорошо подходит для моделирования.

Ожидаемое время выживания или ожидаемая продолжительность жизни определяется по формуле:

после интегрирования по частям и с учетом того, что

мы получаем:

Функция правдоподобия для цензурированных данных

Предположим, у нас есть когорта из N наблюдений, управляемая функцией выживания S(t) с плотностью смерти f(t) и опасностью λ(t) . Если i-й субъект наблюдается в момент времени tᵢ и субъект умер, его вклад в функцию правдоподобия является:

Если субъект жив в момент времени t=tᵢ, все, что мы знаем, это то, что время выживания Tсубъекта i превышает tᵢ. Вероятность этого события:

вклад в вероятность цензурированного наблюдения.

Пусть dᵢ∈ {0,1} будет переменной-индикатором смерти для i-го наблюдения, которое принимает значение 1, если субъект умер, и 0 в противном случае. Тогда мы можем записать функцию правдоподобия как:

Это оставляет для логарифмического правдоподобия:

Пример. Вернемся к нашему старому доброму примеру с экспоненциальным распределением. Функция степени опасности является постоянной λ(t)=λ. Тогда для кумулятивной опасности имеем Λ(t)=λt. Затем логарифмическая вероятность выглядит следующим образом:

if:

это общее количество смертей и

это общее время наблюдения. Мы можем записать логарифмическую вероятность как:

Мы получаем функцию оценки, дифференцируя это выражение по λ:

установка этой оценки на ноль дает оценку максимальной вероятности опасности:

общее количество смертей, деленное на общую продолжительность, т. е. уровень смертности. Ожидаемая информация получена из:

используя оценку максимального правдоподобия λₘₗ, мы получаем оценку дисперсии:

который используется в качестве доверительного интервала для λ.

Подходы к анализу выживания

Оценка функции выживания

У нас может возникнуть соблазн оценить эту величину как:

где nₛ(t) — количество субъектов, выживших после T=t и n — общее количество субъектов. Но учтите, что при наличии цензуры эту оценку использовать нельзя, так как числитель не всегда определен. Чтобы увидеть это, рассмотрите данные в приведенном ниже примере: субъект, индексированный как 13, имеет время выживания T=25 дней, но переменная ответа status имеет значение 0, что означает наблюдение было подвергнуто цензуре на 25-й день. Для этого субъекта у нас нет информации о возникновении события. Единственное, что мы знаем, это то, что событие не происходило в течение 25 первых дней наблюдения, когда мы потеряли предмет из виду. Эта правильная цензура не позволяет нам вычислить реальное количество выживших субъектов, числитель

оценщик.

Чтобы проиллюстрировать оценку функции выживания, мы будем использовать набор данных, относящийся к Исследованию рака легких, проводимому ветеранами. Он состоит из рандомизированного исследования двух схем лечения рака легких. Набор данных содержит следующие переменные:

  • Лечение: обозначает тип лечения рака легких; стандартный и тестируемый препарат.
  • Celltype: обозначает тип задействованной ячейки; плоскоклеточный, мелкоклеточный, аденозный, крупноклеточный.
  • Karnofsky_score: — показатель Карновского.
  • Diag: – время с момента постановки диагноза в месяцах.
  • Возраст: возраст в годах.
  • Prior_Therapy: обозначает любую предшествующую терапию; нет или да.
  • Статус: обозначает статус пациента как живого или мертвого; мертвый или живой.
  • Survival_in_days – это время выживания в днях после лечения.
df = pd.read_csv('veterans.csv')
df.iloc[[11,5,32,13,23]]

Рассмотрим нашу наивную оценку

Из приведенной выше таблицы мы можем вычислить:

Но мы не можем вычислить:

поскольку мы не знаем, жив ли субъект номер 4 (с индексом 13) в t=30, мы знаем только, что субъект был еще жив в t = 25.

Оценщик Каплана-Мейера

Оценка, которую можно использовать для цензурированных справа данных (как в нашем случае), — это так называемая оценка Каплана-Мейера.

Определение: [Каплан-Мейер] Пусть dᵢ будет количеством событий во времени T=tᵢ . Пусть также nᵢ будет количеством субъектов, о которых известно, что они выжили, т. е. число субъектов без каких-либо событий или подвергнутых цензуре время tᵢ. Оценка Каплана-Мейера для функции выживания определяется как:

Пациенты в этом исследовании получали два лечения, эта информация хранится в столбце trt набора данных. Мы можем использовать этот флажок, чтобы исследовать влияние двух различных методов лечения на приведенные выше кривые выживания:

df['trt'].value_counts()
1    69
2    68
Name: trt, dtype: int64

Приведенные выше оценки Каплана-Мейера для кривых выживания позволяют изучить влияние ковариат на функцию выживания. Этот подход разделения набора данных на более мелкие группы в соответствии с переменной полезен, но невозможен, когда мы хотим рассмотреть более пары переменных, потому что подгруппы очень быстро станут очень маленькими. Для таких ситуаций нам нужны модели! БАМ!

Пропорциональные модели опасности Кокса

Это семейство моделей, представленное Коксом (1972), фокусируется на функции риска, а не на функции выживания. Простейшим членом этого семейства моделей является модель пропорциональных рисков. Для этой модели предполагается, что опасность во время t для человека с ковариатами xᵢ составляет:

где λ₀(t) — базовая функция риска, описывающая риск для лиц с xᵢ=0, что служить ссылкой. Мультипликативный термин

это относительный риск, пропорциональное увеличение или уменьшение риска, зависящее от набора ковариат xᵢ. Обратите внимание, что это увеличение или уменьшение риска не меняется со временем. На самом деле, допущение пропорционального риска удовлетворяется двумя выборками x₁ и x₂, когда:

постоянна во времени. Логарифмируя, получаем:

что дает линейную модель логарифма опасности. Возвращаясь к выражению для опасности, мы можем интегрировать его, чтобы получить кумулятивные функции опасности:

возводя в степень, получаем функции выживания:

Таким образом, влияние фактора относительного риска на функцию выживания заключается в возведении ее в степень. Вероятность наблюдаемого события для подмножества i в момент времени tᵢ может быть записана как:

tшляпа представляет вероятность того, что субъект iстолкнется с событием во время tᵢсреди те, кто находится в группе риска во время tᵢ. Обратите внимание, что базовая опасность в приведенном выше выражении аннулируется, поэтому мы можем написать для вероятности

Приведенное выше выражение необходимо скорректировать для связанных событий, то есть когда число смертейd во время tᵢбольше чем 1. Это дает:

Пример. Вернемся к набору данных о раке легких и попробуем оценить показатели риска и кривую выживаемости.

Теперь мы можем вызвать метод predict_survival и построить кривые выживания:

Экспоненциальные модели и модели Вейбулла

Различные виды моделей пропорциональных рисков получаются в зависимости от различных вариантов базовой функции выживания. Например, мы могли бы (по какой-то причине) установить базовую функцию риска λ₀(t)=λ₀ как постоянное значение. Тогда опасность субъекта i в момент времени t выглядит следующим образом:

Другая распространенная модель получается, если предположить, что функция плотности смертей соответствует распределению Вейбулла, из которого функция выживания равна:

и функция опасности

где λ›0 и p›0 — параметры. Если p=1, эта модель сводится к экспоненциальной модели и имеет постоянный риск во времени. Если p›1, то риск со временем увеличивается. Если p‹1, то риск со временем снижается. Это можно увидеть, взяв лог приведенного выше выражения:

Изменяющиеся во времени ковариаты

Характер модели позволяет рассматривать зависимые от времени ковариаты. Пусть xᵢ(t) обозначает значение вектора ковариации для субъекта iв момент времени (или продолжительности) t. Модель пропорциональной опасности теперь гласит:

Зависящие от времени «Эффекты»

Модель также можно обобщить, так что коэффициенты β(t) теперь являются функцией времени. Пример: возможно, что определенные социальные характеристики могут иметь большое влияние на опасность для детей вскоре после рождения, но могут иметь относительно небольшое влияние в более позднем возрасте. Итак, мы хотели бы что-то вроде:

Модель общей степени опасности

Естественно, мы можем положить все на время, чтобы получить наиболее общую версию модели степени опасности:

Или… в более общем смысле:

где fₘₗ обозначает вашу любимую модель машинного обучения или любой другой аппроксиматор функций.

Дискретные модели времени

Бывают случаи, когда временная переменная на самом деле является «дискретным» значением, что в основном означает, что она не непрерывна. Одним из примеров является использование возраста в качестве предиктора риска смерти. Наличие дискретного описания полезно. На практике именно дискретное описание используется для аппроксимации функций с помощью наших компьютеров.

Пусть время выживания T — дискретная случайная величина, принимающая значения t₁ ‹ t₂ ‹ … ‹ tₙ с функция плотности:

Тогда функция выжившего в момент времени tⱼ будет следующей:

Затем нам нужно определить опасность в момент времени tⱼ как условную вероятность смерти в этот момент при условии, что человек дожил до этого момента, так что:

Другим интересным результатом дискретного времени является то, что функция выживания во времени может быть записана в терминах опасности во все предыдущие моменты времени t₁, t₂, …, tⱼ₋₁ как :

что имеет смысл, так как нужно выжить t₁, тогда нужно выжить t₂, учитывая, что выжил t₁ и т. д., в конце концов сохранившись до tⱼ₋₁до этого момента.

Измерение эффективности моделей выживания

Обычные метрики, такие как RMSE или корреляция, могут не подходить для цензурированных переменных. Наиболее часто используемой оценочной метрикой моделей выживания является индекс конкордации (c-index, c-statistic). Это ранговая корреляция между прогнозируемыми показателями риска f и наблюдаемыми моментами времени t, которая тесно связана с τ. Он определяется как отношение правильно упорядоченных (согласованных) пар к сопоставимым парам.

То есть два образца iи jназываются сопоставимыми, если образец с более низким временем выживания tпережил событие, т. е. если tⱼ › tᵢ и dᵢ=1. Сопоставимая пара является конкордантной, если предполагаемый риск f по модели выживания выше для субъектов с меньшим временем выживания, т. е.

в противном случае говорят, что пара дискордантна. По сути, речь идет о следующем количестве:

Харрелл предложил оценку этой величины, реализованную в scikit-survival. Известен как индекс соответствия Харрелла, с-индекс или с-статистика. Интерпретация идентична нашей обычной кривой ROC.

Хотя индекс конкордации Харрелла легко интерпретировать и вычислять, у него есть некоторые недостатки:

  1. было показано, что это слишком оптимистично с увеличением количества цензуры
  2. это бесполезная мера эффективности, если основной интерес представляет конкретный диапазон времени (например, прогнозирование смерти в течение 2 лет).

Первый пункт решается с помощью процедуры взвешивания, известной как IPCW (обратная вероятность цензурирующих весов).

Вторая проблема решается путем распространения хорошо известной кривой ROC на возможно подвергнутое цензуре время выживания. Учитывая момент времени t, мы можем оценить, насколько хорошо прогностическая модель может различать субъектов, которые испытают событие, по времени t(чувствительность) от тех, кто не будет (конкретность). Результирующая AUC для этой зависящей от времени кумулятивной/динамической кривой ROC:

Оценщик для этой величины реализован в scikit-survival и определяется как:

где ωᵢ — обратная вероятность цензурирующих весов (IPCW). Для оценки IPCW требуется доступ к времени выживания из обучающих данных для оценки распределения цензурирования. Обратите внимание, что для этого требуется, чтобы время выживания выживания_тест лежало в диапазоне времени выживания выживания_поезда. Это может быть достигнуто путем соответствующего указания времени, например. установив значение times[-1] немного ниже максимального ожидаемого времени последующего действия. IPCW вычисляются с использованием оценки Каплана-Мейера, которая ограничена ситуациями, когда выполняется предположение о случайной цензуре, а цензура не зависит от признаков. Наконец, функция также предоставляет единую сводную меру, которая относится к среднему значению AUC(t) во временном диапазоне (τ₁,τ₂). :

где S(t) — оценка Каплана-Мейера функции выживания.

Давайте загрузим еще один встроенный набор данных из scikit-survival. Набор данных содержит 7874 образца и 9 признаков:

  • возраст: возраст в годах
  • пол: Ж=женщина, М=мужчина
  • sample.yr: календарный год, в котором был получен образец крови.
  • каппа: бессывороточная легкая цепь, каппа-часть
  • лямбда: бессывороточная легкая цепь, лямбда-часть
  • flc.grp: группа легких цепей, не содержащихся в сыворотке, для субъекта, используемая в исходном анализе.
  • креатинин: креатинин сыворотки
  • mgus: была ли у субъекта диагностирована моноклональная гаммапотия (MGUS)
  • глава: для тех, кто умер, группировка их основной причины смерти по заголовкам глав Международного кодекса болезней МКБ-9
  • Конечной точкой является смерть, которая наступила у 2169 человек (27,5%).

Давайте обработаем некоторые пропущенные значения:

Нам нужно быть немного осторожными при выборе тестовых данных и моментов времени, в которые мы хотим оценить ROC, из-за зависимости оценки от обратной вероятности цензурирующего взвешивания. Во-первых, мы собираемся проверить, находится ли наблюдаемое время тестовых данных в наблюдаемом диапазоне времени обучающих данных.

При выборе временных точек для оценки ROC важно не забывать выбирать последнюю временную точку таким образом, чтобы вероятность цензурирования после последней временной точки была ненулевой. Здесь мы используем более консервативный подход, установив верхнюю границу 80% процентиля наблюдаемых моментов времени, потому что уровень цензурирования довольно велик и составляет 72,5%.

На графике показана расчетная площадь под зависящим от времени ROC в каждый момент времени и среднее значение по всем моментам времени в виде пунктирной линии.

Мы видим, что возраст в целом является наиболее отличительной чертой, за которой следуют κ и λ FLC. Тот факт, что возраст является самым сильным предиктором общей выживаемости населения в целом, вряд ли удивителен (в конце концов, в какой-то момент мы должны умереть). Больше различий становится очевидным при рассмотрении времени: дискриминационная способность FLC снижается в более поздние моменты времени, а с возрастом увеличивается. Наблюдение за возрастом снова следует здравому смыслу. Напротив, СЛЦ кажется хорошим предиктором смерти в ближайшем будущем, но не так сильно, если это произойдет спустя десятилетия.

Оценка прогнозов модели

Модель пропорциональных опасностей Кокса

Подбирайте пропорциональную модель опасности Кокса к данным о поезде:

Используя тестовые данные, мы хотим оценить, насколько хорошо модель может отличать выживших от умерших с недельными интервалами в течение 6 месяцев после зачисления.

График показывает, что модель работает в среднем умеренно хорошо с AUC ≈ 0,72 (пунктирная линия). Тем не менее, есть явная разница в производительности между первой и второй половиной временного диапазона. Производительность на тестовых данных увеличивается до 56 дней после регистрации, остается высокой до 98 дней, а затем быстро падает. Таким образом, мы можем сделать вывод, что модель наиболее эффективна для прогнозирования смертельных случаев в среднесрочной перспективе.

Случайные леса выживания

Недостатком модели пропорциональных опасностей Кокса является то, что она может предсказывать только оценку риска, не зависящую от времени (из-за встроенного предположения о пропорциональных опасностях). Таким образом, одна прогнозируемая оценка риска должна хорошо работать для каждой временной точки. Напротив, случайный лес выживания не имеет этого ограничения. Итак, давайте подгоним такую ​​модель к обучающим данным. Подробности внутренней работы Random Survival Forests можно найти здесь.

Для прогнозирования мы не вызываем функцию «предсказать», которая возвращает независимую от времени оценку риска, а вызываем функцию «предсказание_кумулятивная_опасность_функция», которая возвращает функцию риска с течением времени для каждого пациента. Мы получаем зависящие от времени оценки риска, оценивая каждую кумулятивную функцию риска в интересующие нас моменты времени.

Действительно, Random Survival Forest в среднем работает немного лучше, в основном из-за лучшей производительности в интервалах 25–50 дней и 112–147 дней. Более 147 дней на самом деле дела обстоят хуже. Это показывает, что средняя AUC удобна для оценки общей производительности, но может скрывать интересные характеристики, которые становятся видимыми только при просмотре AUC в отдельные моменты времени.

Выводы и дополнительные ссылки

ИМХО Анализ выживания заключается в том, чтобы задавать правильные вопросы:

  • Имеет ли смысл решать ту или иную задачу как обычную задачу классификации?
  • Что мы пытаемся моделировать?
  • Является ли время важной переменной для предсказания?

Существует значительное количество ресурсов (с открытым исходным кодом), написанных на питоне, таких как pySurvival, scikit-survival и lifelines. Кроме того, люди в Loft умело использовали XGBoost и некоторые другие методы для решения проблемы анализа выживания, известной как XGBoost-survival-embeddings.

Исходный код, использованный для написания этой статьи, можно найти по адресу: https://github.com/pibieta/blog-posts/blob/main/survival-analysis/intro-survival-analysis.ipynb.