Instantly share code, notes, and snippets.

# AustinRochford/Bayesian Survival analysis with PyMC3.ipynb

Last active April 10, 2022 05:02
Show Gist options
• Save AustinRochford/afe6862e622c31494b2f to your computer and use it in GitHub Desktop.
Bayesian Survival analysis with PyMC3
Display the source blob
Display the rendered blob
Raw
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.

### jacoblevine commented Sep 30, 2017

First of all, thank you @AustinRochford for a wonderful demo.

Like others, I also had trouble reproducing the result for the first hazard regression (i.e. not the time-varying coefficient model, which I found to be reproducible). I was able to get similar results with a different model specification that is perhaps more typical for Bayesian regression:

```with pm.Model() as model:
lambda0 = pm.Gamma('lambda0', 0.01, 0.01, shape=n_intervals)

tau = pm.Gamma('tau', 10., 10.)
mu_beta = pm.Normal('mu_beta', 0., 10 ** -2)
beta = pm.Normal('beta', mu_beta, tau)

lambda_ = pm.Deterministic('lambda_', T.outer(T.exp(beta * df.metastized), lambda0))
mu = pm.Deterministic('mu', exposure * lambda_)

obs = pm.Poisson('obs', mu, observed=death)```

The gamma prior on `tau` produces a distribution for `beta` that looks reasonable for a regression model (i.e. centered on zero with a fair amount of density between -2 and 2).

Since `beta` is generated by a Gaussian random walk with fixed `tau=1` in the time-varying model, this can explain why someone running this demo could have a problem with the first but not second example. Incidentally, a gamma distribution with parameters (a=10, b=10) as in the code written above produces a distribution for `tau` with a mean of 1.

Of course, what I can't explain is why the model specification as it appears in the notebook worked in the first place. Anyway, hope this helps anyone else struggling with it.

And FWIW I'm using Python 3.6.1 on Mac OS X, pymc3 3.1, Theano 0.9.0, numpy 1.12.1