Statistical Specifications

Survival Model Specification

This package can only be used to fit proportional hazards models of the form:

$$\log(h_i(t \mid \theta, \psi_i)) = \log(h_0(t \mid \theta)) + X_i \beta + \sum_j G_j(t \mid \psi_i)$$

Where:

• $h_0(.)$ is a parametric baseline hazard function
• $t$ is the event time
• $\theta$ is a vector of parameters that parameterise the baseline hazard
• $\psi_i$ is an arbitrary vector of parameters for subject $i$ specified by the longitudinal model
• $G_j(.)$ is a link function that maps $\psi_i$ to a contribution to the log-hazard function where $j$ simply indexes the given link function
• $X_i$ is the subjects covariate design matrix
• $\beta$ is the corresponding coefficients to scale the design matrix covariates contribution to the log-hazard function

The following sections outline the available distributions to users which can be selected for the baseline hazard $h_0(.)$. Please note that some of these distributions do not have the proportional-hazards property meaning that the resulting survival model corresponding to the hazard $h_i()$ will not be of the same parametric family as the baseline distribution with the hazard $h_0(.)$.

Exponential Distribution

$$h(t \mid \lambda) = \lambda$$

Where: - $\lambda > 0$ is the rate parameter

Weibull Distribution (Proportional Hazard Parameterisation)

$$h(t \mid \lambda, \gamma) = \lambda \gamma t^{\gamma - 1 };$$

Where: - $\lambda > 0$ is the rate parameter - $\gamma > 0$ is the shape parameter Note that with $\gamma = 1$ we obtain the exponential distribution as a special case.

Log-Logistic Distribution

$$h(t \mid a, b) = \frac {(b/a)(t/a)^{(b-1)}} {1 + (t/a)^b}$$

Where:

• $a > 0$ is the scale parameter
• $b > 0$ is the shape parameter

Gamma Distribution

$$h(t \mid k, \theta) = \frac{ t^{k - 1} e^{-t/\theta} } { \theta ^ k ( \Gamma(k) - \gamma(k, t/\theta) ) }$$

Where:

• $k > 0 $ is the shape parameter
• $\theta > 0$ is the scale parameter
• $\Gamma(.)$ is the complete gamma function
• $\gamma(., .)$ is the lower incomplete gamma function

Flat Prior

prior_flat() instead contributes a constant density on the declared parameter scale. It is improper on an unbounded interval and is not invariant under a nonlinear change of scale. Its initialization distribution contributes no density term.

Priors on Transformed Parameters

Let $x$ denote a declared model parameter and $y = g(x)$ the scale on which prior information is available. A transformed prior has density

$$p_x(x) = \frac{p_y(g(x))}{F_y(g(u)) - F_y(g(l))}\, |g'(x)|, \qquad l < x < u,$$

where $(l,u)$ is the admissible parameter interval, including distribution support. The supported transformations are increasing and one-to-one on their domains. jmpost evaluates truncation bounds on the $y$ scale and adds the log absolute Jacobian to the Stan target:

Transformation	Parameter domain	Log absolute Jacobian
prior_transform_exp()	$x \in \mathbb{R}$	$x$
prior_transform_exp(scale = a), $a>0$	$x \in \mathbb{R}$	$\log a + ax$
prior_transform_square()	$x>0$	$\log 2 + \log x$

For example, prior_invgamma(1.5, 0.3, transform = prior_transform_exp(2)) places the inverse-gamma prior on $\exp(2x)$. If $x$ is a log standard deviation, this is a prior on the variance. It differs from placing an inverse-gamma prior directly on $x$. See the quickstart for examples using residual variances and log-scale intercepts.

The transform argument specifies the prior scale independently of the predictor parametrization. It includes a change-of-variable Jacobian; changing sampling coordinates with an affine multiplier preserves the specified prior distribution.

Covariate Selection with Regularized Horseshoe Prior

In the survival model the covariate coefficients $\beta$ typically just have a normal prior distribution, using prior_normal() or prior_normal_vector(). However, if the user wishes to perform covariate selection then they can instead use the regularized horseshoe prior via prior_horseshoe(). Here we explain how the prior is defined and how the hyperparameters can be set. Note that we use the same parametrization and notation as in the brms package (Bürkner 2017).

The horseshoe prior is a global-local shrinkage prior (Carvalho, Polson, and Scott 2010). The global parameter shrinks all coefficients towards zero, while coefficient-specific local parameters allow some coefficients to escape the shrinkage. The regularized horseshoe prior (Piironen and Vehtari 2017) additionally introduces a slab component that regularizes coefficients that are far from zero. This is useful for covariate selection because all candidate covariates can be included in one model fit, with irrelevant coefficients strongly shrunk towards zero.

For $p$ covariate coefficients $\beta_{j}, j = 1, \dotsc, p$, the regularized horseshoe prior applied with prior_horseshoe(), and with positive hyperparameters $\nu_{\lambda}, \nu_{\tau}, \nu_c, s_{\tau}, s_c$, can be written as

$$\begin{align*} z_j &\sim \operatorname{Normal}(0, 1), \\ \lambda_j &\sim t^+_{\nu_{\lambda}}(0, 1), \\ \tau &\sim t^+_{\nu_{\tau}}(0, s_{\tau}), \\ \xi &\sim \text{Inv-Gamma}\!\left(\frac{\nu_c}{2}, \frac{\nu_c}{2}\right), \\ c^2 &= s_c^2 \xi, \\ \tilde{\lambda}_j &= \left( \frac{c^2\lambda_{j}^2} {c^2 + \tau^2\lambda_{j}^2} \right)^{1/2}, \\ \beta_j &= z_j\tau\tilde{\lambda}_j, \qquad j = 1, \dotsc, p . \end{align*}$$

Here $t^+_{\nu}(0, s)$ denotes a half-Student-$t$ distribution with degrees of freedom $\nu$ and scale $s$. The $p$ local shrinkage parameters $\lambda_j$ control the amount of shrinkage for each coefficient, while the global shrinkage parameter $\tau$ controls the overall amount of shrinkage across all coefficients. The slab scale $s_c$ and slab degrees of freedom $\nu_c$ determine how strongly large coefficients are regularized.

Equivalently, conditional on $\lambda_j$, $\tau$, and $c$, this implies

$$\beta_j \mid \lambda_j, \tau, c \sim \operatorname{Normal}\!\left( 0,\, \tau^2 \frac{c^2\lambda_j^2}{c^2 + \tau^2\lambda_j^2} \right).$$

The arguments of prior_horseshoe() correspond to these hyperparameters:

• df: local degrees of freedom, $\nu_{\lambda}$
• df_global: global degrees of freedom, $\nu_{\tau}$
• df_slab: slab degrees of freedom, $\nu_c$
• scale_global: global scale, $s_{\tau}$
• scale_slab: slab scale, $s_c$

As a special case, setting $\nu_{\lambda} = 1$ and $\nu_{\tau} = 1$ gives half-Cauchy priors for the local and global shrinkage parameters. The slab part has the interpretation

$$c^2 \sim \text{Inv-Gamma}\!\left( \frac{\nu_c}{2}, \frac{\nu_c s_c^2}{2} \right),$$

which corresponds to a Student-$t$ slab with $\nu_c$ degrees of freedom and scale $s_c$ for coefficients that escape the main shrinkage. Smaller values of $s_{\tau}$ imply stronger prior shrinkage towards zero, while larger values allow more coefficients to remain non-zero. The slab scale should be chosen large enough to permit clinically or scientifically plausible covariate effects, but not so large that unrealistically large effects are left essentially unregularized.

A useful summary of the amount of shrinkage applied to coefficient $j$ is the shrinkage factor

$$\kappa_j = \frac{1}{1 + \tau^2\lambda_j^2 / c^2}.$$

Values of $\kappa_j$ close to 1 indicate strong shrinkage towards zero, while values close to 0 indicate little shrinkage. Regularized horseshoe priors have been applied in PK/PD covariate selection by Pourzanjani and Davis (2026), where the main practical advantage is the ability to assess all candidate covariates in a single Bayesian model fit.

Longitudinal Model Specification

The longitudinal models in jmpost have each two versions: one without covariates, and one with covariates. In each of the below subsections, we describe both versions. The models without covariates are simpler, and comprise parameters which are study- and treatment-arm-specific and correspondingly indexed. To start with, we describe the commonalities of the longitudinal models having covariates.

Models with Covariates

In the longitudinal models with covariates, the study- and treatment-arm-specific parameters of the longitudinal models without covariates are replaced by subject-specific predictors. For a model parameter $\theta$, let $\boldsymbol{x}_{\theta i}$ be the baseline covariates for subject $i$ and let $\boldsymbol{\beta}_{\theta}$ be their coefficients, while $\boldsymbol{\alpha}_{\theta}$ is the intercept. The model parameter value for subject $i$ may then be derived from the intercept and slope parameters together with the covariates by one of the following forms:

$$\theta_{i} = \begin{cases} \alpha_{\theta} + \boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta} & \text{(linear)}, \\ \alpha_{\theta}\left(1 + \boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta}\right) & \text{(proportional)}, \\ \alpha_{\theta}^{\,\boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta}} & \text{(exponential)}, \\ \exp\!\left(\alpha_{\theta} + \boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta}\right) & \text{(log-linear)}. \end{cases}$$

Standard-deviation parameters can only use the exponential or log-linear form, because these are strictly positive. The exponential form is $\alpha_{\theta}^{\,\boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta}}$ with $\alpha_{\theta}>0$, and the log-linear form is $\exp(\alpha_{\theta}+\boldsymbol{x}_{\theta i}^{\top}\boldsymbol{\beta}_{\theta})$. The exponential form fixes the reference-profile standard deviation at one. When selecting it, the intercept prior must support a positive base; the default scale-intercept priors are chosen for the log-linear form.

A couple of additional points worth calling out:

• Different covariates and parameterisations may be used for each $\theta$.
• Note that we use here location-scale models, where the mean and standard deviation of the response variable are modelled separately.
• Regarding the detailed implementation in the Stan code, we have introduced the options to switch between centred and non-centred parameterizations for the random effects. By default the baseline random effects use the centred parameterization, because it has been found to generally improve the sampling efficiency and convergence of the model for the baseline parameters. On the other hand, the random effects of all other parameters by default use the non-centred parameterization. See the arguments of the constructors for the respective longitudinal models with covariates for details, for example ?LongitudinalSteinFojoCov describes the arguments centred_baseline, centred_shrinkage and centred_growth.

Random-Slope Model

$$y_{ij} = \mu_{l(i)} + s_i t_{ij} + \epsilon_{ij}$$

where:

• $y_{ij}$ is the tumour size for subject $i$ at timepoint $j$
• $\mu_{l(i)}$ is the intercept for subject $i$
• $s_i$ is the random slope for subject $i$ with $s_i \sim N(\mu_{sk(i)}, \sigma_{sk(i)})$
• $\mu_{sk(i)}$ is the mean for the random slope within treatment arm $k(i)$
• $\sigma_{sk(i)}$ is the variance term for the slope in treatment arm $k(i)$
• $\epsilon_{ij}$ is the error term with $\epsilon_{ij} \sim N(0, \sigma)$
• $k(i)$ is the treatment arm index for subject $i$
• $l(i)$ is the study index for subject $i$

Covariate Version

Following models with covariates, the random-slope model then becomes

$$\begin{align*} y_{ij} &= \mu_i + s_i t_{ij} + \epsilon_{ij}, \\ s_i &\sim N\!\left(\mu_{s i},\ \sigma_{s i}\right). \end{align*}$$

Any predictor form above can be used for $\mu_i$. The mean $\mu_{s i}$ and standard deviation $\sigma_{s i}$ of the random slope standard deviation $s_i$ use either the exponential or log-linear parametrization. For example, under the log-linear parametrization, $\log(\sigma_{s i}) = \alpha_{\sigma_s} + \boldsymbol{x}_{\sigma_s i}^{\top}\boldsymbol{\beta}_{\sigma_s}$.

The original model is recovered by using study indicators for the covariates of $\mu_i$ and treatment-arm indicators for the the covariates of $\mu_{s i}$ and $\sigma_{s i}$.

Derivative of the SLD Trajectory Link / Growth Parameter Link

$$G(t_{ij} \mid \mu_0, s_i) = s_i$$

Accessible via linkDSLD() & linkGrowth()

Identity Link

$$\begin{align*} G(t_{ij} \mid \mu_0, s_i) = \mu_0 + s_i t_{ij} \end{align*}$$

Accessible via linkIdentity()

Stein-Fojo Model

$$\begin{align*} y_{ij} &\sim N(SLD_{ij},\ SLD_{ij}^2 \sigma^2) \\ \\ SLD_{ij} &= \begin{cases} b_i[e^{-s_it_{ij}} + e^{g_i t_{ij}} - 1] & \text{if } t_{ij}\geq 0 \\ b_i e^{g_i t_{ij}} & \text{if } t_{ij}\lt 0 \end{cases}\\\\ b_i &\sim \text{LogNormal}(\mu_{bl(i)}, \omega_{b l(i)}) \\ s_i &\sim \text{LogNormal}(\mu_{sk(i)}, \omega_{s k(i)}) \\ g_i &\sim \text{LogNormal}(\mu_{gk(i)}, \omega_{g k(i)}) \\ \end{align*}$$

Where:

• $i$ is the subject index
• $j$ is the visit index
• $y_{ij}$ is the observed tumour measurements
• $SLD_{ij}$ is the expected sum of longest diameter for subject $i$ at time point $j$
• $t_{ij}$ is the time since first treatment for subject $i$ at visit $j$
• $b_i$ is the subject baseline SLD measurement
• $s_i$ is the subject kinetics shrinkage parameter
• $g_i$ is the subject kinetics tumour growth parameter
• $k(i)$ is the treatment arm index for subject $i$
• $l(i)$ is the study index for subject $i$
• $\mu_{\theta k(i)}$ is the population mean for parameter $\theta$ in group $k(i)$
• $\omega_{\theta k(i)}$ is the population variance for parameter $\theta$.

If using the non-centred parameterisation then the following alternative formulation is used:

$$\begin{align*} b_i &= \exp(\mu_{bl(i)} + \omega_{b l(i)} * \eta_{b i}) \\ s_i &= \exp(\mu_{sk(i)} + \omega_{s k(i)} * \eta_{s i}) \\ g_i &= \exp(\mu_{gk(i)} + \omega_{g k(i)} * \eta_{g i}) \\ \\ \eta_{b i} &\sim N(0, 1)\\ \eta_{s i} &\sim N(0, 1) \\ \eta_{g i} &\sim N(0, 1) \\ \end{align*}$$

Where:

• $\eta_{\theta i}$ is a random effects offset on parameter $\theta$ for subject $i$

If using the unscaled variance parameterisation then the following alternative formulation is used: $$y_{ij} \sim N(SLD_{ij},\sigma^2)$$

Covariate Version

Following models with covariates, here we let $\theta_{i}$ be a subject-specific predictor for parameter $\theta$, constructed from parameter-specific baseline covariates. The standard-deviation predictors $\omega_{\theta i}$ can only use the exponential or log-linear form.

That means all population parameters, $\mu_{b i}$, $\omega_{b i}$, $\mu_{s i}$, $\omega_{s i}$, $\mu_{g i}$, and $\omega_{g i}$, are constructed from subject-specific predictors.

The original model is recovered by using study indicators for the predictors of $b_i$ and $\omega_{b i}$ and treatment-arm indicators for those of $s_i$ and $\omega_{s i}$ as well as those of $g_i$ and $\omega_{g i}$.

Derivative of the SLD Trajectory Link

$$G(t_{ij} \mid b_i, s_i, g_i) = \begin{cases} b_i(g_i e^{g_i t_{ij}} -s_ie^{-s_it_{ij}} ) & \text{if } t_{ij}\geq 0 \\ b_i g_i e^{g_i t_{ij}} & \text{if } t_{ij}\lt 0 \end{cases}$$

Accessible via linkDSLD()

Time to Growth Link

$$G(t_{ij} \mid b_i, s_i, g_i) = \max \left( \frac{ \text{log}(s_i) - \text{log}(g_i) }{ s_i + g_i }, 0 \right)$$

Accessible via linkTTG()

Identity Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i) &= SLD_{ij} \end{align*}$$

Accessible via linkIdentity()

Growth Parameter Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i) &= log(g_i) \end{align*}$$

Accessible via linkGrowth()

Shrinkage Parameter Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i) &= log(s_i) \end{align*}$$

Accessible via linkShrinkage()

Population Quantities

When generating population quantities for the $b$, $s$, and $g$ parameters, the median of the distribution is used. This has the same interpretation as using a non-centred parameterisation and setting the random effect to zero.

Generalized Stein-Fojo (GSF) Model

$$\begin{align*} y_{ij} &\sim N(SLD_{ij},\ SLD_{ij}^2 \sigma^2) \\ \\ SLD_{ij} &= \begin{cases} b_i[\phi_i e^{-s_it_{ij}} + (1-\phi_i)e^{g_i t_{ij}}] & \text{if } t_{ij}\geq 0 \\ b_i e^{g_i t_{ij}} & \text{if } t_{ij}\lt 0 \end{cases}\\ \\ b_i &\sim \text{LogNormal}(\mu_{bl(i)}, \omega_{b l(i)}) \\ s_i &\sim \text{LogNormal}(\mu_{sk(i)}, \omega_{s k(i)}) \\ g_i &\sim \text{LogNormal}(\mu_{gk(i)}, \omega_{g k(i)}) \\ \phi_i &\sim \text{LogitNormal}(\mu_{\phi k(i)}, \omega_{\phi k(i)}) \end{align*}$$

The logit-normal distribution is parameterised here as $\operatorname{logit}(\phi_i) \sim N(\mu_{\phi k(i)}, \omega_{\phi k(i)})$.

Where:

• $i$ is the subject index
• $j$ is the visit index
• $y_{ij}$ is the observed tumour measurements
• $SLD_{ij}$ is the expected sum of longest diameter for subject $i$ at time point $j$
• $t_{ij}$ is the time since first treatment for subject $i$ at visit $j$
• $b_i$ is the subject baseline SLD measurement
• $s_i$ is the subject kinetics shrinkage parameter
• $g_i$ is the subject kinetics tumour growth parameter
• $\phi_i$ is the subject proportion of cells affected by the treatment
• $k(i)$ is the treatment arm index for subject $i$
• $l(i)$ is the study index for subject $i$
• $\mu_{\theta k(i)}$ is the population location on the transformed scale for parameter $\theta$ in group $k(i)$
• $\omega_{\theta k(i)}$ is the population standard deviation on the transformed scale for parameter $\theta$ in group $k(i)$.

If using the non-centred parameterisation then the following alternative formulation is used:

$$\begin{align*} b_i &= \exp(\mu_{bl(i)} + \omega_{b l(i)} * \eta_{b i}) \\ s_i &= \exp(\mu_{sk(i)} + \omega_{s k(i)} * \eta_{s i}) \\ g_i &= \exp(\mu_{gk(i)} + \omega_{g k(i)} * \eta_{g i}) \\ \phi_i &= \text{logit}^{-1}(\mu_{\phi k(i)} + \omega_{\phi k(i)} * \eta_{\phi i}) \\ \\ \eta_{b i} &\sim N(0, 1)\\ \eta_{s i} &\sim N(0, 1) \\ \eta_{g i} &\sim N(0, 1) \\ \eta_{\phi i} &\sim N(0, 1) \\ \end{align*}$$

Where:

• $\eta_{\theta i}$ is a random effects offset on parameter $\theta$ for subject $i$

If using the unscaled variance parameterisation then the following alternative formulation is used: $$y_{ij} \sim N(SLD_{ij}, \sigma^2)$$

Covariate Version

Following models with covariates, here we let $\theta_{i}$ be a subject-specific predictor for parameter $\theta$, constructed from parameter-specific baseline covariates. The standard-deviation predictors $\omega_{\theta i}$ can only use the exponential or log-linear form.

That means all population parameters, $\mu_{b i}$, $\omega_{b i}$, $\mu_{s i}$, $\omega_{s i}$, $\mu_{g i}$, $\omega_{g i}$, $\mu_{\phi i}$, and $\omega_{\phi i}$, are constructed from subject-specific predictors.

The original model is recovered by using study indicators for the predictors of $b_i$ and $\omega_{b i}$ and treatment-arm indicators for those of $s_i$ and $\omega_{s i}$, $g_i$ and $\omega_{g i}$, and $\phi_i$ and $\omega_{\phi i}$.

Derivative of the SLD Trajectory Link

$$G(t_{ij} \mid b_i, s_i, g_i, \phi_i) = \begin{cases} b_i(-s_i\phi_ie^{-s_it_{ij}} + (1 - \phi_i)g_i e^{g_i t_{ij}}) & \text{if } t_{ij}\geq 0 \\ b_i g_i e^{g_i t_{ij}} & \text{if } t_{ij}\lt 0 \end{cases}$$

Accessible via linkDSLD()

Time to Growth Link

$$G(t_{ij} \mid b_i, s_i, g_i, \phi_i) = \max \left( \frac{ \log(s_i) + \log(\phi_i) - \log(g_i) - \log(1-\phi_i) }{ s_i + g_i }, 0 \right)$$

Accessible via linkTTG()

Identity Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i, \phi_i) &= SLD_{ij} \end{align*}$$

Accessible via linkIdentity()

Growth Parameter Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i, \phi_i) &= log(g_i) \end{align*}$$

Accessible via linkGrowth()

Shrinkage Parameter Link

$$\begin{align*} G(t_{ij} \mid b_i, s_i, g_i, \phi_i) &= log(s_i) \end{align*}$$

Accessible via linkShrinkage()

Population Quantities

Note that when generating population quantities for the $b$, $s$, $g$ and $\phi$ parameters the median of the distribution is used. This is because this has the same interpretation as using a non-centred parameterisation and setting the “random effects” term to be 0.

Claret-Bruno Model

$$\begin{align*} y_{ij} &\sim N(SLD_{ij},\ SLD_{ij}^2 \sigma^2) \\ \\ SLD_{ij} &= \begin{cases} b_i e^{g_i t_{ij}} & \text{if } t_{ij} < 0, \\ b_i \cdot \exp\left(g_i t_{ij} - \frac{p_i}{c_i} \left(1 - e^{-c_i t_{ij} }\right)\right) & \text{if } t_{ij} \geq 0. \end{cases}\\ \\ b_i &\sim \text{LogNormal}(\mu_{bl(i)}, \omega_{b l(i)}) \\ g_i &\sim \text{LogNormal}(\mu_{gk(i)}, \omega_{g k(i)}) \\ c_i &\sim \text{LogNormal}(\mu_{ck(i)}, \omega_{c k(i)}) \\ p_i &\sim \text{LogNormal}(\mu_{pk(i)}, \omega_{p k(i)}) \\ \end{align*}$$

Where:

• $i$ is the subject index
• $j$ is the visit index
• $y_{ij}$ is the observed tumour measurements
• $SLD_{ij}$ is the expected sum of longest diameter for subject $i$ at time point $j$
• $t_{ij}$ is the time since first treatment for subject $i$ at visit $j$
• $b_i$ is the subject baseline SLD value.
• $g_i$ is the subject tumour growth rate.
• $c_i$ is the subject treatment resistance rate.
• $p_i$ is the subject treatment growth inhibition response.
• $k(i)$ is the treatment arm index for subject $i$
• $l(i)$ is the study index for subject $i$
• $\mu_{\theta k(i)}$ is the population mean for parameter $\theta$ in group $k(i)$
• $\omega_{\theta k(i)}$ is the population variance for parameter $\theta$ in group $k(i)$

If using the non-centred parameterisation then the following alternative formulation is used:

$$\begin{align*} b_i &= \exp(\mu_{b l(i)} + \omega_{b l(i)} * \eta_{b i}) \\ g_i &= \exp(\mu_{g k(i)} + \omega_{g k(i)} * \eta_{g i}) \\ c_i &= \exp(\mu_{c k(i)} + \omega_{c k(i)} * \eta_{c i}) \\ p_i &= \exp(\mu_{p k(i)} + \omega_{p k(i)} * \eta_{p i}) \\ \\ \eta_{b i} &\sim N(0, 1)\\ \eta_{g i} &\sim N(0, 1) \\ \eta_{c i} &\sim N(0, 1) \\ \eta_{p i} &\sim N(0, 1) \\ \end{align*}$$

If using the unscaled variance parameterisation then the following alternative formulation is used: $$y_{ij} \sim N(SLD_{ij}, \sigma^2)$$

Covariate Version

Following models with covariates, here we let $\theta_{i}$ be a subject-specific predictor for parameter $\theta$, constructed from parameter-specific baseline covariates. The standard-deviation predictors $\omega_{\theta i}$ can only use the exponential or log-linear form.

That means all population parameters, $\mu_{b i}$, $\omega_{b i}$, $\mu_{g i}$, $\omega_{g i}$, $\mu_{c i}$, $\omega_{c i}$, $\mu_{p i}$, and $\omega_{p i}$, are constructed from subject-specific predictors.

The original model is recovered by using study indicators for the predictors of $b_i$ and $\omega_{b i}$ and treatment-arm indicators for those of $g_i$ and $\omega_{g i}$, $c_i$ and $\omega_{c i}$, and $p_i$ and $\omega_{p i}$.

Derivative of the SLD Trajectory Link

$$G(t_{ij} \mid b_i, g_i, c_i, p_i) = \begin{cases} b_i g_i e^{g_i t_{ij}} & \text{if } t_{ij}\lt 0 \\ \left( g_i - p_i e^{-c_i t_{ij}} \right) SLD_{ij} & \text{if } t_{ij}\geq 0 \end{cases}$$

Accessible via linkDSLD()

Time to Growth Link

$$G(t_{ij} \mid b_i, g_i, c_i, p_i) = \max \left( \frac{ln(\frac{p_i}{g_i})}{c_i}, 0 \right)$$

Accessible via linkTTG()

Identity Link

$$\begin{align*} G(t_{ij} \mid b_i, g_i, c_i, p_i) &= SLD_{ij} \end{align*}$$

Accessible via linkIdentity()

Growth Parameter Link

$$\begin{align*} G(t_{ij} \mid b_i, g_i, c_i, p_i) &= log(g_i) \end{align*}$$

Accessible via linkGrowth()

Population Quantities

When generating population quantities for the $b$, $g$, $c$ and $p$ parameters the median of the distribution is used. This is because this has the same interpretation as using a non-centred parameterisation and setting the “random effects” term to be 0.

Post-Processing

Brier Score

The Brier Score is used to measure a models predictive performance. In the case of Survival Models it is the weighted squared difference between whether a subject died within a given time interval and the models estimated probability of them dying within said interval. Within jmpost the following formula (as described in Blanche et al. (2015)) has been implemented:

$$\hat{BS}(t) = \frac{1}{n}\sum_{i=1}^{n} \hat{W}_i(t) \Big(D_i(t) - \pi_i(t) \Big)^2$$

$$\hat{W}_i(t) = \frac{\mathbb{1}_{(T_i \gt t)}}{\hat{G}(t)} + \frac{\mathbb{1}_{(T_i \le t)} \cdot \Delta_i}{\hat{G}(T_i)}$$

Where:

• $T_i$ is the observed event/censor time for subject $i$
• $\Delta_i$ is the event indicator which is 1 if $T_i$ is an event and 0 if $T_i$ is censored
• $\mathbb{1}_{(.)}$ is the indicator which is 1 if ${(.)}$ is true else 0
• $D_i(t)$ is the event indicator function for subject $i$ e.g. $D_i(t) = \mathbb{1}_{(T_i \lt t,\ \delta_i = 1)}$
• $\pi_i(t)$ is a model predicted probability of subject $i$ dying before time $t$
• $\hat{G}(u)$ is the Kaplan-Meier estimator of survival function of the censoring time at $u$

Note that by default $\hat{G}(T_i)$ is estimated by $\hat{G}(T_i-)$ and that in the case of ties event times are always considered to have occurred before censored times (this is in contrast to survival::survfit which regards censored times as coming before event times when estimating the censoring distribution). Both of these default options can be changed if required.

References

Blanche P, Proust-Lima C, Loubère L, Berr C, Dartigues J-F, Jacqmin-Gadda H (2015). ‘Quantifying and Comparing Dynamic Predictive Accuracy of Joint Models for Longitudinal Marker and Time-to-Event in Presence of Censoring and Competing Risks’. Biometrics, 71(1), 102–113. https://doi.org/https://doi.org/10.1111/biom.12232.

Bürkner P-C (2017). ‘brms: An R Package for Bayesian Multilevel Models Using Stan’. Journal of Statistical Software, 80(1), 1–28. https://doi.org/10.18637/jss.v080.i01.

Carvalho CM, Polson NG, Scott JG (2010). ‘The Horseshoe Estimator for Sparse Signals’. Biometrika, 97(2), 465–480. https://doi.org/10.1093/biomet/asq017.

Piironen J, Vehtari A (2017). ‘Sparsity Information and Regularization in the Horseshoe and Other Shrinkage Priors’. Electronic Journal of Statistics, 11(2), 5018–5051. https://doi.org/10.1214/17-EJS1337SI.

Pourzanjani A, Davis C (2026). ‘The Regularized Horseshoe for Covariate Selection Improves Convenience and Predictive Performance in Population PK/PD Models’. CPT: Pharmacometrics & Systems Pharmacology, 15(2), e70198. https://doi.org/10.1002/psp4.70198.