This version adds non-parametric delay distributions, both a direct PMF over fixed bins (step CDF) and a discrete-time hazard parameterisation, with support for fitting them via fitdistdoublecens() and pcd_cmdstan_model().
pcens_cdf() gains methods for flexsurv::pgengamma.orig() (Stacy parameterisation, always analytical) and flexsurv::pgengamma() (Prentice parameterisation, analytical for Q > 0 and numeric otherwise).
pprimarycensored() and related functions already worked with any pdist, including these, via numeric integration.
In Stan, dist_lcdf() now supports dist_id = 5 (generalised gamma, parameters [shape, scale, k]) through a new gengamma_lcdf() function and primarycensored_gengamma_uniform_lcdf() provides the analytical solution used by primarycensored_lpmf() and related functions, including in pcd_cmdstan_model().
flexsurv has been added to Suggests.flexsurv::pgengamma.orig()) are now identified by name when passed as pdist or dprimary, so analytical solutions are used without needing add_name_attribute().pdiscretestep(), ddiscretestep(), and rdiscretestep() represent a delay via a direct PMF over fixed bins, and the discrete-time hazard variant pdiscretehazard(), ddiscretehazard(), and rdiscretehazard() parameterises the same family via per-bin hazards, with conversion utilities hazards_to_pmf() and pmf_to_hazards(). See #218.fitdistdoublecens() accepts distr = "discretestep" and distr = "discretehazard" through the same code path as parametric distributions, with K (the number of bins) inferred from start. Hazard priors are user-settable via prior = list(...).fitdistdoublecens() the hazard_model argument picks between "rw" (random walk on the logit hazards, the default) and "re" (IID logit random effects around a mean intercept, logit(h_i) = alpha + sigma * eps_i with eps_i ~ N(0, 1)).discretehazard_start() for building the named start-value list expected by fitdistdoublecens() with distr = "discretehazard".dist_id = 26 is the step CDF with the PMF supplied in params and dist_id = 27 the step CDF with the hazards supplied in params. Both take params = c(boundaries, weights) of length 2 * K + 1, so dist_lcdf() and primarycensored_lpmf() evaluate them from params alone. dist_id = 28 shares the likelihood of 27 and only differs in the prior used by the package model. The primary event censored CDF is analytic for the uniform and exponential growth primaries, for any bin widths and for boundaries that start below zero.pcd_cmdstan_model() fits the non-parametric families directly: 26 with a Dirichlet prior on the PMF, 27 with a Gaussian random walk on the logit hazards and 28 with IID logit random effects on the hazards. pcd_as_stan_data() accepts dist_options = list(K = ..., boundaries = ...) to carry the bin shape, and the existing priors argument carries the prior: priors$scale is the length-K Dirichlet concentration for dist_id = 26, and priors$location / priors$scale are length-2 (mean, sd) for alpha and log_sigma for dist_id = 27 and 28.pdist argument of pprimarycensored(), dprimarycensored(), and qprimarycensored() now accepts a character string, looked up via the pcd_distributions registry; passing a function still works as before.pprimary is now a user-facing argument on pprimarycensored(), dprimarycensored(), qprimarycensored(), and fitdistdoublecens(). It accepts either a function or a character string naming a primary distribution in pcd_primary_distributions. When omitted, the primary CDF is looked up automatically from the registry using the "name" attribute of dprimary. Supplying both dprimary and pprimary with inconsistent names (e.g. dunif + pexpgrowth) now errors clearly.dprimary_args in favour of primary_args on pprimarycensored(), dprimarycensored(), qprimarycensored(), new_pcens(), and fitdistdoublecens(). The new name reflects that the arguments are passed to both dprimary and the matching primary CDF. The deprecation uses lifecycle::deprecate_soft(), so a warning is shown for direct use and during package tests, and calls from other packages keep working quietly.pcens objects: a specific method for a (delay, primary) pair is tried first, then a delay-only general method, then a numerical default. The analytic primary convolution path now works for any primary with a known CDF, with the primary CDF plumbed through pprimarycensored(), dprimarycensored(), and qprimarycensored() via the pprimary argument.dprimarycensored() and fitdistdoublecens() now accept observations whose secondary censoring interval straddles D (left < D <= right). The upper endpoint is internally clipped to D and the likelihood becomes P(X in [left, min(right, D)] | L <= X <= D). This is a no-op when right <= D and removes the need to pad D when fitting non-parametric delays whose support reaches D. Observations with left >= D are still rejected because under truncation at D no event with latent value >= D is observable. See #312.fitting-nonparametric-delays demonstrating end-to-end non-parametric delay estimation.pprimarycensored(), dprimarycensored() and qprimarycensored() gain a check argument. It defaults to TRUE, which keeps the existing validation of pdist via check_pdist() and dprimary via check_dprimary(). Setting check = FALSE skips both, for callers that have already validated their functions. Because check_pdist() evaluates pdist at four points drawn with runif(), skipping it also leaves the random number stream untouched, so seeded code no longer depends on how many times these functions were called. check follows ... so it must be given by its full name and cannot capture an argument intended for pdist. See #330.fitdistdoublecens() gains a matching check argument.NaN gradient in the Stan primarycensored_lcdf() deep in the lower tail of a narrow lognormal delay. lognormal_lcdf() underflows to -inf once the standardised value falls below about -38.6, and its autodiff partial is then 0 / 0. Stan's reverse pass chains that NaN into mu and sigma even where the term carries zero weight, so the log density came back finite while the gradient did not. Models that evaluate the delay on a grid starting at zero reached this as soon as a proposal was narrow, and saw only Gradient evaluated at the initial value is not finite. Affected terms are now dropped before the underflowing call, in both the analytic uniform primary solution and the ODE path. Where the density really is zero the result is reported as log(0). See #333.dprimarycensored() validated four times, once directly and once inside each of its three internal pprimarycensored() calls. Inside fitdistdoublecens() validation ran once per observation per likelihood evaluation, so a 100 row fit validated several hundred times. It now runs once per fit.This patch release fixes a performance regression introduced in 1.5.0 that slowed the Stan likelihood for positive-support delays.
primarycensored_lpmf and primarycensored_lcdf functions for positive-support delays.
After the lower-truncation guard was relaxed in 1.5.0 to support negative L, every likelihood evaluation for positive-support delays entered a truncation-normalisation block that cancels to a no-op but still adds gradient calculations (an exp/log_diff_exp per call), giving a roughly 30% slowdown on the likelihood block.
The guard now skips this block unless there is a finite D, a strictly positive L, or a finite L on a real-support distribution, restoring 1.4.0 performance without changing results or losing the negative-support behaviour.
Thanks to @sbfnk for reporting (#323).This minor release extends the L (lower truncation) parameter to accept negative and -Inf values in both the R and Stan code, letting delay distributions with support below zero (e.g. normal, logistic, Cauchy, Gumbel) be used with primary censoring and fitted via fitdistdoublecens() and pcd_cmdstan_model(). The default value of L has changed from 0 to -Inf, which leaves results unchanged for distributions with non-negative support (e.g. lognormal, gamma, Weibull). It also fixes a normalisation bug in the exponential growth primary distribution and rewrites the analytical CDFs in a CDF-direct form, dropping the pracma dependency.
L in pprimarycensored(), dprimarycensored(), qprimarycensored(), rprimarycensored(), and pcens_quantile() has changed from 0 to -Inf. fitdistdoublecens() now treats a missing L column as L = -Inf to match. For delay distributions with support on the non-negative reals (e.g. lognormal, gamma, Weibull) F_cens(0) = 0, so the new default leaves results unchanged. Callers that relied on the implicit left truncation at 0 (for example to truncate signed-support delays such as pnorm) must now pass L = 0 explicitly (or add an L column in fitdistdoublecens()). (#267)L may now be negative or -Inf in pprimarycensored(), dprimarycensored(), qprimarycensored(), and rprimarycensored(). This lets delay distributions with support below zero (e.g. normal, logistic, Cauchy) be used with primary censoring. L = -Inf is the sentinel for "no left truncation"; any finite L left-truncates the distribution at L. (#267)pcd_as_stan_data() mirror the R-side handling of L: negative and -Inf values are accepted, and a missing start_relative_obs_time column defaults to -Inf. pcd_cmdstan_model() now accepts negative observed delays and fully-negative truncation windows, letting distributions with support on the reals (e.g. logistic, Cauchy, Gumbel) be fitted. (#313)dist_id upper bound in pcens_model.stan has been raised from 17 to 25, exposing every delay distribution that dist_lcdf already dispatches (Normal, Double Exponential, Pareto, scaled inverse chi-square, Student's t, Uniform, von Mises) through pcd_cmdstan_model(). (#314)fitdistdoublecens() and pcd_cmdstan_model() from doubly-censored, right-truncated samples that include negative observed delays.dexpgrowth(), pexpgrowth(), and their Stan equivalents (expgrowth_pdf, expgrowth_lpdf, expgrowth_cdf, expgrowth_rng) when min is non-zero.
The PDF and CDF formulas contained a stray exp(-r * min) factor from using exp(r * (x - min)) instead of exp(r * x).
The Stan RNG had a compensating xmin + offset.
The bug did not affect results when min = 0 (the default and only value used within the package's primary censoring functions).
Thanks to @TimTaylor for reporting (#290).q = 0 and q > 0 code paths (single algebraic expression, better for NUTS), the outer log_diff_exp ordering is now mathematically guaranteed, and the Gamma case uses the incomplete-gamma recursion $P(k{+}1, y) = P(k, y) - y^k e^{-y}/\Gamma(k{+}1)$ to halve gamma_lcdf evaluations. Behaviour and tests are unchanged.pracma::gammainc with stats::pgamma in the Weibull g() helper, dropping the pracma dependency. The previous pwindow > 3 fallback to numeric integration (and the internal overflow guard) is no longer needed — the base R implementation is stable across the full parameter range. Closes #127.This major release adds left-truncation support via the L parameter, enabling distributions to be truncated over [L, D] rather than just [0, D]. It also removes deprecated functionality that was soft-deprecated in version 1.1.0 and removes the lifecycle and rlang packages from dependencies.
pdist_name and dprimary_name arguments from pprimarycensored(), dprimarycensored(), new_pcens(), and fitdistdoublecens(). Use add_name_attribute() on the pdist and dprimary functions instead to enable analytical solutions.pwindow and D arguments in fitdistdoublecens(). These must now be column names in censdata.lifecycle and rlang packages from dependencies.L parameter before D in their signatures. If you were using positional arguments for D (e.g., pprimarycensored(x, pdist, pwindow, 10)), you must now use named arguments (e.g., pprimarycensored(x, pdist, pwindow, D = 10)). The affected functions are dprimarycensored(), pprimarycensored(), rprimarycensored(), and qprimarycensored(). (#63)L parameter for lower truncation, placed before D. This affects primarycensored_lpmf, primarycensored_lcdf, primarycensored_cdf, primarycensored_pmf, and all vectorized variants. Update your Stan code to include the new parameter: primarycensored_lpmf(d | dist_id, params, pwindow, d_upper, L, D, primary_id, primary_params). (#63)dependencies argument to pcd_load_stan_functions() that automatically resolves and includes all functions that the requested functions depend on. When TRUE, dependencies are included in topological order (dependencies before functions that use them). (#171)pcd_stan_function_deps() to query the dependency graph of Stan functions, returning all dependencies for a given function in topological order. (#171)pcd_distributions that were previously only available in Stan: normal, inverse chi-square, double exponential, pareto, scaled inverse chi-square, student t, uniform, and von Mises (IDs 18–25). These are now accessible via pcd_stan_dist_id(). (#277)L parameter to all primary censored distribution functions (dprimarycensored(), pprimarycensored(), rprimarycensored(), qprimarycensored()) and Stan functions. The L parameter specifies the minimum delay (lower truncation point), enabling distributions to be truncated over [L, D] rather than just [0, D]. This is useful for generation intervals and other settings where delays below a threshold cannot occur. Defaults to L = 0 for backward compatibility. (#63)expgrowth_lpdf Stan function to support negative growth rates. Previously, log(r) returned NaN for r < 0. Now uses log(abs(r)) and log(abs(...)) for the denominator. (#276)dist_lcdf distribution ID mapping to match R's pcd_distributions table. Previously, most distribution IDs were mismatched between R and Stan (e.g. pcd_stan_dist_id("weibull") returned 3 but Stan's dist_lcdf used ID 3 for the Normal distribution). Only lognormal (1), gamma (2), and exponential (4) were correct. This affected the ODE numerical integration path for all mismatched distributions. (#277)min and max parameters to xmin and xmax in Stan functions (expgrowth_pdf, expgrowth_lpdf, expgrowth_cdf, expgrowth_lcdf, expgrowth_rng, primary_lpdf) to avoid conflicts with Stan built-in functions. CmdStan 2.38.0 now strictly enforces reserved keyword restrictions when exposing Stan functions to R. (#258)This minor release improves documentation for fitdistdoublecens() and adds learning objective sections to vignettes. It also fixes floating-point precision issues in dprimarycensored() and adds bounds checking to CDF methods.
fitdistdoublecens():
@details sections explaining how distribution names are resolved (e.g., "gamma" uses dgamma() and pgamma()), and what the function does internally.distr parameter documentation with examples and guidance on custom distributions.pcd_as_stan_data() and fitdistdoublecens() correctly handle data frames with additional columns beyond those required. (#213)dprimarycensored() could return very small negative values (e.g., -2.2e-16) due to floating-point precision when computing PMF as CDF differences. PMF values are now clamped to be non-negative. (#238)pcens_cdf() methods to ensure CDF values are always in [0, 1], complementing the existing upper bound check with a lower bound check.This minor release adds quantile functions for primary censored distributions and enhances the fitdistdoublecens() function to support varying primary censoring windows and truncation times. The release also improves documentation, particularly the Stan reference, making it easier for users to navigate and work with the Stan code.
pcens_quantile(): Core implementation using numerical optimisation to find
quantiles by inverting the CDF. The implementation allows for analytical
solutions to be added in future versions using the same interface as pcens_cdf().qpcens() and qprimarycensored(): Convenient wrapper functions that provide
alternative interfaces to pcens_quantile().fitdistdoublecens() to allow for varying primary censoring windows and truncation times. As part of this refactor the interface has also been improved to allow for more flexible data input.fitdistrplus vignette to use more complex data where the relative observation time is not constant. Also removed the note that the fitdistdoublecens() function is limited to a single primary censoring windows and truncation time as this is no longer the case.This minor release improves the documentation and the internal system used to automatically discover analytical solutions. It also adds lookups for supported supported distributions and tooling for working with these lookups. This makes it easier for package developers using primarycensored and also makes it easier to work with the Stan likelihood functions by improving the discoverability of the distribution functions. Minor bugs were also fixed.
epidist.pdist_name and dprimary_name arguments throughout. Users wishing to pass distribution names (i.e. to potentially leverage analytical solutions) are advised to use the newly introduced add_name_attribute() function. Adds transient dependency on lifecycle and rlang packages. See #188 by @pearsonca.pcd_stan_dist_id() to allow for discovery of distribution IDs for Stan models.pcd_dist_name() to allow for discovery of distribution names for R functions as needed for add_name_attribute().methods(pcens_cdf) to find analytical solutions.pcd_distributions and pcd_primary_distributions to document the distributions supported by primarycensored.primarycensored_lcdf().@family tag to the pcens functions. This omission resulted in the Weibull analytical solution not being visible in the package documentation.primarycensored_cdf() stan function, avoiding errors on some platforms due to narrowing conversions in aggregate initialisation.D to be of type real in pcens_model.stan in order to support infinite relative_obs_time.num_elements() over size() in all stan code to resolve compilation issues on some platforms as reported by @sbfnk.This is the first major release of primarycensored and has been submitted to CRAN.
primary_lpdf.rhub checks to the Github Actions workflow.dependencies: "hard" to the R-CMD-check workflow to ensure checks pass without optional dependencies.This release renames the package to primarycensored from primarycensoredist and also renames many of the functions to remove the dist in their name. This was done to make the package name and the functions more consistent and to remove the need to use the dist suffix. It also aligns it with the new PrimaryCensored.jl package in our Julia ecosystem.
Aside from name changes, this release also adds an analytical solution for the weibull distribution with uniform primary censoring, removes the need to assign functions to the global environment for fitdistdoublecens() by using withr, and adds a check_truncation() function to check if the truncation time is larger than the maximum observed delay. This is used in fitdistdoublecens() and pcd_as_stan_data() to ensure that the truncation time is appropriate to maximise computational efficiency.
fitdistdoublecens() by using withr.check_truncation() function to check if the truncation time is larger than the maximum observed delay. This is used in fitdistdoublecens() and pcd_as_stan_data() to ensure that the truncation time is appropriate to maximise computational efficiency.pcd_as_cmdstan_data() has been renamed to pcd_as_stan_data() to better reflect that it is used for Stan models in general rather than just the CmdStan models.stan folder and the stan model has been moved into the stan folder. All paths to the stan code have been updated to reflect this.primarycensored as have all functions that use "dist" in their name.pdist, dprimary, rdist, and rprimary arguments in the getting started vignette to make it easier to link to mathematical details.This release adds a new {touchstone} based benchmark suite to the package. It also adds a new "How it works" vignette which aims to give the reader more details into how the primary censored distributions work.
As part of the "How it works" we (@SamuelBrand1) found analytical solutions for the gamma, lognormal, and weibull distributions with uniform primary censoring. These are now implemented for the lognormal and gamma distributions in the R and stan code providing significant speedups to the fitting process (~10-20 times faster). The Weibull will be added in the next release.
{touchstone} based benchmarks for benchmarking R utility functions, and fitting the stan and fitdistplus models.primarycensored S3 class.R and stan code.In this release, we have added a new package stan model for fitting distributions using the cmdstanr package. We have also added a new function fitdistdoublecens() to allow for fitting of double censored and truncated data using the fitdistrplus package. As well as these functionality improvements this release focuses on improving the stability of the stan model and improving the speed of the primarycensored_ode function.
fitdistdoublecens() to allow for fitting of double censored and truncated data using the fitdistrplus package.primarycensored_ode function.CmdStan model for fitting distributions using the cmdstanr package.CmdStan model and added an example to the vignette.CmdStan model which tests the primarycensored_lpmf function when used with NUTS based fitting.This release fixes and improves truncation handling across the code base. It also adds a new vignette showcasing how to use the primarycensored and fitdistrplus packages together to fit distributions.
primarycensored_lpmf when used for NUTS based fitting (i.e. in Stan).primarycensored and fitdistrplus packages together to fit distributions.This release puts in place initial documentation and vignettes. It also includes a new primary censored distribution interface to allow for non-secondary event censored distributions. Development of this release as identified some numerical issues in the
gradient evaluations for the primary censored distributions which may lead to breaking
interface changes in 0.3.0 for the Stan code.
swindow = 0 to rprimarycensored to allow for non-secondary event censored distributions.rprimarycensored so that truncation is based on the primary censored distribution before secondary events are censored. This better matches the generative process.cmdstanr.cmdstanr package.This is the initial primarycensored release and includes R and stan tools for dealing with potentially truncated primary event censored delay distributions. We expect all current features to work but the UI may change as the package matures over the next few versions.
0.2.0 release.