Logistische Regression mit Regressionssplines in R

12

Ich habe ein logistisches Regressionsmodell entwickelt, das auf retrospektiven Daten aus einer nationalen Traumadatenbank für Kopfverletzungen in Großbritannien basiert. Das Hauptergebnis ist die 30-Tage-Mortalität (als "Überlebensmaß" bezeichnet). Andere Maßnahmen mit veröffentlichten Hinweisen auf signifikante Auswirkungen auf das Ergebnis in früheren Studien umfassen:

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)

Unter Verwendung dieser Modelle habe ich angesichts der dichotomen abhängigen Variablen eine logistische Regression mit lrm erstellt.

Die Methode zur Auswahl der Modellvariablen basierte auf vorhandener klinischer Literatur, die dieselbe Diagnose modellierte. Alle wurden mit einer linearen Anpassung modelliert, mit Ausnahme der ISS, die traditionell durch fraktionierte Polynome modelliert wurde. Keine Veröffentlichung hat bekannte signifikante Wechselwirkungen zwischen den obigen Variablen identifiziert.

Auf Anraten von Frank Harrell habe ich Regressionssplines zum Modellieren der ISS verwendet (in den Kommentaren unten werden die Vorteile dieses Ansatzes hervorgehoben). Das Modell wurde daher wie folgt vorab spezifiziert:

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

Ergebnisse des Modells waren:

> 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  

Ich habe dann die Kalibrierungsfunktion im RMS-Paket verwendet, um die Genauigkeit der Vorhersagen aus dem Modell zu bewerten. Die folgenden Ergebnisse wurden erhalten:

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

Bootstrap-Kalibrierungskurven werden wegen Überanpassung bestraft

Nach Abschluss des Modellentwurfs habe ich das folgende Diagramm erstellt, um die Auswirkung des Jahres des Vorfalls auf das Überleben zu veranschaulichen, wobei die Werte des Medians in kontinuierlichen Variablen und der Modus in kategorialen Variablen zugrunde gelegt wurden:

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)

Der obige Code ergab die folgende Ausgabe:

Jahrestrend mit unterer und oberer

Meine verbleibenden Fragen sind die folgenden:

1. Spline-Interpretation - Wie kann ich den p-Wert für die für die Gesamtvariable kombinierten Splines berechnen?

Dan Brunnen
quelle
4
Gute Arbeit. Um den Effekt von Year anzuzeigen, empfehle ich, die anderen Variablen auf Standardwerte (Median für Continuous, Modus für Categorical) zu setzen und Year auf der x-Achse zu variieren, z plot(Predict(rcs.ASDH, Year)). Sie können andere Variablen variieren lassen und verschiedene Kurven bilden, indem Sie Dinge wie tun plot(Predict(rcs.ASDH, Year, age=c(25, 35))).
Frank Harrell
1
Ich weiß nicht warum - aber ich habe in der Literatur nicht viele Beispiele für vorspannungskorrigierte Kalibrierungskurven gesehen. Scheint eine gute Idee zu sein
Charles
1
Um die Gesamtassoziation mit mehreren df-Tests zu testen, verwenden Sie anova(rcs.ASDH).
Frank Harrell

Antworten:

8

χ2

Zwei empfohlene Methoden zur Bewertung der Modellanpassung sind:

  1. Bootstrap-Überanpassung-korrigierte glatte nichtparametrische (z. B. * Löss) Kalibrierungskurve zur Überprüfung der absoluten Genauigkeit von Vorhersagen
  2. χ2

Es gibt einige Vorteile von Regressionssplines gegenüber fraktionellen Polynomen, einschließlich:

  1. 0
  2. Sie brauchen sich keine Gedanken über die Herkunft eines Prädiktors zu machen. FPs gehen davon aus, dass Null ein "magischer" Ursprung für Prädiktoren ist, wohingegen Regressionssplines nicht dazu geeignet sind, einen Prädiktor um eine Konstante zu verschieben.

Weitere Informationen zu Regressionssplines und zur Bewertung von Linearität und Additivität finden Sie in meinen Handouts unter http://biostat.mc.vanderbilt.edu/CourseBios330 sowie in der R- rmsPaketfunktion rcs. Für Bootstrap-Kalibrierungskurven, die wegen Überanpassung bestraft werden, siehe die rms calibrateFunktion.

Frank Harrell
quelle
Ich habe versucht, das Modell vorab zu spezifizieren, indem ich ein vollständiges Modell mit allen klinisch verfügbaren Variablen einschließlich der anfänglichen und bekannten Relevanz für die Diagnose- und Ergebnisvariable erstellt habe. Sie waren alle linear kontinuierliche oder dichotome Variablen mit Ausnahme von ISS, die in früheren Studien identifiziert wurden und mit fraktionellen Polynomen modelliert werden können. Die Methode, mit der ich das Modell entwickelt habe, steht im Einklang mit Willi Sauerbreys "Multivariate Model-Building". Vielleicht könnte ich Ihr RMS-Paket in R verwenden, um die globale Anpassungsgüte zu bewerten? Wenn ja, welche Formel würden Sie empfehlen?
Dan Fountain
Ich habe meine Antwort erweitert, um einige davon anzusprechen.
Frank Harrell
Könnten Sie Pakete für die Ausführung von 1 und 2 empfehlen, um die Modellanpassung zu bewerten? Könnte ich die Validierungsfunktion in Effektivwert verwenden? In Bezug auf Regressionssplines lese ich Ihre Vorlesungsunterlagen und versuche derzeit, eine Kopie Ihres Buches zu finden (ich würde es kaufen, wenn ich könnte!). Können Sie ein weiteres R-Paket / eine weitere R-Funktion für die Konstruktion von Regressionssplines empfehlen? für die nichtlinearen stetigen Variablen in diesem Modell? Vielen Dank für all Ihre Hilfe.
Dan Fountain
Vielleicht das Erdpaket, das auf Friedmans multivariaten adaptiven Regressionssplines basiert?
Dan Fountain
Eine letzte Frage. Die Formel für den Schweregrad von Verletzungen (ISS) lautet A ^ 2 + B ^ 2 + C ^ 2, wobei A, B und C unterschiedliche Körperteile mit einem unabhängigen Schweregrad von 1 bis 5 sind. Somit ist das Maximum 75 und das Minimum 1 innerhalb dieses Datensatzes. Wären nach dieser Formel Bruchpolynome die nähere Darstellung zur tatsächlichen Berechnung des Scores im Vergleich zu Regressionssplines?
Dan Fountain