nbstat {lrstat}R Documentation

Negative binomial rate ratio

Description

Obtains the number of subjects accrued, number of events, number of dropouts, number of subjects reaching the maximum follow-up, total exposure, and variance for log rate in each group, rate ratio, variance, and Wald test statistic of log rate ratio at given calendar times.

Usage

nbstat(
  time = NA_real_,
  rateRatioH0 = 1,
  allocationRatioPlanned = 1,
  accrualTime = 0L,
  accrualIntensity = NA_real_,
  piecewiseSurvivalTime = 0L,
  stratumFraction = 1L,
  kappa1 = NA_real_,
  kappa2 = NA_real_,
  lambda1 = NA_real_,
  lambda2 = NA_real_,
  gamma1 = 0L,
  gamma2 = 0L,
  accrualDuration = NA_real_,
  followupTime = NA_real_,
  fixedFollowup = 0L,
  nullVariance = 0L
)

Arguments

time

A vector of calendar times for data cut.

rateRatioH0

Rate ratio under the null hypothesis.

allocationRatioPlanned

Allocation ratio for the active treatment versus control. Defaults to 1 for equal randomization.

accrualTime

A vector that specifies the starting time of piecewise Poisson enrollment time intervals. Must start with 0, e.g., c(0, 3) breaks the time axis into 2 accrual intervals: [0, 3) and [3, Inf).

accrualIntensity

A vector of accrual intensities. One for each accrual time interval.

piecewiseSurvivalTime

A vector that specifies the starting time of piecewise exponential survival time intervals. Must start with 0, e.g., c(0, 6) breaks the time axis into 2 event intervals: [0, 6) and [6, Inf). Defaults to 0 for exponential distribution.

stratumFraction

A vector of stratum fractions that sum to 1. Defaults to 1 for no stratification.

kappa1

The dispersion parameter (reciprocal of the shape parameter of the gamma mixing distribution) for the active treatment group by stratum.

kappa2

The dispersion parameter (reciprocal of the shape parameter of the gamma mixing distribution) for the control group by stratum.

lambda1

The rate parameter of the negative binomial distribution for the active treatment group by stratum.

lambda2

The rate parameter of the negative binomial distribution for the control group by stratum.

gamma1

The hazard rate for exponential dropout, a vector of hazard rates for piecewise exponential dropout applicable for all strata, or a vector of hazard rates for dropout in each analysis time interval by stratum for the active treatment group.

gamma2

The hazard rate for exponential dropout, a vector of hazard rates for piecewise exponential dropout applicable for all strata, or a vector of hazard rates for dropout in each analysis time interval by stratum for the control group.

accrualDuration

Duration of the enrollment period.

followupTime

Follow-up time for the last enrolled subject.

fixedFollowup

Whether a fixed follow-up design is used. Defaults to 0 for variable follow-up.

nullVariance

Whether to calculate the variance for log rate ratio under the null hypothesis.

Details

The probability mass function for a negative binomial distribution with dispersion parameter \kappa_i and rate parameter \lambda_i is given by

P(Y_{ij} = y) = \frac{\Gamma(y+1/\kappa_i)}{\Gamma(1/\kappa_i) y!} \left(\frac{1}{1 + \kappa_i \lambda_i t_{ij}}\right)^{1/\kappa_i} \left(\frac{\kappa_i \lambda_i t_{ij}} {1 + \kappa_i \lambda_i t_{ij}}\right)^{y},

where Y_{ij} is the event count for subject j in treatment group i, and t_{ij} is the exposure time for the subject. If \kappa_i=0, the negative binomial distribution reduces to the Poisson distribution.

For treatment group i, let \beta_i = \log(\lambda_i). The likelihood for \{(\kappa_i, \beta_i):i=1,2\} can be written as

l = \sum_{i=1}^{2}\sum_{j=1}^{n_{i}} \{\log \Gamma(y_{ij} + 1/\kappa_i) - \log \Gamma(1/\kappa_i) + y_{ij} (\log(\kappa_i) + \beta_i) - (y_{ij} + 1/\kappa_i) \log(1+ \kappa_i \exp(\beta_i) t_{ij})\}.

It follows that

\frac{\partial l}{\partial \beta_i} = \sum_{j=1}^{n_i} \left\{y_{ij} - (y_{ij} + 1/\kappa_i) \frac{\kappa_i \exp(\beta_i) t_{ij}} {1 + \kappa_i \exp(\beta_i)t_{ij}}\right\},

and

-\frac{\partial^2 l}{\partial \beta_i^2} = \sum_{j=1}^{n_i} (y_{ij} + 1/\kappa_i) \frac{\kappa_i \lambda_i t_{ij}} {(1 + \kappa_i \lambda_i t_{ij})^2}.

The Fisher information for \beta_i is

E\left(-\frac{\partial^2 l}{\partial \beta_i^2}\right) = n_i E\left(\frac{\lambda_i t_{ij}} {1 + \kappa_i \lambda_i t_{ij}}\right).

In addition, we can show that

E\left(-\frac{\partial^2 l} {\partial \beta_i \partial \kappa_i}\right) = 0.

Therefore, the variance of \hat{\beta}_i is

Var(\hat{\beta}_i) = \frac{1}{n_i} \left\{ E\left(\frac{\lambda_i t_{ij}}{1 + \kappa_i \lambda_i t_{ij}}\right) \right\}^{-1}.

To evaluate the integral, we need to obtain the distribution of the exposure time,

t_{ij} = \min(\tau - W_{ij}, C_{ij}, T_{fmax}),

where \tau denotes the calendar time since trial start, W_{ij} denotes the enrollment time for subject j in treatment group i, C_{ij} denotes the time to dropout after enrollment for subject j in treatment group i, and T_{fmax} denotes the maximum follow-up time for all subjects. Therefore,

P(t_{ij} \geq t) = P(W_{ij} \leq \tau - t)P(C_{ij} \geq t) I(t\leq T_{fmax}).

Let H denote the distribution function of the enrollment time, and G_i denote the survival function of the dropout time for treatment group i. By the change of variables, we have

E\left(\frac{\lambda_i t_{ij}}{1 + \kappa_i \lambda_i t_{ij}} \right) = \int_{0}^{\tau \wedge T_{fmax}} \frac{\lambda_i}{(1 + \kappa_i \lambda_i t)^2} H(\tau - t) G_i(t) dt.

A numerical integration algorithm for a univariate function can be used to evaluate the above integral.

For the restricted maximum likelihood (reml) estimate of (\beta_1,\beta_2) subject to the constraint that \beta_1 - \beta_2 = \Delta, we express the log-likelihood in terms of (\beta_2,\Delta,\kappa_1,\kappa_2), and takes the derivative of the log-likelihood function with respect to \beta_2. The resulting score equation has asymptotic limit

E\left(\frac{\partial l}{\partial \beta_2}\right) = s_1 + s_2,

where

s_1 = n r E\left\{\lambda1_1 t_{1j} - \left(\lambda_1t_{1j} + \frac{1}{\kappa_1}\right) \frac{\kappa_1 e^{\tilde{\beta}_2 + \Delta}t_{1j}}{1 + \kappa_1 e^{\tilde{\beta}_2 +\Delta}t_{1j}}\right\},

and

s_2 = n (1-r) E\left\{\lambda_2 t_{2j} - \left(\lambda_2 t_{2j} + \frac{1}{\kappa_2}\right) \frac{\kappa_2 e^{\tilde{\beta}_2} t_{2j}} {1 + \kappa_2 e^{\tilde{\beta}_2}t_{2j}}\right\}.

Here r is the randomization probability for the active treatment group. The asymptotic limit of the reml of \beta_2 is the solution \tilde{\beta}_2 to E\left(\frac{\partial l}{\partial \beta_2}\right) = 0.

Value

A list with two components:

Author(s)

Kaifeng Lu, kaifenglu@gmail.com

Examples

# Example 1: Variable follow-up design

nbstat(time = c(1, 1.25, 2, 3, 4),
       accrualIntensity = 1956/1.25,
       kappa1 = 5,
       kappa2 = 5,
       lambda1 = 0.7*0.125,
       lambda2 = 0.125,
       gamma1 = 0,
       gamma2 = 0,
       accrualDuration = 1.25,
       followupTime = 2.75)

# Example 2: Fixed follow-up design

nbstat(time = c(0.5, 1, 1.5, 2),
       accrualIntensity = 220/1.5,
       stratumFraction = c(0.2, 0.8),
       kappa1 = 3,
       kappa2 = 3,
       lambda1 = c(0.5*8.4, 0.6*10.5),
       lambda2 = c(8.4, 10.5),
       gamma1 = 0,
       gamma2 = 0,
       accrualDuration = 1.5,
       followupTime = 0.5,
       fixedFollowup = 1,
       nullVariance = 1)
[Package lrstat version 0.2.9 Index]