The PHREG Procedure

Proportional Rates/Means Models for Recurrent Events

Let be the number of events experienced by a subject over the time interval . Let be the increment of the counting process N over . The rate function is given by

$d\mu _{\bZ }(t)= E[dN(t) | Z(t)] = \mr {e}^{\bbeta \bZ (t)}d\mu _0(t)$

where $\mu _0(.)$ is an unknown continuous function. If the $\bZ$ are time independent, the rate model is reduced to the mean model

$\mu _ Z(t)= \mr {e}^{\bbeta Z} \mu _0(t)$

The partial likelihood for n independent triplets $(N_{i},Y_{i},\mb {Z}_{i}), i=1, \ldots , n$ , of counting, at-risk, and covariate processes is the same as that of the multiplicative hazards model. However, a robust sandwich estimate is used for the covariance matrix of the parameter estimator instead of the model-based estimate.

Let $T_{ki}$ be the kth event time of the ith subject. Let be the censoring time of the ith subject. The at-risk indicator and the failure indicator are, respectively,

$Y_ i(t) = I(C_ i\ge t) ~ ~ \mr {and}~ ~ \Delta _{ki}= I(T_{ki} \le C_ i)$

Denote

$S^{(0)}(\bbeta ,t) = \sum _{i=1}^ n Y_ i(t)\mr {e}^{\bbeta \bZ _ i(t)} ~ ~ \mr {and}~ ~ \bar{\bZ }(\bbeta ,t) = \frac{\sum _{i=1}^ n Y_ i(t)\mr {e}^{\bbeta \bZ _ i(t)} \bZ _ i(t)}{S^{(0)}(\bbeta ,t)}$

Let $\hat{\bbeta }$ be the maximum likelihood estimate of $\bbeta$ , and let $\mc {I}(\hat{\bbeta })$ be the observed information matrix. The robust sandwich covariance matrix estimate is given by

$\mc {I}^{-1}(\hat{\bbeta })\sum _{i=1}^ n \biggl [ W_ i(\hat{\bbeta }) W’_ i(\hat{\bbeta }) \biggl ] \mc {I}^{-1}(\hat{\bbeta })$

where

$\displaystyle W_ i(\bbeta )$	$\displaystyle =$	$\displaystyle \sum _ k \Delta _{ki}\biggl \{ Z_ i(T_{ki}) - \bar{\bZ }(\bbeta ,T_{ki}) \biggr \} -$
$\displaystyle$	$\displaystyle$	$\displaystyle \sum _{i=1}^ n \sum _ l \frac{\Delta _{lj}Y_ i(T_{lj}) \mr {e}^{\bbeta \bZ _ i(T_{lj})}}{S^{0}(\bbeta ,T_{lj}) } \biggl \{ Z_ i(T_{lj}) - \bar{\bZ }(\bbeta ,T_{lj}) \biggr \}$

For a given realization of the covariates $\bxi$ , the Nelson estimator is used to predict the mean function

$\hat{\mu }_{\bxi }(t) = \mr {e}^{\hat{\bbeta }\bxi } \sum _{i=1}^ n \sum _ k \frac{I(T_{ki} \le t) \Delta _{ki}}{S^{(0)}(\hat{\bbeta },T_{ki})}$

with standard error estimate given by

$\hat{\sigma }^2(\hat{\mu }_{\bxi }(t)) = \sum _{i=1}^ n \biggl (\frac{1}{n}\hat{\Psi }_ i(t,\bxi ) \biggr )^2$

where

$\displaystyle \frac{1}{n}\hat{\Psi }_ i(t,\bxi )$	$\displaystyle =$	$\displaystyle \mr {e}^{\hat{\bbeta }\bxi } \biggl \{ \sum _ k \frac{I(T_{ki}\le t)\Delta _{ik}}{S^{(0)}(\hat{\bbeta },T_{ki})} - \sum _{j=1}^ n\sum _ k \frac{Y_ i(T_{kj}) \mr {e}^{\hat{\bbeta }\bZ _ i(T_{kj})}I(T_{kj} \le t)\Delta _{kj}}{[S^{(0)}(\hat{\bbeta },T_{kj})]^2} -$
$\displaystyle$	$\displaystyle$	$\displaystyle \biggl [ \sum _{i=1}^ n \sum _ k \frac{I(T_{ki}\le t)\Delta _{ik} [\bar{\bZ }(\hat{\bbeta },T_{ki}) - \bxi ]}{S^{(0)}(\hat{\bbeta },T_{ki})} \biggr ]$
$\displaystyle$	$\displaystyle$	$\displaystyle \times \mc {I}^{-1}(\hat{\bbeta })\int _0^{\tau } [\bZ _ i(u) - \bar{\bZ }(\hat{\bbeta },u)]d\hat{M}_ i(u)\biggl \}$

Since the cumulative mean function is always nonnegative, the log transform is used to compute confidence intervals. The $100(1-\alpha )$ % pointwise confidence limits for $\mu _{\bxi }(t)$ are

$\hat{\mu }_{\bxi }(t) \mr {e}^{\pm z_{\alpha /2} \frac{\hat{\sigma }(\hat{\mu }_{\bxi }(t)) }{\hat{\mu }{\bxi }(t)}}$

where $z_{\alpha /2}$ is the upper $100\alpha /2$ percentage point of the standard normal distribution.