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?
plot(Predict(rcs.ASDH, Year))
. Możesz zmieniać inne zmienne, tworząc różne krzywe, wykonując takie czynnościplot(Predict(rcs.ASDH, Year, age=c(25, 35)))
.anova(rcs.ASDH)
.Odpowiedzi:
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 .χ2
Dwa zalecane sposoby oceny dopasowania modelu to:
Istnieją pewne zalety splajnów regresji w porównaniu z wielomianami ułamkowymi, w tym:
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
rms
pakietrcs
funkcji. Aby zapoznać się z krzywymi kalibracji bootstrap karanymi za przeregulowanie, patrzrms
calibrate
funkcja.źródło