情节生存并使用曲线survreg的风险函数（）(Plot survival and hazard f

我有以下survreg模型：

Call:
survreg(formula = Surv(time = (ev.time), event = ev) ~ age, 
    data = my.data, dist = "weib")
             Value Std. Error    z        p
(Intercept) 4.0961     0.5566 7.36 1.86e-13
age         0.0388     0.0133 2.91 3.60e-03
Log(scale)  0.1421     0.1208 1.18 2.39e-01
Scale= 1.15 

Weibull distribution

我想绘制风险函数，并根据上述估计生存函数。
我不希望使用predict()或pweibull()这里所提出参数生存或这里SO问题 。

我想用curve()函数。 任何想法如何，我可以做到这一点？ 看来survreg的威布尔函数使用比例和形状比通常的（和不同之处在于例如rweibull）的其他定义。

更新：我想我真的需要它来表达危害/存活的估计的函数Intercept ， age (+ other potential covariates) ， Scale不使用任何现成的*weilbull功能。

Answer 1:

您的参数是：

scale=exp(Intercept+beta*x)在你的榜样，让我们说年龄= 40

scale=283.7

您的形状参数为尺度的倒数，该模型输出

shape=1/1.15

那么危害是：

curve((shape/scale)*(x/scale)^(shape-1), from=0,to=12,ylab=expression(hat(h)(t)), col="darkblue",xlab="t", lwd=5)

累积风险函数为：

curve((x/scale)^(shape), from=0,to=12,ylab=expression(hat(F)(t)), col="darkgreen",xlab="t", lwd=5)

生存函数为1-累积风险函数，因此：

curve(1-((x/scale)^(shape)), from=0,to=12,ylab=expression(hat(S)(t)), col="darkred",xlab="t", lwd=5, ylim=c(0,1))

还检查了eha封装，功能hweibull和Hweibull

Answer 2:

在你提供的第一个链接实际上是对如何工作，有一个可爱例如，沿着理论解释清楚。 （感谢你为这个，这是一个很好的资源，我会在自己的工作中使用。）

要使用curve的功能，你就需要通过一些函数作为参数。 这是事实，该*weibull系列函数使用不同的参数化韦伯比survreg ，但它可以很容易地转化，为您解释第一个链接。 另外，从文档survreg ：

有多种方式来参数威布尔分布。 所述survreg函数它体内埋下在一般位置规模的家庭内，这是一个不同的参数比rweibull功能，并经常导致混乱。

  survreg's scale = 1/(rweibull shape) survreg's intercept = log(rweibull scale) 

这里是一个简单的转换的实现：

# The parameters
intercept<-4.0961
scale<-1.15

par(mfrow=c(1,2),mar=c(5.1,5.1,4.1,2.1)) # Make room for the hat.
# S(t), the survival function
curve(pweibull(x, scale=exp(intercept), shape=1/scale, lower.tail=FALSE), 
      from=0, to=100, col='red', lwd=2, ylab=expression(hat(S)(t)), xlab='t',bty='n',ylim=c(0,1))
# h(t), the hazard function
curve(dweibull(x, scale=exp(intercept), shape=1/scale)
      /pweibull(x, scale=exp(intercept), shape=1/scale, lower.tail=FALSE), 
      from=0, to=100, col='blue', lwd=2, ylab=expression(hat(h)(t)), xlab='t',bty='n')
par(mfrow=c(1,1),mar=c(5.1,4.1,4.1,2.1))

据我所知，你在你的答案，你不想使用提到pweibull功能，但我猜你不想使用它，因为它使用了不同的参数。 否则，你可以简单地写自己的版本pweibull使用该survreg的参数：

my.weibull.surv<-function(x,intercept,scale) pweibull(x,scale=exp(intercept),shape=1/scale,lower.tail=FALSE)
my.weibull.haz<-function(x,intercept,scale) dweibull(x, scale=exp(intercept), shape=1/scale) / pweibull(x,scale=exp(intercept),shape=1/scale,lower.tail=FALSE)

curve(my.weibull.surv(x,intercept,scale),1,100,lwd=2,col='red',ylim=c(0,1),bty='n')
curve(my.weibull.haz(x,intercept,scale),1,100,lwd=2,col='blue',bty='n')

---

正如我在评论中提到的，我不知道你为什么会这么做（除非这是作业），但你可以handcode pweibull和dweibull如果你喜欢：

my.dweibull <- function(x,shape,scale) (shape/scale) * (x/scale)^(shape-1) * exp(- (x/scale)^shape)
my.pweibull <- function(x,shape,scale) exp(- (x/scale)^shape)

这些定义来直出的?dweibull 。 现在只是包装的，慢的，未经测试的功能，而不是pweibull和dweibull直接。

文章来源: Plot survival and hazard function of survreg using curve()