diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 8c52264e..854bca90 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -14,7 +14,7 @@ repos: ##### # Python - repo: https://github.com/astral-sh/ruff-pre-commit - rev: v0.16.1 + rev: v0.16.5 hooks: # Run the linter - id: ruff-check @@ -27,7 +27,7 @@ repos: ##### # Quarto - repo: https://github.com/jolars/panache-pre-commit - rev: v3.0.2 + rev: v3.7.0 hooks: - id: panache-format - id: panache-lint diff --git a/README.md b/README.md index cf0d576a..c692aa96 100755 --- a/README.md +++ b/README.md @@ -8,7 +8,7 @@ From this, it infers $\mathcal{R}(t)$, the time-varying reproduction number, whi The core renewal equation is: $$ - I(t) = \mathcal{R}(t) \sum_{s} I(t-s) \, w(s) + I(t) = \mathcal{R}(t) \sum_{s} I(t - s) \, w(s) $$ where $w(s)$ is the generation interval distribution: the probability that $s$ time units separate infection in an index case and a secondary case. @@ -17,7 +17,7 @@ However inference is complicated by the fact that observational data require the The observation equation links infections to expected observations: $$ - \mu(t) = \alpha \sum_{s} I(t-s) \, \pi(s) + \mu(t) = \alpha \sum_{s} I(t - s) \, \pi(s) $$ where $\alpha$ is the ascertainment rate and $\pi(s)$ is the delay distribution from infection to observation. diff --git a/docs/tutorials/ascertainment.qmd b/docs/tutorials/ascertainment.qmd index bba71461..8fcd3ce1 100644 --- a/docs/tutorials/ascertainment.qmd +++ b/docs/tutorials/ascertainment.qmd @@ -41,9 +41,7 @@ For emergency department visits, it can be called the infection-to-ED-visit rate For a count observation process, the ascertainment rate scales delayed latent infections to obtain the expected number of observed events: $$ - \mu_s(t) - = - \alpha_s \sum_{d=0}^{D_s} I(t-d)\,\pi_{s,d}, + \mu_s(t) = \alpha_s \sum_{d=0}^{D_s} I(t - d)\,\pi_{s,d}, $$ where $s$ identifies the signal, $\alpha_s$ is its ascertainment rate, and $\pi_{s,d}$ is its delay distribution. @@ -132,14 +130,8 @@ Use `RatioLinkedAscertainment` when the relative quantity, such as IHR / IEDR, i $$ \boldsymbol{\eta} - \sim - \operatorname{MultivariateNormal} - \left( - \operatorname{logit}(\boldsymbol{\alpha}_0), - \boldsymbol{\Sigma} - \right), - \qquad - \boldsymbol{\alpha} = \operatorname{logit}^{-1}(\boldsymbol{\eta}). + \sim \operatorname{MultivariateNormal} \left( \operatorname{logit}(\boldsymbol{\alpha}_0), \boldsymbol{\Sigma} \right), \qquad \boldsymbol{\alpha} + = \operatorname{logit}^{-1}(\boldsymbol{\eta}). $$ The baseline rates $\boldsymbol{\alpha}_0$ set the prior location. @@ -170,13 +162,8 @@ During model execution, this component creates one vector sample site, `he_ascer `RatioLinkedAscertainment` samples a base rate and a ratio of the linked rate to the base rate: $$ - \alpha_{\mathrm{linked}} - = - \alpha_{\mathrm{base}}\rho, - \qquad - \rho - = - \frac{\alpha_{\mathrm{linked}}}{\alpha_{\mathrm{base}}}. + \alpha_{\mathrm{linked}} = \alpha_{\mathrm{base}}\rho, \qquad \rho + = \frac{\alpha_{\mathrm{linked}}}{\alpha_{\mathrm{base}}}. $$ For a hospital and ED model, the IEDR can be the base rate and the ratio can represent IHR / IEDR. @@ -203,7 +190,7 @@ The priors are assigned to the base rate and the ratio, and together they induce The ratio may be any positive value. Unlike `JointAscertainment`, this component does not automatically constrain the derived linked rate to the unit interval. -If both rates should be $<1$ (for instance because both are probabilities), make sure to choose priors that make $\alpha_{\mathrm{base}}\rho > 1$ negligibly likely. +If both rates should be $< 1$ (for instance because both are probabilities), make sure to choose priors that make $\alpha_{\mathrm{base}}\rho > 1$ negligibly likely. The component does not clip an invalid linked rate, because clipping would change the specified statistical model. ## Building a model an ascertainment component diff --git a/docs/tutorials/day_of_week_effects.qmd b/docs/tutorials/day_of_week_effects.qmd index 357f4afe..7322ec42 100644 --- a/docs/tutorials/day_of_week_effects.qmd +++ b/docs/tutorials/day_of_week_effects.qmd @@ -50,7 +50,7 @@ Ignoring this weekly periodicity forces the noise model to absorb systematic var PyRenew models day-of-week effects as a **multiplicative adjustment** applied to predicted counts after the delay convolution and ascertainment scaling: $$ - \lambda(t) = d_{w(t)} \cdot \alpha \sum_{s} I(t-s)\,\pi(s) + \lambda(t) = d_{w(t)} \cdot \alpha \sum_{s} I(t - s)\,\pi(s) $$ where $d_{w(t)}$ is the day-of-week multiplier for the weekday of timepoint $t$, $\alpha$ is the ascertainment rate, and $\pi(s)$ is the delay PMF. @@ -335,7 +335,7 @@ Day-of-week effects and right-truncation are independent adjustments that compos Day-of-week is applied first (adjusting the expected counts for reporting patterns), then right-truncation scales down recent counts for incomplete reporting: $$ - \lambda(t) = F(k_t) \cdot d_{w(t)} \cdot \alpha \sum_s I(t-s)\,\pi(s) + \lambda(t) = F(k_t) \cdot d_{w(t)} \cdot \alpha \sum_s I(t - s)\,\pi(s) $$ ```{python} diff --git a/docs/tutorials/latent_infections.qmd b/docs/tutorials/latent_infections.qmd index 8a774f60..c11d6918 100644 --- a/docs/tutorials/latent_infections.qmd +++ b/docs/tutorials/latent_infections.qmd @@ -156,17 +156,17 @@ Let $\tau = 0, 1, \ldots, n_{\text{init}} - 1$ index positions in the initializa Then: $$ - I_{\text{init}}(\tau) = I_0 \cdot e^{r \cdot \tau}, \quad \tau = 0, 1, - \ldots, n_{\text{init}} - 1 + I_{\text{init}}(\tau) = I_0 \cdot e^{r \cdot \tau}, \quad \tau + = 0, 1, \ldots, n_{\text{init}} - 1 $$ -where $r$ is the asymptotic growth rate implied by the reproduction number at the start of the observation period, $\mathcal{R}(t=0) = e^{\text{log\_rt\_time\_0}}$, and the generation interval. -The function `r_approx_from_R` converts $\mathcal{R}(t=0)$ and the generation interval into $r$ using Newton's method. +where $r$ is the asymptotic growth rate implied by the reproduction number at the start of the observation period, $\mathcal{R}(t = 0) = e^{\text{log\_rt\_time\_0}}$, and the generation interval. +The function `r_approx_from_R` converts $\mathcal{R}(t = 0)$ and the generation interval into $r$ using Newton's method. - **The level is set by `I0`**.
`I0` is the infection prevalence at the earliest point in the initialization period, $n_{\text{init}} - 1$ time points before $t = 0$. It sets the scale of the entire initialization vector: $I_{\text{init}}(0) = I_0$, with subsequent entries growing or declining exponentially toward $t = 0$. -- **The shape is set by `log_rt_time_0`**.
`log_rt_time_0` enters the model in two places: it is the starting point of the $\mathcal{R}(t)$ trajectory ($\mathcal{R}(t=0) = e^{\text{log\_rt\_time\_0}}$), and it determines the exponential growth rate $r$ used to construct the initialization vector. +- **The shape is set by `log_rt_time_0`**.
`log_rt_time_0` enters the model in two places: it is the starting point of the $\mathcal{R}(t)$ trajectory ($\mathcal{R}(t = 0) = e^{\text{log\_rt\_time\_0}}$), and it determines the exponential growth rate $r$ used to construct the initialization vector. When `log_rt_time_0 = 0`, $r = 0$ and the initialization vector is flat at level `I0`. When `log_rt_time_0 > 0`, infections are growing exponentially at $t = 0$; when `log_rt_time_0 < 0`, they are declining. @@ -232,7 +232,7 @@ init_df = pd.DataFrame(init_data) ) ``` -The initialization vector matters because the renewal equation is a convolution: infections on day 0 depend on infections from days $-1$ through $-(K-1)$, weighted by the generation interval. +The initialization vector matters because the renewal equation is a convolution: infections on day 0 depend on infections from days $-1$ through $-(K - 1)$, weighted by the generation interval. A flat initialization (stable) means the renewal equation starts with uniform recent history. A growing initialization means the most recent days have disproportionately more infections, which amplifies the effect of the generation interval's short-lag weights. @@ -362,7 +362,7 @@ The variance of $x_t$ grows linearly with time: $\text{Var}(x_t) = \sigma^2 t$. The further into the future, the less constrained the process is. **Hyperparameter:** `innovation_sd` ($\sigma$) is the standard deviation of each daily step on the log scale. -With `innovation_sd = 0.05`, each day's $\log \mathcal{R}$ changes by roughly $\pm 0.05$, which corresponds to roughly $\pm 5\%$ multiplicative change in $\mathcal{R}$. +With `innovation_sd = 0.05`, each day's $\log \mathcal{R}$ changes by roughly $\pm0.05$, which corresponds to roughly $\pm5\%$ multiplicative change in $\mathcal{R}$. ```{python} #| label: rw-sample @@ -408,7 +408,7 @@ If $\mathcal{R}(t)$ drifts above 1, the $\phi$ coefficient pulls it back; if it - `innovation_sd` ($\sigma$): standard deviation of daily noise. The two hyperparameters jointly determine the **stationary standard deviation** $\sigma_{\text{stat}} = \sigma / \sqrt{1 - \phi^2}$, which is the long-run spread of $\log \mathcal{R}(t)$. -For example, `autoreg = 0.9` and `innovation_sd = 0.05` give $\sigma_{\text{stat}} \approx 0.115$, meaning 95% of long-run $\log \mathcal{R}$ values fall within $\pm 0.23$ of zero, or equivalently $\mathcal{R} \in [0.79, 1.26]$. +For example, `autoreg = 0.9` and `innovation_sd = 0.05` give $\sigma_{\text{stat}} \approx 0.115$, meaning 95% of long-run $\log \mathcal{R}$ values fall within $\pm0.23$ of zero, or equivalently $\mathcal{R} \in [0.79, 1.26]$. ```{python} #| label: ar1-sample @@ -441,7 +441,8 @@ DifferencedAR1 models *changes* in $\mathcal{R}(t)$ as an autoregressive process In terms of increments, $$ - \Delta x_t = \texttt{autoreg} \cdot \Delta x_{t-1} + \varepsilon_t, \quad x_t = x_{t-1} + \Delta x_t, + \Delta x_t = \texttt{autoreg} \cdot \Delta x_{t-1} + \varepsilon_t, \quad x_t + = x_{t-1} + \Delta x_t, $$ so the autoregressive structure applies to the *rate of change* rather than the level. diff --git a/docs/tutorials/latent_subpopulation_infections.qmd b/docs/tutorials/latent_subpopulation_infections.qmd index 39c2e099..23174e4d 100644 --- a/docs/tutorials/latent_subpopulation_infections.qmd +++ b/docs/tutorials/latent_subpopulation_infections.qmd @@ -146,10 +146,9 @@ Two temporal processes define the evolution of $\log \mathcal{R}_k(t)$: Together, these define the full set of reproduction numbers: $$ - \log \mathcal{R}_k(t) - = \log \mathcal{R}_{\text{baseline}}(t) + \delta_k(t), - \qquad - \mathcal{R}_k(t) = \exp\big(\log \mathcal{R}_k(t)\big). + \log \mathcal{R}_k(t) = \log \mathcal{R}_{\text{baseline}}(t) + + \delta_k(t), \qquad \mathcal{R}_k(t) + = \exp\big(\log \mathcal{R}_k(t)\big). $$ ### Population structure diff --git a/docs/tutorials/observation_processes_counts.qmd b/docs/tutorials/observation_processes_counts.qmd index 392e79fc..b0aaefff 100644 --- a/docs/tutorials/observation_processes_counts.qmd +++ b/docs/tutorials/observation_processes_counts.qmd @@ -59,12 +59,12 @@ Observed data can be aggregated or available as subpopulation-level counts, whic The deterministic transformation is given by the observation equation: $$ - \mu(t) = \alpha \sum_{d=0}^{D} I(t-d)\, \pi_d + \mu(t) = \alpha \sum_{d=0}^{D} I(t - d)\, \pi_d $$ where: -- $I(t-d)$ is the number of incident (new) infections on day $t-d$ +- $I(t - d)$ is the number of incident (new) infections on day $t - d$ - $\alpha$ is the **ascertainment rate**, the probability that an infection results in an observed event (e.g., hospitalization) - $\pi_d$ is the delay distribution from infection to observation, conditional on an infection leading to an observed event @@ -103,7 +103,7 @@ In this tutorial, we simulate the observation process by specifying infections d For hospital admissions data, we construct a `PopulationCounts` observation process. -The delay is the key mechanism: infections from $d$ days ago ($I(t-d)$) contribute to today's predicted hospital admissions ($\mu(t)$), weighted by the probability $\pi_d$ that an infection leads to hospitalization after exactly $d$ days. +The delay is the key mechanism: infections from $d$ days ago ($I(t - d)$) contribute to today's predicted hospital admissions ($\mu(t)$), weighted by the probability $\pi_d$ that an infection leads to hospitalization after exactly $d$ days. The convolution sums these contributions across all past days. Observed hospital admissions are then generated by sampling from a negative binomial distribution: @@ -223,8 +223,8 @@ hosp_process = PopulationCounts( The observation process convolves infections with a delay distribution, maintaining alignment between input and output: day $t$ in the output corresponds to day $t$ in the input. Hospital admissions depend on infections from prior days. -A delay PMF of length $L$ covers delays 0 to $L-1$, requiring $L-1$ days of prior infection history. -The method `lookback_days()` returns $L-1$; the first valid observation day is at index `lookback_days()`. +A delay PMF of length $L$ covers delays 0 to $L - 1$, requiring $L - 1$ days of prior infection history. +The method `lookback_days()` returns $L - 1$; the first valid observation day is at index `lookback_days()`. Earlier days are marked invalid. ```{python} diff --git a/docs/tutorials/observation_processes_measurements.qmd b/docs/tutorials/observation_processes_measurements.qmd index a8316b8c..3c013ff7 100644 --- a/docs/tutorials/observation_processes_measurements.qmd +++ b/docs/tutorials/observation_processes_measurements.qmd @@ -198,12 +198,12 @@ The concentration is determined by: The predicted concentration at time $t$ is given by $$ - \mu(t) = \frac{G}{V} \sum_{d=0}^{D} I(t-d)\, \pi_d, + \mu(t) = \frac{G}{V} \sum_{d=0}^{D} I(t - d)\, \pi_d, $$ where: -- $I(t-d)$ is the number of infections on day $t-d$ +- $I(t - d)$ is the number of infections on day $t - d$ - $\pi_d$ is the shedding kinetics PMF, giving the fraction of total shedding occurring $d$ days after infection, analogous to the delay distribution in count observation models - $G$ is the number of genome copies shed per infection - $V$ is the wastewater volume per person per day @@ -213,7 +213,8 @@ This has the same convolution structure as the observation equation for count da We model the observed **log-concentration** as $$ - Y(t) \sim \text{Normal}\big(\log(\mu(t)) + \text{sensor\_mode}, \ \text{sensor\_sd}\big), + Y(t) + \sim \text{Normal}\big(\log(\mu(t)) + \text{sensor\_mode}, \ \text{sensor\_sd}\big), $$ where: @@ -437,8 +438,8 @@ print(f"Required lookback: {ww_process.lookback_days()} days") ### Timeline alignment The observation process maintains alignment: day $t$ in output corresponds to day $t$ in input. -A temporal PMF of length $L$ covers lags 0 to $L-1$, requiring $L-1$ days of prior history. -The method `lookback_days()` returns $L-1$; the first valid observation day is at index `lookback_days()`. +A temporal PMF of length $L$ covers lags 0 to $L - 1$, requiring $L - 1$ days of prior history. +The method `lookback_days()` returns $L - 1$; the first valid observation day is at index `lookback_days()`. ```{python} #| label: helper-function diff --git a/docs/tutorials/right_truncation.qmd b/docs/tutorials/right_truncation.qmd index c3809a4e..95429db8 100644 --- a/docs/tutorials/right_truncation.qmd +++ b/docs/tutorials/right_truncation.qmd @@ -47,7 +47,7 @@ Ignoring this produces a spurious decline in recent counts. PyRenew's observation equation defines the expected observation count as: $$ - \mu(t) = \alpha \sum_{s} I(t-s) \, \pi(s) + \mu(t) = \alpha \sum_{s} I(t - s) \, \pi(s) $$ where $\alpha$ is the ascertainment rate and $\pi(s)$ is the infection-to-observation delay distribution. @@ -65,7 +65,8 @@ Concretely, let $T$ denote the last observation day and let $\text{offset} = \te Then: $$ - k_t = (T - t) + \text{offset} = (T - t) + (\text{data pull date} - T) = \text{data pull date} - t + k_t = (T - t) + \text{offset} = (T - t) + (\text{data pull date} - T) + = \text{data pull date} - t $$ Because $k_t$ depends only on the data pull date and the timepoint $t$, its behavior is straightforward: timepoints far in the past (small $t$) have large $k_t$, so $F(k_t) \approx 1$ and counts are fully reported.