Fix #534: Add marginal model to model outline document (#594)

Co-authored-by: Sam Abbott <contact@samabbott.co.uk>

Author: seabbs-bot

README.Rmd

@@ -142,9 +142,11 @@ citation("epidist")
 If using our methodology, or the methodology on which ours is based, please cite the relevant papers.
 This may include:
 
-* [Estimating epidemiological delay distributions for infectious diseases](https://www.medrxiv.org/content/10.1101/2024.01.12.24301247v1) by Park *et al.* (2024)
+* [Estimating epidemiological delay distributions for infectious diseases](https://www.medrxiv.org/content/10.1101/2024.01.12.24301247v1) by Park *et al.* (2024) -- if using the latent model
 * [Best practices for estimating and reporting epidemiological delay distributions of infectious diseases using public health surveillance and healthcare data](https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1012520) by Charniga *et al.* (2024)
 
+If using the marginal model, please also cite the [`primarycensored`](https://primarycensored.epinowcast.org/) package using `citation("primarycensored")`.
+
 ## Contributors
 
 

vignettes/model.Rmd

@@ -154,6 +154,7 @@ This then allows the modelling of the continuous distribution, adjusted for the
 Whilst this is an approximation [@park2024estimating] showed good recovery of simulated distributions in a range of settings.
 However, the use of two latent variables per observed delay means that this approach may scale poorly with larger datasets.
 That being said this approach has been used successfully in multiple real-world outbreak settings ([@ward2022transmission]).
+If using the latent model, please cite @park2024estimating in addition to `epidist`.
 
 Mathematically this model is described as follows.
 We look at the conditional probability that the secondary event $S$ falls between $S_L$ and $S_R$, given that the primary event $P$ falls between $P_L$ and $P_R$ and that the secondary event $S$ occurs before the truncation time $T$:
@@ -177,4 +178,28 @@ y_i &\sim \text{Unif}(s_{L, i}, s_{R, i}) \\
 $$
 As before, $g_P(z \, | \, p_{L, i}, p_{R, i})$ represents the conditional distribution of the primary event given lower $P_L$ and upper $P_R$ bounds; this is equivalent to modelling the incidence in primary events.
 
+# The marginal model
+
+The marginal model corrects for the same biases as the latent model but integrates out the exact event times numerically, or analytically where closed-form solutions exist, rather than sampling latent variables.
+This approach uses the primary event censored distribution implemented in the [`primarycensored`](https://primarycensored.epinowcast.org/) package [@primarycensored].
+If using the marginal model, please cite `primarycensored` in addition to `epidist`.
+
+Under the assumption that the forward distribution does not change within the censoring interval (i.e. $f_x = f$ for $x \in [P_L, P_R]$), the double censoring probability from Section \@ref(interval-censoring) simplifies to
+$$
+\mathbb{P}(S_L < S < S_R \mid P_L < P < P_R) = \int_{P_L}^{P_R} g_P(x \mid P_L, P_R) \left[F(S_R - x) - F(S_L - x)\right] \text{d}x.
+$$
+For common delay and primary event distributions, such as gamma or lognormal delays with uniform primary events, `primarycensored` provides closed-form analytical solutions to this integral.
+For other combinations, numerical integration is used.
+
+Right truncation at time $T$ is handled by normalising the likelihood as in the latent model:
+$$
+\mathcal{L}(\mathbf{Y} \mid \mathbf{\theta}) = \prod_i \frac{\mathbb{P}(S_{L,i} < S_i < S_{R,i} \mid P_{L,i} < P_i < P_{R,i})}{\int_{P_{L,i}}^{P_{R,i}} g_P(z \mid p_{L,i}, p_{R,i}) F(T - z) \, \text{d}z}.
+$$
+
+Removing the latent variables reduces the number of parameters that must be sampled, and where analytical solutions exist the likelihood can be evaluated without numerical integration.
+In addition, identical observations can be aggregated and the likelihood computed once per unique combination of delay, censoring windows, and covariates.
+Together these make the marginal model substantially more efficient than the latent model, particularly for larger datasets with daily-censored data where many observations share the same structure.
+
+For the mathematical details of primary event censored distributions, including the survival function derivation and closed-form solutions for specific distributions, see `vignette("why-it-works", package = "primarycensored")` and `vignette("analytic-solutions", package = "primarycensored")`.
+
 ## References {-}

vignettes/references.bib

@@ -10,6 +10,13 @@ @article {park2024estimating
 	journal = {medRxiv}
 }
 
+@Manual{primarycensored,
+  title = {primarycensored: Primary Event Censored Distributions},
+  author = {Sam Abbott and Sam Brand and James Mba Azam and Carl Pearson and Sebastian Funk and Kelly Charniga},
+  year = {2025},
+  doi = {10.5281/zenodo.13632839},
+}
+
 @article{charniga2024best,
     doi = {10.1371/journal.pcbi.1012520},
     author = {Charniga, Kelly and Park, Sang Woo and Akhmetzhanov, Andrei R. and Cori, Anne and Dushoff, Jonathan and Funk, Sebastian and Gostic, Katelyn M. and Linton, Natalie M. and Lison, Adrian and Overton, Christopher E. and Pulliam, Juliet R. C. and Ward, Thomas and Cauchemez, Simon and Abbott, Sam},