# Fitting time-to-event data with Weibull Hazard using brm function

**URL:** <https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638>\
**Category:** brms\
**Created:** [June 26, 2018, 10:24am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638 "2018-06-26T10:24:55Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![Ivan\_Azarov](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/ivan_azarov/32/2954_2.png) [@Ivan\_Azarov](https://discourse.mc-stan.org/u/Ivan_Azarov)\
**Post date:** [June 26, 2018, 10:24am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/1 "2018-06-26T10:24:55Z")

</div>

Good afternoon, colleagues!

I have recently been the user of R package ‘brms’.  
I’m trying to use it to fit time-to-event data about the survival of patients.  
I planned to use a Hazard function of this kind: the Weibull Hazard multiplied by the factor responsible for the effect of covariates.  
My question is this: did I understand correctly that in this package this fitting in brm will be recorded as follows?

> fit2 ← brm(time (time) | cens (CENSOR) ~ WHOSTATN + BASESLD + NLine,  
> data = df1, family = weibull(), inits = “0”,  
> control = list(adapt\_delta = 0.95, max\_treedepth = 12))

Is it enough, that the type of hazard is “contained” in the likelihood function?

* * *

- Operating System: Windows 7
- brms Version: 2.3.1

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 26, 2018, 10:30am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/2 "2018-06-26T10:30:11Z")

</div>

I don’t understand the `time (time)` specification. I would rather write:

```
fit2 <- brm(time | cens(CENSOR) ~ WHOSTATN + BASESLD + NLine,
            data = df1, family = weibull(), inits = “0”,
            control = list(adapt_delta = 0.95, max_treedepth = 12))

```

The rest seems sensible to me.

---

<div class="post-metadata">

**Author:** ![Ivan\_Azarov](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/ivan_azarov/32/2954_2.png) [@Ivan\_Azarov](https://discourse.mc-stan.org/u/Ivan_Azarov)\
**Post date:** [June 26, 2018, 10:45am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/3 "2018-06-26T10:45:24Z")

</div>

> [@paul.buerkner](#):
>
> The rest seems sensible to me.

Thank you, Paul!  
With such time specification it works well and reasonably. Compilation and calculations — all is ok!  
But without this specification I obtained such error message:

> Error: Family ‘weibull’ requires response greater than 0.

Maybe, such specification is the bad way, but without this trick compilation can not start.

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 26, 2018, 10:51am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/4 "2018-06-26T10:51:03Z")

</div>

Depending on the input, the `time()` transformation in R does something I am not sure you want to do. Check the output of `time(df1$time)`.

Rather, you should take the error message of brms seriously. Perhaps there is something in the `time` variable which is incorrect or you have `0` in it, which cannot be handled by such models. Alternatively, you may need to transform `time` to a numeric variable first via `as.numeric(time)`.

---

<div class="post-metadata">

**Author:** ![Ivan\_Azarov](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/ivan_azarov/32/2954_2.png) [@Ivan\_Azarov](https://discourse.mc-stan.org/u/Ivan_Azarov)\
**Post date:** [June 26, 2018, 11:54am UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/5 "2018-06-26T11:54:16Z")

</div>

Paul, you are absolutely right!  
I didn’t think, that time = 0 can not be used.  
Much grateful!

But I still have a question: why brm doesn’t give the results of estimating the two parameters in the Weibull risk function?  
Since there is not only the shape parameter k, but also the “lambda” is a scale parameter:  
(k/lambda)\*(x/lambda)^(k-1).

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 26, 2018, 12:00pm UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/6 "2018-06-26T12:00:22Z")

</div>

Please see [https://cran.r-project.org/web/packages/brms/vignettes/brms\_families.html](https://cran.r-project.org/web/packages/brms/vignettes/brms_families.html) for the parameterization of the weibull family (and other families) in brms.

---

<div class="post-metadata">

**Author:** ![Ivan\_Azarov](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/ivan_azarov/32/2954_2.png) [@Ivan\_Azarov](https://discourse.mc-stan.org/u/Ivan_Azarov)\
**Post date:** [June 26, 2018, 2:52pm UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/7 "2018-06-26T14:52:06Z")

</div>

Paul, this is precisely the parametrization of the Weibul distribution, to which I am accustomed. And in the summary we can see:

> Family: weibull  
> Links: mu = log; shape = identity  
> Formula: time | cens(CENSOR) ~ AGE  
> Data: df1 (Number of observations: 250)  
> Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;  
> total post-warmup samples = 4000
> 
> Population-Level Effects:  
> Estimate Est.Error l-95% CI u-95% CI Eff.Sample Rhat  
> Intercept 14.99 3.37 8.92 22.22 2021 1.00  
> AGE -0.04 0.05 -0.14 0.06 2143 1.00
> 
> Family Specific Parameters:  
> Estimate Est.Error l-95% CI u-95% CI Eff.Sample Rhat  
> shape 0.28 0.03 0.22 0.35 1813 1.00
> 
> Samples were drawn using sampling(NUTS). For each parameter, Eff.Sample  
> is a crude measure of effective sample size, and Rhat is the potential  
> scale reduction factor on split chains (at convergence, Rhat = 1).

We can see value of **shape parameter**.  
But, if I’m interested in looking at the resulting value of the **scale parameter** , how can it be derived?

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 26, 2018, 4:50pm UTC](https://discourse.mc-stan.org/t/fitting-time-to-event-data-with-weibull-hazard-using-brm-function/4638/8 "2018-06-26T16:50:01Z")

</div>

The parameterization is in terms of mean and shape not in terms of scale and shape.

You can extract the mean fitted values via `fitted(fit2)`. If you know the transformation formula to compute the scale based on mean and shape, you can use

```
ps_mean <- fitted(fit2, summary = FALSE)
ps_shape <- as.vector(as.matrix(fit2, pars = "^shape$"))
ps_scale <- trans(ps_mean, ps_shape)
posterior_summary(ps_scale)

```

where `trans()` is a vectorized function computing the scale based on mean and shape.