*! version 1.0 January 15, 2020 @ 13:36:49

*! version 1.1 September 9, 2018 @ 18:13:10

Title

wtdttt -- Estimate maximum likelihood estimate for
parametric Waiting Time Distribution (WTD) based on observed
prescription redemptions with adjustment for covariates. Reports
estimates of prevalence fraction and specified percentile of
inter-arrival density together with regression coefficients.

Syntax

wtdttt
varname
, disttype(recurrence_dens) start(time) end(time) [options]

Recurrence distribution
disttype(string) Parametric distribution for Forward or
Backward Recurrence Density (FRD/BRD)

Time window
reverse Estimate reverse WTD

start(string) Start of observation window

end(string) End of observation window

conttime Estimate durations for continuous data

Covariates
logitpcovar(fvvarlist) Covariates for logit(p)

mucovar(fvvarlist) Covariates for mu (lnorm)

lnsigmacovar(fvvarlist) Covariates for
log(sigma) (lnorm)

lnbetacovar(fvvarlist) Covariates for
log(beta) (exp | wei)

lnalphacovar(fvvarlist) Covariates for
log(alpha) (exp | wei)

allcovar(fvvarlist) Covariates for all parameters

Reporting
iadpercentile(#) Percentile to estimate in the
Inter-Arrival Distribution (IAD); default is 0.8 (80th percentile)

eform(string) Report exponentiated regression coefficients

Maximum likelihood options (all options available with ml)
niter(#) perform maximum of # iterations; default is niter(50)

Description

wtdttt estimates a parametric Waiting Time Distribution (WTD)
to varname and then computes an estimate of the specified percentile
together with an estimate of the proportion of prevalent users in the
sample. All parameters of the reverse WTD can be allowed to depend on
covariates.

Note: To use this command you first need to create a dataset
which contains only the first prescription redemption of each
individual within an observation window (ordinary WTD), or the last
(reverse WTD), respectively. You can typically achieve
this with the following two lines of code:

. keep if rxdate >= startdate & rxdate <= enddate

. bysort pid (rxdate): keep if _n == 1

or for the reverse WTD replace the last line with

. bysort pid (rxdate): keep if _n == _N

Here pid is a variable containing a person identifier
and rxdate is a variable containing the time of observed
prescription redemptions, typically dates. Before estimation it is
advisable to make a histogram of the variable rxdate to ensure
the distribution has the shape of a WTD, ordinary or reverse.

To assess the fit, the command wtdtttdiag can be used to
obtain diagnostic plots.

After estimation, the command wtdtttpredprob allows
prediction of treatment probabilities on times of interest based on
observed prescription redemptions. Similarly, the command
wtdtttpreddur allows prediction of prescription durations
based on a specified percentile of the inter-arrival distribution.

A key requirement for the WTD is that the index date is chosen
independently of the individual redemption processes. Instead of using
a fixed index date for the WTD and identical observation windows for
each individual, it is possible to use randomly sampled index dates.
This is implemented with the command ranwtdttt, which takes a
dataset of all prescription redemption dates of all users as input.
It then uniformly samples unique index dates for each individual, such
that only the first prescription subsequent to (ordinary WTD) - or the last
prescription prior to (reverse WTD) - the random index date for each
individual are considered in the WTD estimation.

Options

Recurrence distribution

disttype specifies the forward recurrence density to use.
Possible choices are named after their corresponding interarrival
density and there are three different choices implemented: exp
means Exponential, lnorm means Log-Normal, and wei means
Weibull. See Remarks below for a description of these and their
parametrization.

Time window

reverse indicates that observations represent the last
prescription redemption observed in the interval for each patient and
a reverse WTD is estimated. If not specified (default), observations
are assumed to be first prescription redemptions and the ordinary WTD
is estimated.

start is either a string such as "1Jan2014" (a date for discrete data)
or "0" (a number for continuous data) which gives the start of the observation
window. In the discrete case (default) the string must conform to requirements
for the date function td().

end is either a string such as "31Dec2014" (a date for discrete data)
or "1" (a number for continuous data) which gives the end of the observation
window. In the discrete case (default) the string must conform to requirements
for the date function td().

conttime indicates that the data is continuous. If not
specified (default), observations are assumed to be discrete.

Covariates

logitpcovar(fvvarlist) specifies covariates included in
the regression equation for the parameter logit(p) (log-odds of
prevalence).

mucovar(fvvarlist) specifies covariates included in the
regression equation for mu, when using a Log-Normal recurrence
distribution (lnorm).

lnsigmacovar(fvvarlist) Covariates for log(sigma) (lnorm)

lnbetacovar(fvvarlist) Covariates for log(beta) (exp | wei)

lnalphacovar(fvvarlist) Covariates for log(alpha) (exp | wei)

allcovar(fvvarlist) Covariates included in all regression
equations for the parameters.

iadpercentile The percentile of the IAD, which is to be
estimated specified as a fraction between 0 and 1 (default is 0.8).

Remarks - Methods and Formulas

The WTD is parametrized as a two-component mixture distribution with
one density component for prevalent users, g(t), and a uniform
distribution over the observation window for incident users, i.e. the
likelihood contribution for one patient is:

l(t) = p * g(t) + (1 - p) / delta

where p is the proportion of prevalent users in the
sample, and delta is the width of the observation window. Only
patients with at least one prescription redemption in the observation
window are considered. When a reverse WTD is estimated, t is time
from last prescription in the time window to its end.

The density g(t) is known as the forward (ordinary WTD) or
backward (reverse WTD) recurrence density corresponding to the
interarrival density, f(t), which governs the distribution of
time from one prescription redemption of a patient to the subsequent
one. g(t) is given by

g(t) = (1 - F(t)) / M

where F(t) is the cumulative distribution function for
f(t) and M is the mean for f(t).

The actual parametrizations used are:

Exponenential:

f(t) = exp(-(exp(lnbeta) * t))

where lnbeta = ln(beta).

Weibull:

f(t) = exp(-(exp(lnbeta) * t)^exp(lnalpha))

where lnalpha = ln(alpha) and lnbeta = ln(beta).

Log-Normal:

f(t) = normprob(-(ln(x) - mu) / exp(lnsigma))

where lnsigma = ln(sigma).

The ML procedure reports estimates of lnbeta
(Exponential), (lnalpha, lnbeta) (Weibull) or (mu, lnsigma)
(Log-Normal) together with an estimate of the log-odds of prevalent
users, logitp. The latter is also reported as the estimated proportion
of prevalent users in the sample after an
inverse-logit-transformation, i.e. exp(logitp)/(1 + exp(logitp))
accompanied by a 95% confidence interval.

Note that for the percentile reported on the ordinary time-scale (not log time), the
reported confidence interval is z-based using the delta-method. However, for
smaller sample sizes the assumed asymptotic normality of this estimate is
questionable and instead the confidence interval on the log-scale should be
used and then be exponentially transformed to the original time scale (see the
example do-file for a worked example of this).

Examples

. wtdttt rx1time, disttype(lnorm) start(1jan2014) end(31dec2014)

. wtdttt last_rxtime, disttype(lnorm) start(0) end(1) conttime reverse mucovar(i.packsize) logitpcovar(i.packsize)

To get bootstrap confidence intervals we can do the following - notice
the use of eform in the second statement to obtain the
percentile itself and not its logarithm:

. bootstrap logtimeperc = r(logtimeperc), reps(50): ///

wtdttt rx1time, disttype(lnorm) start(1jan2014) end(31dec2014)

. bootstrap, eform

Further examples are provided in the example do-file
wtdttt_ex.do, which contains analyses based on the datafile
wtddat.dta - a simulated dataset, which is also enclosed. Be sure
to read comments in the do-file for further explanations.

Stored results

wtdttt stores the following scalars in r():

Scalars

r(logtimeperc) Logarithm of estimated IAD percentile

r(selogtimeperc) Standard error of estimated logarithm of IAD percentile

r(timepercentile) Estimated IAD percentile

r(setimepercentile) Standard error of estimated IAD percentile

r(prevprop) Estimated proportion of prevalent users

r(seprev) Standard error of estimated proportion of
prevalent users

r(disttype) Model type (backward or forward recurrence distribution)

r(reverse) If undefined: Ordinary WTD. If defined and
equal to "reverse": Reverse WTD.

r(delta) Length of observation window

r(start) Time value at start of time window

r(end) Time value at end of time window

Apart from the above, all results obtained by the maximum likelihood
estimation are stored by ml in the usual e() macros, see
help ml.

References

Støvring H, Pottegård A, Hallas J. Refining estimates of
prescription durations by using observed covariates in
pharmacoepidemiological databases: an application of the reverse
waiting time distribution. Pharmacoepidemiol Drug Saf. 2017.
doi:10.1002/pds.4216. http://dx.doi.org/10.1002/pds.4216

Støvring H, Pottegård A, Hallas J. Estimating medication
stopping fraction and real-time prevalence of drug use in
pharmaco-epidemiologic databases. An application of the reverse
waiting time distribution. Pharmacoepidemiol Drug Saf. 2017.
doi:10.1002/pds.4217. http://dx.doi.org/10.1002/pds.4217

Støvring H, Pottegård A, Hallas J. Determining prescription
durations based on the parametric waiting time distribution.
Pharmacoepidemiol Drug Saf. 2016. doi:10.1002/pds.4114. http://dx.doi.org/10.1002/pds.4114

Hallas J, Gaist D, Bjerrum L. The Waiting Time Distribution
as a Graphical Approach to Epidemiologic Measures of Drug Utilization.
Epidemiology. 1997;8:666-670. http://insights.ovid.com/pubmed?pmid=9345667

Stovring H, Vach W. Estimation of prevalence and incidence
based on occurrence of health-related events. Stat Med.
2005;24(20):3139-3154. https://doi.org/10.1002/sim.2142

Author

Katrine Bødkergaard Nielsen, Aarhus University, kani@clin.au.dk.

Henrik Støvring, Aarhus University, stovring@ph.au.dk.