2017-09-11 4 views
0

Je suis nouveau à l'aide de R, et j'essaie d'utiliser l'analyse de survie afin de trouver une corrélation dans les données censurées. La donnée x est la masse de l'enveloppe des protostars. Les données y sont l'intensité d'une ligne moléculaire observée, et certaines valeurs sont des limites supérieures. Les données sont:R survie survreg ne produisant pas un bon ajustement

x <- c(17.299, 4.309, 7.368, 29.382, 1.407, 3.404, 0.450, 0.815, 1.027, 0.549, 0.018) 
y <- c(2.37, 0.91, 1.70, 1.97, 0.60, 1.45, 0.25, 0.16, 0.36, 0.88, 0.42) 
censor <- c(0, 1, 1, 0, 0, 1, 0, 1, 0, 0, 1) 

J'utilise la fonction survreg de la bibliothèque R survie

modeldata<-survreg(formula=Surv(y,censor)~x, dist="exponential", control = list(maxiter=90)) 

Ce qui donne le résultat suivant:

summary(modeldata) 

Call: 
survreg(formula = Surv(y, censor) ~ x, dist = "exponential", 
control = list(maxiter = 90)) 
Value Std. Error  z  p 
(Intercept) -0.114  0.568 -0.20 0.841 
x   0.153  0.110 1.39 0.163 

Scale fixed at 1 

Exponential distribution 
Loglik(model)= -6.9 Loglik(intercept only)= -9 
Chisq= 4.21 on 1 degrees of freedom, p= 0.04 
Number of Newton-Raphson Iterations: 5 
n= 11 

Cependant, quand je tracer les données et le modèle utilisant la méthode suivante:

plot(x,y,pch=(censor+1)) 
xnew<-seq(0,30) 
model<-predict(modeldata,list(x=xnew)) 
lines(xnew,model,col="red") 

Je reçois ce plot of x and y data; triangles are censored data

Je ne sais pas où je vais mal. J'ai essayé différentes distributions, mais toutes produisent des résultats similaires. La même chose est vrai quand j'utiliser d'autres données, par exemple:

x <- c(1.14, 1.14, 1.19, 0.78, 0.43, 0.24, 0.19, 0.16, 0.17, 0.66, 0.40) 

Je suis pas sûr si j'interprète correctement les résultats.

J'ai essayé d'autres exemples en utilisant la même méthode (par exemple https://stats.idre.ucla.edu/r/examples/asa/r-applied-survival-analysis-ch-1/), et cela fonctionne bien, pour autant que je sache.

Mes questions sont les suivantes:

  1. Suis-je utiliser la fonction correcte pour le montage des données? Si non, lequel serait le plus approprié?

  2. Si la fonction est correcte, pourquoi le modèle ne correspond-il pas très étroitement aux données? Est-ce que cela a à voir avec le complot?

Nous vous remercions de votre aide.

+0

Je pense qu'il est plus juste de dire que vous avez des données tronquées. Et la distribution semble plus logarithmique qu'exponentielle. –

Répondre

0

La "forme" de la relation semble concave vers le bas, donc je l'aurais deviné un ajustement ~ log(x) pourrait être plus approprié:

dfrm <- data.frame(x = c(17.299, 4.309, 7.368, 29.382, 1.407, 3.404, 0.450, 0.815, 1.027, 0.549, 0.018), 
y = c(2.37, 0.91, 1.70, 1.97, 0.60, 1.45, 0.25, 0.16, 0.36, 0.88, 0.42), 
censor= c(0, 1, 1, 0, 0, 1, 0, 1, 0, 0, 1)) 

modeldata<-survreg(formula=Surv(y,censor)~log(x), data=dfrm, dist="loggaussian", control = list(maxiter=90)) 

code Vous semblait approprié:

png(); plot(y~x,pch=(censor+1),data=dfrm) 
xnew<-seq(0,30) 
model<-predict(modeldata,list(x=xnew)) 
lines(xnew,model,col="red"); dev.off() 

enter image description here

modeldata 
Call: 
survreg(formula = Surv(y, censor) ~ log(x), data = dfrm, dist = "loggaussian", 
    control = list(maxiter = 90)) 

Coefficients: 
(Intercept)  log(x) 
0.02092589 0.32536509 

Scale= 0.7861798 

Loglik(model)= -6.6 Loglik(intercept only)= -8.8 
    Chisq= 4.31 on 1 degrees of freedom, p= 0.038 
n= 11 
+0

Merci pour la réponse, cela a fonctionné. J'ai remarqué que dans mon cas using formule = Surv (y, censor) ~ x a fonctionné mieux que formule = Surv (y, censure) ~ log (x), j'essaye toujours de comprendre pourquoi. Dans les deux cas, votre réponse a résolu mon problème. – dragon