Regresja logistyczna z splajnami regresji w R.

Opracowuję model regresji logistycznej oparty na danych retrospektywnych z krajowej bazy danych dotyczących urazów głowy w Wielkiej Brytanii. Kluczowym rezultatem jest 30-dniowa śmiertelność (oznaczona jako miara „przetrwania”). Inne miary z opublikowanymi dowodami znaczącego wpływu na wyniki poprzednich badań obejmują:

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)

Używając tych modeli, biorąc pod uwagę dychotomiczną zmienną zależną, zbudowałem regresję logistyczną za pomocą lrm.

Metodę selekcji zmiennych modelowych oparto na istniejącej literaturze klinicznej modelującej tę samą diagnozę. Wszystkie zostały modelowane z dopasowaniem liniowym, z wyjątkiem ISS, który został tradycyjnie modelowany za pomocą wielomianów ułamkowych. Żadna publikacja nie zidentyfikowała znanych znaczących interakcji między powyższymi zmiennymi.

Zgodnie z radą Franka Harrella zastosowałem splajny regresji do modelu ISS (zalety tego podejścia podkreślono w komentarzach poniżej). Model został więc wstępnie określony w następujący sposób:

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

Wyniki modelu to:

> 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  

Następnie użyłem funkcji kalibracji w pakiecie rms w celu oceny dokładności prognoz z modelu. Następujące wyniki zostały osiągnięte:

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

Po zakończeniu projektowania modelu stworzyłem następujący wykres, aby zademonstrować wpływ roku zdarzenia na przeżycie, opierając wartości mediany w zmiennych ciągłych i tryb w zmiennych kategorialnych:

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)

Powyższy kod zaowocował następującymi danymi wyjściowymi:

Moje pozostałe pytania to:

1. Interpretacja splajnu - Jak obliczyć wartość p dla splajnów połączonych dla zmiennej ogólnej?

Bardzo trudno jest zinterpretować wyniki, jeśli nie określono wcześniej modelu, ale przeprowadzono znaczące testy i wynikające z niego modyfikacje modelu. I nie polecam zastosowanego globalnego testu dopasowania, ponieważ ma on zdegenerowany rozkład, a nie .$\chi^2$

Dwa zalecane sposoby oceny dopasowania modelu to:

1. Łagodna, nieparametryczna (np. * Lessowa) krzywa kalibracji bootstrapu korygująca przepełnienie w celu sprawdzenia absolutnej dokładności prognoz
2. Zbadaj różne uogólnienia modelu, sprawdzając, czy bardziej elastyczna specyfikacja modelu działa lepiej. Na przykład wykonaj testy prawdopodobieństwa lub testy Wald dodatkowych warunków nieliniowych lub interakcji.$\chi^2$

Istnieją pewne zalety splajnów regresji w porównaniu z wielomianami ułamkowymi, w tym:

1. Predyktor może mieć wartości ; nie możesz brać dzienników takich wartości, jakie są wymagane przez FP$\leq 0$
2. Nie musisz się martwić o pochodzenie predyktora. FP zakładają, że zero jest „magicznym” źródłem dla predyktorów, podczas gdy splajny regresji są niezmienne dla przesunięcia predyktora o stałą.

Więcej informacji na temat wypustami regresji i liniowość i oceny addytywności zobaczyć moje materiały informacyjne na http://biostat.mc.vanderbilt.edu/CourseBios330 jak również R rmspakiet rcsfunkcji. Aby zapoznać się z krzywymi kalibracji bootstrap karanymi za przeregulowanie, patrz rms calibratefunkcja.