Логистическая регрессия со сплайнами регрессии в R

12

Я разрабатывал модель логистической регрессии на основе ретроспективных данных из национальной базы данных о травмах головы в Великобритании. Ключевым результатом является 30-дневная смертность (обозначается как «выживаемая» мера). Другие меры с опубликованным доказательством существенного влияния на результат в предыдущих исследованиях включают в себя:

Year - Year of procedure = 1994-2013
Age - Age of patient = 16.0-101.5
ISS - Injury Severity Score = 0-75
Sex - Gender of patient = Male or Female
inctoCran - Time from head injury to craniotomy in minutes = 0-2880 (After 2880 minutes is defined as a separate diagnosis)

Используя эти модели, учитывая дихотомическую зависимую переменную, я построил логистическую регрессию, используя lrm.

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

Следуя совету Фрэнка Харрелла, я приступил к использованию сплайнов регрессии для моделирования МКС (преимущества этого подхода выделены в комментариях ниже). Таким образом, модель была предварительно определена следующим образом:

rcs.ASDH<-lrm(formula = Survive ~ Age + GCS + rcs(ISS) +
    Year + inctoCran + oth, data = ASDH_Paper1.1, x=TRUE, y=TRUE)

Результаты модели были:

> rcs.ASDH

Logistic Regression Model

lrm(formula = Survive ~ Age + GCS + rcs(ISS) + Year + inctoCran + 
    oth, data = ASDH_Paper1.1, x = TRUE, y = TRUE)

                      Model Likelihood     Discrimination    Rank Discrim.    
                         Ratio Test            Indexes          Indexes       
Obs          2135    LR chi2     342.48    R2       0.211    C       0.743    
 0            629    d.f.             8    g        1.195    Dxy     0.486    
 1           1506    Pr(> chi2) <0.0001    gr       3.303    gamma   0.487    
max |deriv| 5e-05                          gp       0.202    tau-a   0.202    
                                           Brier    0.176                     

          Coef     S.E.    Wald Z Pr(>|Z|)
Intercept -62.1040 18.8611 -3.29  0.0010  
Age        -0.0266  0.0030 -8.83  <0.0001 
GCS         0.1423  0.0135 10.56  <0.0001 
ISS        -0.2125  0.0393 -5.40  <0.0001 
ISS'        0.3706  0.1948  1.90  0.0572  
ISS''      -0.9544  0.7409 -1.29  0.1976  
Year        0.0339  0.0094  3.60  0.0003  
inctoCran   0.0003  0.0001  2.78  0.0054  
oth=1       0.3577  0.2009  1.78  0.0750  

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

plot(calibrate(rcs.ASDH, B=1000), main="rcs.ASDH")

Кривые начальной загрузки оштрафованы за переоснащение

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

ASDH <- Predict(rcs.ASDH, Year=seq(1994,2013,by=1),Age=48.7,ISS=25,inctoCran=356,Other=0,GCS=8,Sex="Male",neuroYN=1,neuroFirst=1)
Probabilities <- data.frame(cbind(ASDH$yhat,exp(ASDH$yhat)/(1+exp(ASDH$yhat)),exp(ASDH$lower)/(1+exp(ASDH$lower)),exp(ASDH$upper)/(1+exp(ASDH$upper))))
names(Probabilities) <- c("yhat","p.yhat","p.lower","p.upper")
ASDH<-merge(ASDH,Probabilities,by="yhat")
plot(ASDH$Year,ASDH$p.yhat,xlab="Year",ylab="Probability of Survival",main="30 Day Outcome Following Craniotomy for Acute SDH by Year", ylim=range(c(ASDH$p.lower,ASDH$p.upper)),pch=19)
arrows(ASDH$Year,ASDH$p.lower,ASDH$Year,ASDH$p.upper,length=0.05,angle=90,code=3)

Код выше привел к следующему выводу:

Годовой тренд с нижним и верхним

Мои оставшиеся вопросы следующие:

1. Интерпретация сплайнов. Как рассчитать значение p для сплайнов, объединенных для общей переменной?

Дан Фонтан
источник
4
Хорошо сделано. Для отображения эффекта года я предлагаю дать другие переменной быть установлено на значения по умолчанию (медиана для непрерывного режима для категориального) и изменяются год на оси х, например, plot(Predict(rcs.ASDH, Year)). Вы можете позволить другим переменным изменяться, составляя различные кривые, делая такие вещи, как plot(Predict(rcs.ASDH, Year, age=c(25, 35))).
Фрэнк Харрелл
1
Я не знаю почему - но я не видел много примеров калибровочных кривых с поправкой на смещение в литературе. Похоже, хорошая идея
Чарльз
1
Для проверки общей связи с несколькими тестами DF используйте anova(rcs.ASDH).
Фрэнк Харрелл

Ответы:

8

Очень трудно интерпретировать ваши результаты, когда вы не указали модель заранее, но провели значительное тестирование и в результате доработали модель. И я не рекомендую использовать глобальный критерий соответствия, который вы использовали, потому что он имеет вырожденное распределение, а не .χ2

Два рекомендуемых способа оценки соответствия модели:

  1. Начальная корректировка гладкой непараметрической (например, * лёсс) калибровочной кривой с корректировкой при загрузке для проверки абсолютной точности прогнозов
  2. Изучите различные обобщения модели, проверяя, работает ли более гибкая спецификация модели. Например, сделайте отношение правдоподобия или Wald тесты дополнительных нелинейных или взаимодействующих членов.χ2

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

  1. Предиктор может иметь значения ; вы не можете вести журналы таких значений, как этого требуют FP0
  2. Вам не нужно беспокоиться о происхождении предиктора. ФП предполагают, что ноль является «магическим» источником для предикторов, тогда как сплайны регрессии инвариантны к смещению предиктора на константу.

Для получения дополнительной информации о сплайнах регрессии и оценке линейности и аддитивности см. Мои раздаточные материалы по адресу http://biostat.mc.vanderbilt.edu/CourseBios330, а также функцию rmsпакета R. rcsДля калибровочных кривых начальной загрузки, штрафуемых за переоснащение, смотрите rms calibrateфункцию.

Фрэнк Харрелл
источник
Я попытался заранее указать модель при создании полной модели со всеми клинически доступными переменными, включая первоначально и с известной релевантностью к диагнозу и переменной результата. Все они были линейными непрерывными или дихотомическими переменными, за исключением МКС, которую, как показали предыдущие исследования, можно моделировать дробными полиномами. Метод, который я использовал для разработки модели, я считаю, согласуется с «Многофакторным построением модели» Вилли Сауэрбреи. Возможно, я мог бы использовать ваш RMS-пакет в R, чтобы оценить глобальную пригодность? Если да, то какую формулу вы бы порекомендовали?
Дан Фонтан
Я расширил свой ответ, чтобы решить некоторые из них.
Фрэнк Харрелл
Не могли бы вы порекомендовать пакеты для выполнения 1 и 2 для оценки соответствия модели? Могу ли я использовать функцию проверки в RMS? Что касается сплайнов регрессии, я читаю ваши лекционные заметки и в настоящее время пытаюсь найти копию вашей книги (я бы купил ее, если бы мог!). Можете ли вы порекомендовать еще один пакет / функцию R для построения сплайнов регрессии? для нелинейных непрерывных переменных в этой модели? Большое спасибо за вашу помощь.
Дан Фонтан
Возможно, земной пакет основан на многомерных адаптивных сплайнах регрессии Фридмана?
Дан Фонтан
И последний вопрос. Формула для оценки тяжести травмы (ISS): A ^ 2 + B ^ 2 + C ^ 2, где A, B и C - разные части тела с независимой оценкой тяжести травмы 1-5. Таким образом, максимум составляет 75 и минимум 1 в этом наборе данных. Учитывая эту формулу, будут ли дробные полиномы ближе к тому, как на самом деле рассчитывается оценка по сравнению со сплайнами регрессии?
Дан Фонтан