34

ワイブル モデルを生存データに適合させてプロットしようとしています。データには、2006 年から 2010 年まで実行されるコホートという 1 つの共変量しかありません。では、2010 年のコホートの生存曲線をプロットするために続く 2 行のコードに何を追加すればよいでしょうか?

library(survival)
s <- Surv(subSetCdm$dur,subSetCdm$event)
sWei <- survreg(s ~ cohort,dist='weibull',data=subSetCdm)

Cox PH モデルで同じことを達成するのは、次の行でかなり簡単です。問題は、survfit() が survreg 型のオブジェクトを受け入れないことです。

sCox <- coxph(s ~ cohort,data=subSetCdm)
cohort <- factor(c(2010),levels=2006:2010)
sfCox <- survfit(sCox,newdata=data.frame(cohort))
plot(sfCox,col='green')

（生存パッケージからの）データ肺を使用して、これが私が達成しようとしていることです。

#create a Surv object
s <- with(lung,Surv(time,status))

#plot kaplan-meier estimate, per sex
fKM <- survfit(s ~ sex,data=lung)
plot(fKM)

#plot Cox PH survival curves, per sex
sCox <- coxph(s ~ as.factor(sex),data=lung)
lines(survfit(sCox,newdata=data.frame(sex=1)),col='green')
lines(survfit(sCox,newdata=data.frame(sex=2)),col='green')

#plot weibull survival curves, per sex, DOES NOT RUN
sWei <- survreg(s ~ as.factor(sex),dist='weibull',data=lung)
lines(survfit(sWei,newdata=data.frame(sex=1)),col='red')
lines(survfit(sWei,newdata=data.frame(sex=2)),col='red')
4

5 に答える 5

25

これが役に立ち、誤解を招くような間違いを犯していないことを願っています:

上からコピー：

    #create a Surv object
    s <- with(lung,Surv(time,status))

    #plot kaplan-meier estimate, per sex
    fKM <- survfit(s ~ sex,data=lung)
    plot(fKM)

    #plot Cox PH survival curves, per sex
    sCox <- coxph(s ~ as.factor(sex),data=lung)
    lines(survfit(sCox,newdata=data.frame(sex=1)),col='green')
    lines(survfit(sCox,newdata=data.frame(sex=2)),col='green')

Weibull については、予測を使用して、Vincent からのコメントを参照してください。

    #plot weibull survival curves, per sex,
    sWei <- survreg(s ~ as.factor(sex),dist='weibull',data=lung)

    lines(predict(sWei, newdata=list(sex=1),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red")
    lines(predict(sWei, newdata=list(sex=2),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red")

プロット出力

ここでのトリックは、プロットと予測の分位数の順序を逆にすることでした。これを行うためのより良い方法がある可能性がありますが、ここでは機能します。幸運を！

于 2012-02-06T20:48:35.977 に答える
18

別のオプションは、パッケージを利用することですflexsurv。これにより、パッケージにいくつかの追加機能が提供されsurvivalます-パラメトリック回帰関数flexsurvreg()には、あなたが求めることを行う素敵なプロットメソッドがあることを含みます.

上記のように肺を使用します。

#create a Surv object
s <- with(lung,Surv(time,status))

require(flexsurv)
sWei  <- flexsurvreg(s ~ as.factor(sex),dist='weibull',data=lung)
sLno  <- flexsurvreg(s ~ as.factor(sex),dist='lnorm',data=lung)   

plot(sWei)
lines(sLno, col="blue")

plot.flexsurvreg からの出力

引数を使用して累積ハザードまたはハザード スケールにプロットし、type引数を使用して信頼区間を追加できciます。

于 2014-03-14T13:55:21.257 に答える
9

これは、次のコードを使用するTim Riffe の answerを明確にするメモです。

lines(predict(sWei, newdata=list(sex=1),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red")
lines(predict(sWei, newdata=list(sex=2),type="quantile",p=seq(.01,.99,by=.01)),seq(.99,.01,by=-.01),col="red")

2 つの鏡像シーケンスseq(.01,.99,by=.01)およびのseq(.99,.01,by=-.01)理由は、predict() メソッドがイベント分布 f(t) の分位数、つまり f(t) の逆 CDF の値を与えるためです。一方、生存曲線は1-(f の CDF) 対 t をプロットします。つまり、p と predict(p) をプロットすると CDF が得られ、1-p と predict(p) をプロットすると生存曲線 (1-CDF) が得られます。次のコードはより透過的で、p 値の任意のベクトルに一般化されています。

pct <- seq(.01,.99,by=.01)
lines(predict(sWei, newdata=list(sex=1),type="quantile",p=pct),1-pct,col="red")
lines(predict(sWei, newdata=list(sex=2),type="quantile",p=pct),1-pct,col="red")
于 2015-05-13T12:51:08.100 に答える