intervals: $$\begin{aligned}
L\left(\mathbf{t}_{1:J}\vert\boldsymbol{\beta}_{1:J}\right) & = & \prod_{j=1}^{J}L_{j}\left(\mathbf{t}_{j}|\boldsymbol{\beta}_{j}\right),\label{eq:interval likelihood}\end{aligned}$$ where $$\begin{aligned}
L_{j}\left(\mathbf{t}_{j}|\boldsymbol{\beta}_{j}\right) & =\prod_{i=1}^{n_{j}}\lambda_{ij}^{d_{ij}}\exp\left(-\lambda_{ij}t_{ij}\right),\end{aligned}$$ $\mathbf{t}_{j}$ is the vector of exposures for interval $I_{j}$, $\mathbf{t}_{1:J}=(\mathbf{t}_{1},\ldots,\mathbf{t}_{J})$, $\boldsymbol{\beta}_{1:J}=(\beta_{1},\ldots,\beta_{J})$, and $n_{j}$ is the number of individuals who experienced the event during interval $I_{j}$. For more details and justification of the likelihood given in (\[eq:interval likelihood\]), see @gamerman1991dynamic.

Alternative likelihood expressions can be obtained via data augmentation [@wagner2011bayesian] or partial likelihood [@sargent1997flexible] approaches. The data augmentation approach completes the censored exposure times (corresponding to $d_{ij}=0$) by a latent exponentially distributed residual time. @wagner2011bayesian approximates the log of the augmented survival times by a mixture of ten normal components with known mean and variance parameters, resulting in a dynamic linear model. This method is used in Section \[subsec:Simulations\] as a benchmark for assessing the performance of the inference methodology proposed in this paper.

The PEH model relies on the method of partitioning the time into $J$ intervals. The partition can be based on the event times, placing $\tau_{j}$ at each observed event time [@gamerman1991dynamic], which may be computationally expensive especially when the number of events is large. Alternatively, intervals could be equidistant as in @hemming2002parametric, or they could contain equal number of events per interval (in this 