diff --git a/docs/quickstart.md b/docs/quickstart.md
index c19edff..e82e5a5 100644
--- a/docs/quickstart.md
+++ b/docs/quickstart.md
@@ -2,7 +2,7 @@
## Set up a dimensionless ODE model
-Load pySODM's [`ODE` class](models.md) and define your model. As an example, we'll set up a simple Susceptible-Infectious-Removed (SIR) disease model, schematically represented as follows,
+To set up a simple Susceptible-Infectious-Removed (SIR) disease model, schematically represented as follows,
@@ -16,6 +16,8 @@ N &=& S + I + R, \\
\end{eqnarray}
```
+Load pySODM's [`ODE`](models.md) class, and define a new class inheriting the `ODE` object to define your model. Minimally, you must provide: 1) A list containing the names of your model's states, 2) a list containing the names of your model's parameters, 3) an `integrate` function integrating your model's states. To learn more about the `ODE` class' formatting, check out the [`ODE`](models.md) class description.
+
```python
# Import the ODE class
from pySODM.models.base import ODE
@@ -45,7 +47,7 @@ To initialize the model, provide a dictionary containing the initial values of t
model = SIR(states={'S': 1000, 'I': 1}, parameters={'beta': 0.35, 'gamma': 5})
```
-Simulate the model using the `sim()` method. The solver method and tolerance of `scipy.solve_ivp()`, as well as the timesteps included in the output can be tailored. pySODM supports the use of `datetime` types as timesteps.
+Simulate the model using its `sim()` method. pySODM supports the use of `datetime` types as timesteps.
```python
# Timesteps
@@ -63,12 +65,9 @@ In some situations, the use of a discrete timestep with a fixed length may be pr
```python
# Use a discrete timestepper with step size 1
out = model.sim(121, tau=1)
-
-# Is equivalent to (`tau` overwrites `method`!):
-out = model.sim(121, method='RK45', rtol='1e-4', tau=1)
```
-Results are sent to an `xarray.Dataset`, for more information on indexing and selecting data using `xarray`, [see](https://docs.xarray.dev/en/stable/user-guide/indexing.html).
+Simulations are sent to an `xarray.Dataset`, for more information on indexing and selecting data using `xarray`, [see](https://docs.xarray.dev/en/stable/user-guide/indexing.html).
```bash
Dimensions: (time: 122)
@@ -81,11 +80,12 @@ Data variables:
```
The simulation above results in the following trajectories for the number of susceptibles, infected and recovered individuals.
+
## Set up an ODE model with a labeled dimension
-To transform all our SIR model's states in 1D vectors, referred to as a *dimension* throughout the documentation, add the `dimensions` keyword to the model declaration. Here, we add a dimension representing the age groups individuals belong to.
+To transform all our SIR model's states in 1D vectors, referred to as a *dimension* throughout the documentation, add the `dimensions` keyword to the model declaration. In the example below, we add a dimension representing the age groups individuals belong to. This turns all model states in the `integrate` function into 1D vectors.
```python
# Import the ODE class
@@ -95,7 +95,7 @@ from pySODM.models.base import ODE
class stratified_SIR(ODE):
states = ['S','I','R']
- parameters = ['beta', 'gamma]
+ parameters = ['beta', 'gamma']
dimensions = ['age_groups']
@staticmethod
@@ -111,7 +111,7 @@ class stratified_SIR(ODE):
return dS, dI, dR
```
-When initializing the model, additionally provide a dictionary containing coordinates for every dimension declared previously. In our example, we'll declare four age groups: 0-5, 5-15, 15-65 and 65-120 year olds. **All** model states are now 1D vectors of shape `(4,)`.
+When initializing the model, provide a dictionary containing coordinates for every dimension declared previously. In the example below, we'll declare four age groups: 0-5, 5-15, 15-65 and 65-120 year olds. **All** model states are now 1D vectors of shape `(4,)`.
```python
model = stratified_SIR(states={'S': 1000*np.ones(4), 'I': np.ones(4)},
@@ -139,6 +139,8 @@ pySODM allows model states to have different coordinates and thus different size
+In the example below, the states S, I and R will be 1D vectors in the `integrate` function, while the states S_v and I_v will be scalars.
+
```python
class ODE_SIR_SI(ODE):
"""
@@ -182,7 +184,7 @@ out = model.sim(120)
print(out)
```
-results in the following `xarray.Dataset`. Here it can be seen that the S, I and R state have `age_group` and `time` as dimensions and their shape is thus `(4,)`, while S_v and I_v only have `time` as dimension and their shape is thus `(1,)`.
+results in the following `xarray.Dataset`,
```bash
@@ -197,10 +199,11 @@ Data variables:
S_v (time) float64 1e+06 1e+06 1e+06 ... 8.635e+05 8.562e+05
I_v (time) float64 2.0 1.798 1.725 ... 1.292e+05 1.365e+05 1.438e+05
```
+Here, S, I and R state have `age_group` and `time` as dimensions and their shape is thus `(4,)`, while S_v and I_v only have `time` as dimension and their shape is thus `(1,)`.
## Set up a dimensionless stochastic jump process Model
-To stochastically simulate the simple SIR model, the `JumpProcess` class is loaded and two functions `compute_rates` and `apply_transitionings` must be defined. For a detailed mathematical description of implenting models using Gillespie's tau-leaping method (example of a jump proces), we refer to the [manuscript](https://arxiv.org/abs/2301.10664). The rates dictionary defined in `compute_rates` contains the rates of the possible transitionings in the system. These are contained in a list because a state may have multiple possible transitionings. Further, the transitioning for state `I` is a float and this should be a `numpy.float`. I'm still figuring out if I want to introduce the overhead to check and correct this behavior (likely not).
+To stochastically simulate the simple SIR model, pySODM's `JumpProcess` class is loaded and two functions, `compute_rates` and `apply_transitionings`, must be defined instead of one. For a detailed mathematical description of implenting models using Gillespie's tau-leaping method (example of a jump proces), we refer to the [peer-reviewed paper](https://www.sciencedirect.com/science/article/pii/S1877750323002089). The rates dictionary defined in `compute_rates` contains the rates of the possible transitionings in the system. These are contained within a list so that a state may have multiple transitionings. To learn more about the `JumpProcess` class' formatting, check out the [`JumpProcess`](models.md) class description.
```python
# Import the ODE class
@@ -238,7 +241,7 @@ class SIR(JumpProcess):
### Draw functions
-The simulation functions of the `ODE` and `JumpProcess` classes (`sim()`) can be used to perform {math}`N` repeated simulations by using the optional argument `N`. A *draw function* can be used to instruct the model on how to alter the value of a model parameters during consecutive model runs, thereby offering a powerful tool for sensitivity analysis and uncertainty propagation.
+The simulation functions of the `ODE` and `JumpProcess` classes can be used to perform {math}`N` repeated simulations by using the optional argument `N`. A *draw function* can be used to instruct the model on how to alter the value of a model parameters during consecutive model runs, thereby offering a powerful tool for sensitivity analysis and uncertainty propagation.
Draw functions **always** take the dictionary of model parameters, `parameters` as their first argument, input checks are used to ensure you provide the correct name and type. This can be followed an arbitrary number of user-defined parameters, which are supplied to the `sim()` function by using its `draw_function_kwargs` argument. The output of a draw function is **always** the dictionary of model parameters, without alterations to the dictionaries keys (no new parameters introduced or parameters deleted). In the example below, we attempt to draw a model parameter `gamma` randomly from 0 to 5,
@@ -293,7 +296,7 @@ Parameters can also be varied through the course of one simulation using a *time
```python
def vary_my_parameter(t, states, param, an_additional_parameter):
- """A function to vary a model parameter during a simulation"""
+ """A function to implement a time-dependency on a parameter"""
# Do any computation
param = ...
@@ -305,6 +308,6 @@ When initialising the model, all we need to do is use the `time_dependent_parame
```python
model = SIR(states={'S': 1000, 'I': 1},
- parameters={'beta': 0.35, 'gamma': 5, 'an_additional_parameter': anything_you_want},
+ parameters={'beta': 0.35, 'gamma': 5, 'an_additional_parameter': any_datatype_you_want},
time_dependent_parameters={'beta': vary_my_parameter})
```
\ No newline at end of file