Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Renaming variables throughout for clarity #179

Merged
merged 10 commits into from
Jun 12, 2024
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 27 additions & 25 deletions model/docs/example_with_datasets.qmd
Original file line number Diff line number Diff line change
@@ -127,16 +127,18 @@ from pyrenew import latent, deterministic, metaclass
import jax.numpy as jnp
import numpyro.distributions as dist

inf_hosp_int = deterministic.DeterministicPMF(inf_hosp_int, name="inf_hosp_int")
inf_hosp_int = deterministic.DeterministicPMF(
inf_hosp_int, name="inf_hosp_int"
)

hosp_rate = metaclass.DistributionalRV(
dist=dist.LogNormal(jnp.log(0.05), 0.1),
name="IHR",
)

latent_hosp = latent.HospitalAdmissions(
infection_to_admission_interval=inf_hosp_int,
infect_hosp_rate_dist=hosp_rate,
infection_to_admission_interval_rv=inf_hosp_int,
infect_hosp_rate_rv=hosp_rate,
)
```

@@ -173,12 +175,12 @@ Notice all the components are `RandomVariable` instances. We can now build the m
```{python}
# | label: init-model
hosp_model = model.HospitalAdmissionsModel(
latent_infections=latent_inf,
latent_admissions=latent_hosp,
I0=I0,
gen_int=gen_int,
Rt_process=rtproc,
observation_process=obs,
latent_infections_rv=latent_inf,
latent_hosp_admissions_rv=latent_hosp,
I0_rv=I0,
gen_int_rv=gen_int,
Rt_process_rv=rtproc,
hosp_admission_obs_process_rv=obs,
)
```

@@ -208,7 +210,7 @@ axs[0].plot(sim_data.Rt)
axs[0].set_ylabel("Rt")

# Infections plot
axs[1].plot(sim_data.sampled_admissions)
axs[1].plot(sim_data.sampled_observed_hosp_admissions)
axs[1].set_ylabel("Infections")
axs[1].set_yscale("log")

@@ -230,7 +232,7 @@ import jax
hosp_model.run(
num_samples=2000,
num_warmup=2000,
observed_admissions=dat["daily_hosp_admits"].to_numpy(),
observed_hosp_admissions=dat["daily_hosp_admits"].to_numpy(),
rng_key=jax.random.PRNGKey(54),
mcmc_args=dict(progress_bar=False),
)
@@ -244,7 +246,7 @@ We can use the `plot_posterior` method to visualize the results[^capture]:
# | label: fig-output-hospital-admissions
# | fig-cap: Hospital Admissions posterior distribution
out = hosp_model.plot_posterior(
var="predicted_admissions",
var="observed_hosp_admissions",
ylab="Hospital Admissions",
obs_signal=dat["daily_hosp_admits"].to_numpy(),
)
@@ -268,7 +270,7 @@ dat_w_padding = np.hstack((np.repeat(np.nan, days_to_impute), dat_w_padding))
hosp_model.run(
num_samples=2000,
num_warmup=2000,
observed_admissions=dat_w_padding,
observed_hosp_admissions=dat_w_padding,
rng_key=jax.random.PRNGKey(54),
mcmc_args=dict(progress_bar=False),
padding=days_to_impute, # Padding the model
@@ -281,7 +283,7 @@ And plotting the results:
# | label: fig-output-admissions-with-padding
# | fig-cap: Hospital Admissions posterior distribution
out = hosp_model.plot_posterior(
var="predicted_admissions",
var="observed_hosp_admissions",
ylab="Hospital Admissions",
obs_signal=dat_w_padding,
)
@@ -343,18 +345,18 @@ Notice that the instance's `nweeks` and `len` members are passed during construc
```{python}
# | label: latent-hosp-weekday
latent_hosp_wday_effect = latent.HospitalAdmissions(
infection_to_admission_interval=inf_hosp_int,
infect_hosp_rate_dist=hosp_rate,
weekday_effect_dist=weekday_effect,
infection_to_admission_interval_rv=inf_hosp_int,
infect_hosp_rate_rv=hosp_rate,
weekday_effect_rv=weekday_effect,
)

hosp_model_weekday = model.HospitalAdmissionsModel(
latent_infections=latent_inf,
latent_admissions=latent_hosp_wday_effect,
I0=I0,
gen_int=gen_int,
Rt_process=rtproc,
observation_process=obs,
latent_infections_rv=latent_inf,
latent_hosp_admissions_rv=latent_hosp_wday_effect,
I0_rv=I0,
gen_int_rv=gen_int,
Rt_process_rv=rtproc,
hosp_admission_obs_process_rv=obs,
)
```

@@ -365,7 +367,7 @@ Running the model (with the same padding as before):
hosp_model_weekday.run(
num_samples=2000,
num_warmup=2000,
observed_admissions=dat_w_padding,
observed_hosp_admissions=dat_w_padding,
rng_key=jax.random.PRNGKey(54),
mcmc_args=dict(progress_bar=False),
padding=days_to_impute,
@@ -378,7 +380,7 @@ And plotting the results:
# | label: fig-output-admissions-padding-and-weekday
# | fig-cap: Hospital Admissions posterior distribution
out = hosp_model_weekday.plot_posterior(
var="predicted_admissions",
var="observed_hosp_admissions",
ylab="Hospital Admissions",
obs_signal=dat_w_padding,
)
20 changes: 10 additions & 10 deletions model/docs/extending_pyrenew.qmd
Original file line number Diff line number Diff line change
@@ -63,11 +63,11 @@ With all the components defined, we can build the model:
```{python}
# | label: build1
model0 = RtInfectionsRenewalModel(
gen_int=gen_int,
I0=I0,
latent_infections=latent_infections,
Rt_process=rt,
observation_process=None,
gen_int_rv=gen_int,
I0_rv=I0,
latent_infections_rv=latent_infections,
Rt_process_rv=rt,
infection_obs_process_rv=None,
)
```

@@ -239,11 +239,11 @@ latent_infections2 = InfFeedback(
)

model1 = RtInfectionsRenewalModel(
gen_int=gen_int,
I0=I0,
latent_infections=latent_infections2,
Rt_process=rt,
observation_process=None,
gen_int_rv=gen_int,
I0_rv=I0,
latent_infections_rv=latent_infections2,
Rt_process_rv=rt,
infection_obs_process_rv=None,
)

# Sampling and fitting model 0 (with no obs for infections)
14 changes: 7 additions & 7 deletions model/docs/getting_started.qmd
Original file line number Diff line number Diff line change
@@ -113,11 +113,11 @@ With these five pieces, we can build the basic renewal model as an instance of t
```{python}
# | label: model-creation
model1 = RtInfectionsRenewalModel(
gen_int=gen_int,
I0=I0,
Rt_process=rt_proc,
latent_infections=latent_infections,
observation_process=observation_process,
gen_int_rv=gen_int,
I0_rv=I0,
Rt_process_rv=rt_proc,
latent_infections_rv=latent_infections,
infection_obs_process_rv=observation_process,
)
```

@@ -167,7 +167,7 @@ axs[0].plot(sim_data.Rt)
axs[0].set_ylabel("Rt")

# Infections plot
axs[1].plot(sim_data.sampled_infections)
axs[1].plot(sim_data.sampled_observed_infections)
axs[1].set_ylabel("Infections")

fig.suptitle("Basic renewal model")
@@ -185,7 +185,7 @@ import jax
model1.run(
num_warmup=2000,
num_samples=1000,
observed_infections=sim_data.sampled_infections,
observed_infections=sim_data.sampled_observed_infections,
rng_key=jax.random.PRNGKey(54),
mcmc_args=dict(progress_bar=False),
)
32 changes: 16 additions & 16 deletions model/docs/pyrenew_demo.qmd
Original file line number Diff line number Diff line change
@@ -37,15 +37,15 @@ import numpyro.distributions as dist
from pyrenew.process import SimpleRandomWalkProcess
```

To understand the simple random walk process underlying the sampling within the renewal process model, we first examine a single random walk path. Using the `sample` method from an instance of the `SimpleRandomWalkProcess` class, we first create an instance of the `SimpleRandomWalkProcess` class with a normal distribution of mean = 0 and standard deviation = 0.0001 as its input. Next, the `with` statement sets the seed for the random number generator for the duration of the block that follows. Inside the `with` block, the `q_samp = q.sample(duration=100)` generates the sample instance over a duration of 100 time units. Finally, this single random walk process is visualized using `matplot.pyplot` to plot the exponential of the sample instance.
To understand the simple random walk process underlying the sampling within the renewal process model, we first examine a single random walk path. Using the `sample` method from an instance of the `SimpleRandomWalkProcess` class, we first create an instance of the `SimpleRandomWalkProcess` class with a normal distribution of mean = 0 and standard deviation = 0.0001 as its input. Next, the `with` statement sets the seed for the random number generator for the n_timepoints of the block that follows. Inside the `with` block, the `q_samp = q.sample(n_timepoints=100)` generates the sample instance over a n_timepoints of 100 time units. Finally, this single random walk process is visualized using `matplot.pyplot` to plot the exponential of the sample instance.

```{python}
# | label: fig-randwalk
# | fig-cap: Random walk example
np.random.seed(3312)
q = SimpleRandomWalkProcess(dist.Normal(0, 0.001))
with seed(rng_seed=np.random.randint(0, 1000)):
q_samp = q.sample(duration=100)
q_samp = q.sample(n_timepoints=100)

plt.plot(np.exp(q_samp[0]))
```
@@ -112,8 +112,8 @@ inf_hosp_int = DeterministicPMF(
)

latent_admissions = HospitalAdmissions(
infection_to_admission_interval=inf_hosp_int,
infect_hosp_rate_dist=DistributionalRV(
infection_to_admission_interval_rv=inf_hosp_int,
infect_hosp_rate_rv=DistributionalRV(
dist=dist.LogNormal(jnp.log(0.05), 0.05), name="IHR"
),
)
@@ -131,12 +131,12 @@ The `HospitalAdmissionsModel` is then initialized using the initial conditions j
```{python}
# Initializing the model
hospmodel = HospitalAdmissionsModel(
gen_int=gen_int,
I0=I0,
latent_admissions=latent_admissions,
observation_process=admissions_process,
latent_infections=latent_infections,
Rt_process=Rt_process,
gen_int_rv=gen_int,
I0_rv=I0,
latent_hosp_admissions_rv=latent_admissions,
hosp_admission_obs_process_rv=admissions_process,
latent_infections_rv=latent_infections,
Rt_process_rv=Rt_process,
)
```

@@ -151,13 +151,13 @@ x
Visualizations of the single model output show (top) infections over the 30 time steps, (middle) hospital admissions over the 30 time steps, and (bottom)

```{python}
#| label: fig-hosp
#| fig-cap: Infections
# | label: fig-hosp
# | fig-cap: Infections
fig, ax = plt.subplots(nrows=3, sharex=True)
ax[0].plot(x.latent_infections)
ax[0].set_ylim([1/5, 5])
ax[1].plot(x.latent_admissions)
ax[2].plot(x.sampled_admissions, 'o')
ax[0].set_ylim([1 / 5, 5])
ax[1].plot(x.latent_hosp_admissions)
ax[2].plot(x.sampled_observed_hosp_admissions, "o")
for axis in ax[:-1]:
axis.set_yscale("log")
```
@@ -169,7 +169,7 @@ To fit the `hospmodel` to the simulated data, we call `hospmodel.run()`, an MCMC
hospmodel.run(
num_warmup=1000,
num_samples=1000,
observed_admissions=x.sampled_admissions,
observed_hosp_admissions=x.sampled_observed_hosp_admissions,
rng_key=jax.random.PRNGKey(54),
mcmc_args=dict(progress_bar=False),
)
22 changes: 22 additions & 0 deletions model/pyproject.toml
Original file line number Diff line number Diff line change
@@ -32,6 +32,28 @@ pytest-cov = "^5.0.0"
pytest-mpl = "^0.17.0"
numpydoc = "^1.7.0"

[tool.numpydoc_validation]
checks = [
"GL03",
"GL08",
"SS01",
"PR03",
"PR04",
"PR07",
"RT01"
]
ignore = [
"ES01",
"SA01",
"EX01",
"SS06",
"RT05"
]
exclude = [ # don't report on objects that match any of these regex
'\.undocumented_method$',
'\.__repr__$',
'\.__call__$'
]

[build-system]
requires = ["poetry-core"]
6 changes: 3 additions & 3 deletions model/src/pyrenew/deterministic/nullrv.py
Original file line number Diff line number Diff line change
@@ -129,7 +129,7 @@ def validate() -> None:

def sample(
self,
predicted: ArrayLike,
mu: ArrayLike,
obs: ArrayLike | None = None,
name: str | None = None,
**kwargs,
@@ -139,8 +139,8 @@ def sample(

Parameters
----------
predicted : ArrayLike
Rate parameter of the Poisson distribution.
mu : ArrayLike
Unused parameter, represents mean of non-null distributions
obs : ArrayLike, optional
Observed data. Defaults to None.
name : str, optional
Loading