Применение анализа выживаемости

Анализ выживаемости (survival analysis, SA) — область статистики, в которой основное внимание уделяется оценке времени до наступления каких-либо событий. Обычно анализ выживаемости применяют в клинических исследованиях, чтобы определить эффективность конкретных лекарств (время до смерти пациента), в исследованиях надежности программных систем (время до сбоя) и анализе кредитования (вероятность неуплаты кредитной задолженности).

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

Я опишу базовые концепции SA наряду с C#-реализацией наиболее часто используемого метода оценки Каплана-Майера (Kaplan-Meier estimator). Я воспользуюсь реальным примером для оценки вероятности выживания мобильных приложений.

Вообразите фирму, разрабатывающую ПО, которая выпускает два мобильных приложения: X и Y. Каждое из них разрабатывается разными группами. Фирма стремится выяснить, насколько отказоустойчивы эти мобильные приложения, и определить, является ли одно из приложений значительно менее надежным, чем другое, и не требует ли оно больше усилий для повышения его надежности.

В любой момент на мобильных устройствах может выполняться множество экземпляров X и Y и оставаться в рабочем состоянии. Таким образом, крах мобильного приложения — это как раз то, что наиболее интересно в данном случае. Более длительные периоды времени до события (краха) в этом примере будут указывать на то, что одно из приложений более отказоустойчивое, или имеет большую выживаемость.

Самая базовая концепция SA — это функция выживания.

В демонстрационной программе я буду сначала отображать данные по выживаемости мобильных приложений X и Y (рис. 1). Эти данные показывают, что оба приложения выполнялись десятью пользователями с идентификаторами от 0 до 9. В моем примере приложение либо рухнет (на экранном снимке это описывается как event = app crash), либо будет закрыто пользователем (описывается как event = app off). Также регистрируется день, в которой происходит это событие.

Рис. 1. Демонстрация анализа выживаемости, показывающая жизненные циклы мобильных приложений

Базовые концепции SA

Самая базовая концепция SA — это функция выживания (survival function). Ее принято обозначать как S(t). Если X — случайная переменная (переменная, значения которой присваиваются случайным образом) представляет время события, то S(t) = Pr (X > t). Иначе говоря, S(t) — это вероятность выживания спустя время t. S(0) определяется равной 1. Функция выживания связана с функцией распределения времени жизни, которая обычно обозначается как F(t) и определяется так: F(t) = Pr (X<=t); т. е. она описывает вероятность наступления события до времени t. Следовательно, F(t) = 1 – S(t). Тогда функция плотности распределения f(t) определяется как dF(t)/dt — первая производная F(t), если F(t) является дифференцируемой. Таким образом, f(t) можно рассматривать как частоту события (в единицу времени).

Функция риска (hazard function), или функция интенсивности отказов, — еще одна базовая концепция; она равна f(t)/S(t) и является частотой события ко времени t для индивидуальных сущностей, которые остаются живыми к моменту t.

Вы можете указывать функции выживания параметрически, используя явную функцию или некое семейство функций. Вы также можете вывести их непараметрически из существующих данных, не представляя их в параметризованном формульном виде (parametrized closed form). Полупараметрическая спецификация — нечто среднее между параметрической и непараметрической спецификациями — тоже возможна. Экспоненциальное распределение является простым и распространенным семейством параметризованных функций для описания функций выживания благодаря своим удобным с математической точки зрения свойствам.

Например, S(t) = exp(–0.05t) — это функция выживания из параметризованного экспоненциального распределения (рис. 2). Функции выживания вида S(t) = exp(–at) (где a — параметр, контролирующий интенсивность отказов [уровень риска]) могут описать это распределение. Функция распределения времени жизни имеет вид F(t) = 1 – S(t) = 1 – exp(–at). Рис. 2 помогает нам наглядно представить, как ведут себя функции выживания в течение времени.

Рис. 2. Как ведут себя функции выживания в течение времени

Работая с этой параметрической моделью (parametric model) (в ней для обозначения уровней фактора используются параметры), вы можете вычислять параметры модели на основе имеющихся данных. В случае экспоненциального распределения это параметр a. Один из способов оценки — применение методов оценки по максимуму правдоподобия (Maximum Likelihood Estimation, MLE), но это уже совсем другая история.

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

Цензурирование происходит, когда некоторые наблюдения в наборе данных неполные. К некоему моменту вы потеряли из виду наблюдаемый элемент. В моем примере это означало бы, что мобильное приложение закончило свое выполнение без краха (генерирующего фатальное исключение). Это приложение было корректно закрыто пользователем. Хотя причины завершения приложения без краха могут быть и другими, я исхожу из того, что приложение либо рушится, либо закрывается пользователем.

Метод оценки Каплана-Майера является непараметрическим алгоритмом, который вычисляет функцию выживания.

Существует две основные разновидности цензурирования: справа (right censoring) и слева (left censoring). Цензурирование справа выполняется, когда известно стартовое время, а время события отсутствует. Цензурирование слева происходит, когда время события имеется, а стартового времени нет. В моем примере осуществляется цензурирование справа.

Использование метода оценки Каплана-Майера для вычисления функции выживания

Метод оценки Каплана-Майера (Kaplan-Meier [KM] estimator) является непараметрическим алгоритмом, который вычисления функцию выживания. Подгонка метода оценки KM включает применение высшей математики, в том числе теории мартингалов (martingale theory) и процессов учета, и выходит за рамки этой статьи. Однако реализация самого метода оценки KM прямолинейна и основана на счетчиках.

Подумаем о вычислении метода оценки KM для выживания мобильного приложения X. Этому методу нужно отслеживать три счетчика.

1. Сколько экземпляров мобильного приложения X все еще работают. Это значение представляется переменной atRisk в моей реализации.
2. Количество рухнувших экземпляров. Отслеживается в переменной Crashed.
3. Число экземпляров, корректно завершивших свою работу. Подсчитывается с помощью переменной Censored.

В следующих строках кода (для мобильного приложения X) используется класс CrashMetaData для кодирования данных по выживанию, представленных в табл. 1:

var appX = new CrashMetaData[] {new CrashMetaData{UserID = 0,
  CrashTime = 1, Crashed = false},
           new CrashMetaData{UserID = 1,
              CrashTime = 5, Crashed = true},
           new CrashMetaData{UserID = 2,
              CrashTime = 5, Crashed = false},
           new CrashMetaData{UserID = 3,
              CrashTime = 8, Crashed = false},
           new CrashMetaData{UserID = 4,
              CrashTime = 10, Crashed = false},
           new CrashMetaData{UserID = 5,
              CrashTime = 12, Crashed = true},
           new CrashMetaData{UserID = 6,
              CrashTime = 15, Crashed = false},
           new CrashMetaData{UserID = 7,
              CrashTime = 18, Crashed = true},
           new CrashMetaData{UserID = 8,
              CrashTime = 21, Crashed = false},
           new CrashMetaData{UserID = 9,
              CrashTime = 22, Crashed = true}};

Табл. 1. Данные по выживанию мобильного приложения X

Данные по выживанию содержат время события в сутках (кодированного классом CrashTime) и информацию о том, относится событие к краху приложения или его цензурированию. Если Crashed равно true, приложение рухнуло. Иначе приложение было корректно закрыто (другими словами, подверглось цензурированию). Кроме того, в поле UserID отслеживается экземпляр приложения.

Метод оценки KM реализуется в методе EstimateKaplanMeier. Он делит данные на различные, неперекрывающиеся временные интервалы на основе периодов до событий (в моем примере это крах приложения). В каждом интервале ведутся счетчики.

Важно обратить внимание на счетчик, показывающий, сколько приложений закрывается перед самым событием (это связано с математическим описанием процессов учета). Поэтому за первые пять интервалов в моем примере, охватывающим сутки с 0 до 5, девять из десяти экземпляров были работоспособны и выполнялись до суток номер 5 (один экземпляр был завершен к моменту 1). В интервале далее вплоть до суток номер 5 включительно я наблюдал один крах (который определяет интервал) и завершение двух экземпляров (в сутки под номерами 1 и 5). См. рис. 3.

Рис. 3. Суточные интервалы, созданные методом оценки KM

В таком случае оценка KM для функции выживания является произведением всех интервалов выживания, полученных от счетчиков в разделах (partitions):

1 – (рухнуло в интервале) / (численность группы риска
перед окончанием интервала)

Метод EstimateKapalanMeier возвращает объект класса SurvivalCurve. Он представляет вычисленную функцию выживания. В выводе — ступенчатая функция (step function). Каждая ступень — это значение функции выживания в соответствующем интервале (по оценке KM). Рис. 4 включает часть вывода демонстрационной программы Survival Analysis, соответствующего объекту SurvivalCurve (для приложений X и Y).

Рис. 4. Вывод демонстрационной программы Survival Analysis по оценкам KM для приложений X и Y

Другой вопрос, на который можно ответить, используя оценки KM, — есть ли разница в выживаемости у двух (или более) разных приложений.

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

Рис. 5. Функция выживания, вычисленная по оценке KM для мобильного приложения X

Затем оценку можно использовать для логического определения среднего времени выживания (median survival time), или времени, к которому половина экземпляров будет оставаться работоспособной. Это должно произойти в некий момент между сутками номер 12 (где оценка вероятности выживания составляет 0.711 > 0.5) и номер 18 (где та же оценка дает 0.474 < 0.5). В литераторе по SA существует несколько подходов к описанию того, как точно вычислить эту величину, поскольку она обычно лежит между двумя ступенями.

Я буду определять среднее время выживания как минимальное время выживания, для которого функция выживания равна менее 0.5, что в случае мобильного приложения X дает среднее время выживания в 18 суток. Как можно интерпретировать эту величину? К моменту восемнадцатых суток половина экземпляров мобильного приложения X рушится, а половина остается работоспособной. Эта реализация вычисляет среднее время выживания в методе GetMedianSurvivalTime.

Другой вопрос, на который можно ответить, используя оценки KM, — есть ли разница в выживаемости у двух (или более) разных приложений. Один из подходов к ответу на этот вопрос — визуализация в виде графика оценок KM, соответствующих каждому приложению. Этот тип графика приведен на рис. 6, где сравниваются вычисленные функции выживания приложений X и Y.

Рис. 6. Оценки KM для мобильных приложений X и Y

{Для верстки: здесь придется сделать зеленый график светло-серым, а синий – черным}

Светло-серая кривая представляет функцию выживания приложения X, а черная кривая — функцию выживания приложения Y.

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

Хотя визуализация функций выживания помогает определять различия в выживаемости, в некоторых случаях разница выражена не столь отчетливо. К счастью, есть статистический подход для проверки таких различий строгим формальным способом, который называют Log Rank Test (логарифмический ранговый критерий). Это алгоритм, который проверяет, существует ли значимая разница между двумя (или более) распределениями выживания без параметризации. В литературе по SA можно найти подробные обсуждения этой темы, а большинство статистических библиотек SA включают реализации логарифмического рангового критерия.

Стоит отметить, что имеется другой популярный непараметрический алгоритм для вычисления функции выживания, который называют методом оценки Нельсона-Алена (Nelson-Aalen, NA). NA вычисляет функцию совокупного риска (cumulative hazard function) на основе данных по выживанию. Потом вы можете вывести функцию выживания из этой оценки, используя математическую формулу, которая связывает ее с функцией совокупного риска. Подробнее об этом методе оценки см. в литературе по SA.

Заключение

Я ознакомил вас с базовыми концепциями и терминологией из статистической области анализа выживания. Я показал, как реализовать непараметрический метод оценки Каплана-Майера и применить его к примеру, где сравнивается отказоустойчивость мобильных приложений. Этот метод оценки может помочь определить, есть ли разница в выживаемости между двумя приложениями. Кроме того, я упомянул о строгом статистическом критерии для проверки наличия разницы, который называют логарифмическим ранговым критерием. Другая величина, получаемая мной с помощью метода оценки Каплана-Майера, — среднее время выживания, которое также указывает на различия в выживаемости между приложениями X и Y.

Наконец, я упомянул о методе оценки Нельсона-Алена как об альтернативном непараметрическом методе для вычисления функции выживания. Он вычисляет не саму функцию выживания, как метод Каплана-Майера, а функцию совокупного риска. После этого вы можете вывести функцию выживания из функции совокупного риска.

Хотя визуализация функций выживания помогает определять различия в выживаемости, в некоторых случаях разница выражена не столь отчетливо.

Все это было лишь поверхностным обзором богатейшего поля SA. Области его применения охватывают едва ли не все — от медицины до машиностроения; его методы и алгоритмы реализованы во многих статистических пакетах. С распространением мобильных приложений и корпоративный систем на основе Software as a Service я предвижу, что методы SA могут сыграть роль в мониторинге и повышении качества таких систем.