Starting the monty section

2009fluUK.qmd

@@ -0,0 +1 @@
+# 2009 pandemic in the UK

_quarto.yml

@@ -21,10 +21,17 @@ book:
     - part: "Monty"
       chapters:
       - monty.qmd
+      - model.qmd
+      - montyDSL.qmd
+      - composability.qmd
+      - samplers.qmd
+      - runners.qmd
+      - advanced_monty.qmd
     - part: "Odin & Monty"
       chapters:
       - inference.qmd
       - differentiability.qmd
+      - 2009fluUK.qmd
     - references.qmd
   appendices:
     - installation.qmd

advanced_monty.qmd

@@ -0,0 +1,14 @@
+# Advanced monty
+
+```{r}
+#| include: false
+source("common.R")
+```
+
+## Observers
+
+## Working with other packages
+
+Exporting chains to work with other packages
+
+## Nested models

composability.qmd

@@ -0,0 +1,29 @@
+# Composability {#sec-monty-combine}
+
+```{r}
+#| include: false
+source("common.R")
+```
+
+```{r}
+library(monty)
+```
+
+```{r}
+# A simple example; a model that contains something of interest,
+# and a simple prior from monty_dsl
+likelihood <- monty_example("banana")
+prior <- monty_dsl({
+  alpha ~ Normal(0, 1)
+  beta ~ Normal(0, 10)
+})
+posterior <- likelihood + prior
+posterior
+
+# # The same thing, more explicitly:
+# monty_model_combine(likelihood, prior)
+# 
+# # Control properties of the combined model:
+# monty_model_combine(likelihood, prior,
+#                     monty_model_properties(has_gradient = FALSE))
+```

model.qmd

@@ -0,0 +1,114 @@
+# About (monty) Models {#sec-monty-models}
+
+```{r}
+#| include: false
+source("common.R")
+```
+
+In this book, we explore two fundamentally different yet interconnected approaches to modelling: **dynamical systems** and **statistical models**. To illustrate these differences, let’s begin with a familiar (but maybe dated!) metaphor: the game [*SimCity*](https://en.wikipedia.org/wiki/SimCity).
+
+## Models vs. Models: dynamical and statistical perspectives
+
+### SimCity and dynamical systems
+
+Imagine a city model, as in *SimCity*, where thousands of virtual inhabitants follow routines, interact, and respond to changes like new buildings or natural disasters. These elements evolve according to predefined rules and parameters, similar to how **dynamical systems** simulate real-world processes over time.
+
+In a [**dynamical model**](https://en.wikipedia.org/wiki/Dynamical_system), we track changes in a system's "state" over time. The state is a summary of key information at a specific time point; for example, in an epidemiological model, the state might include the numbers of susceptible, infected, and recovered individuals. More detailed models may add variables like age, location, health risks, and symptoms, allowing us to simulate interventions, such as vaccination campaigns, and explore their potential impacts. There exist many formalisms to describe dynamical systems, incorporating e.g. specific input and output functions in control theory but the central element of these systems is this notion of the "state", a mathematical object summarising the system at a timepoint. 
+
+The `odin` package, which supports differential and difference equations, is an effective tool for building dynamical models, enabling users to simulate time-evolving systems. Dynamical models are well-suited for exploring how specific scenarios or interventions impact a system over time, making them particularly useful for modelling real-world phenomena in a structured way.
+
+TODO: simple odin deterministic model with discussion about results being a time indexed suite of "states"
+
+### Statistical models and uncertainty
+
+[**Statistical models**](https://en.wikipedia.org/wiki/Statistical_model), on the other hand, focus on capturing uncertainty in outcomes. Instead of simulating the system’s evolution over time, they represent possible states of the world probabilistically. Statistical models summarise data patterns or help make predictions based on uncertainties, such as estimating the range of future cases of a disease or identifying risk factors in a population.
+
+The `monty` package facilitates statistical modelling, providing tools for working with probability distributions and Monte Carlo methods. Where `odin` enables us to explore the dynamics of changing systems, `monty` is ideal for making probabilistic inferences and assessing variability across potential outcomes without necessarily modelling time-evolving dynamics.
+
+TODO: simple gaussian model - discussion about parameter being actually the sampling space and explain that it connects with Bayesian inference, show that the important component of the model is the parameters and domain support + density function, maybe shows that the model can be built in different way (R function, monty_model(), DSL).
+
+## Bridging dynamical systems and statistical models
+
+While dynamical systems and statistical models serve distinct purposes, they are not strictly separate. In fact, they can be connected through two powerful approaches that bring probabilistic reasoning into dynamical systems: **stochastic processes** and **Bayesian inference**.
+
+### Stochastic processes: adding uncertainty to dynamical models
+
+A natural way to bridge dynamical and statistical models is by introducing [**stochastic processes**](https://en.wikipedia.org/wiki/Stochastic_process). Traditional dynamical models use fixed parameters and deterministic rules to evolve a system’s state over time. However, many real-world systems have inherent randomness or uncertainty that deterministic models cannot capture.
+
+In a **stochastic process**, the system’s state is no longer a single deterministic value but a collection of potential outcomes, each weighted by a probability. This probabilistic view enables us to model fluctuations and uncertainties within the dynamical framework, capturing both the system’s evolution and the uncertainty in each state. Stochastic processes are thus a natural extension of dynamical models, adding an extra layer of realism by treating system evolution as inherently uncertain.
+
+The `odin` package provides an intuitive framework for writing a class of stochastic systems, making it easier to define models that incorporate randomness in their evolution over time. The `dust` package complements `odin` by enabling efficient large-scale simulations, allowing users to capture and analyse the uncertainty inherent to these systems through repeated runs. Together, `odin` and `dust` offer a powerful toolkit for developing and exploring stochastic models that reflect the complexity and variability of real-world dynamics. 
+
+TODO simple dust example or just link with relevant dust section in the book
+
+### Bayesian inference: statistical modelling of model parameters
+
+**Bayesian inference** is another approach to linking dynamical and statistical models by treating model parameters as random variables rather than fixed values. This introduces a probability distribution over possible parameter values, making parameter estimation a statistical problem.
+
+Using Bayes’ theorem, Bayesian inference combines:
+
+  - **The likelihood**: the probability of observing the data given specific parameter values, and
+  - **The prior distribution**: our initial assumptions about the parameter values before observing the data.
+
+The result is the **posterior distribution** of a parameter $\theta$ given data $y$:
+
+$$
+p(\theta | y) = \frac{p(y|\theta) p(\theta)}{p(y)}
+$$
+
+where:
+
+  - $p(y|\theta)$ is the likelihood,
+  - $p(\theta)$ is the prior distribution, and
+  - $p(y)$ is a normalising constant to ensure the posterior integrates to 1 over all $\theta$ values.
+
+Bayesian inference updates our beliefs about parameters based on observed data, yielding a posterior distribution that reflects both model assumptions and data influence. This probabilistic framework around a dynamical system’s parameters creates a flexible statistical model that adapts as new data becomes available.
+
+In this way, **stochastic processes** add an uncertainty dimension to the state of dynamical systems, while **Bayesian inference** applies probabilistic modelling to the parameters, merging dynamical and statistical perspectives. Together, these methods provide a rich toolkit for modelling complex systems that evolve over time and accounting for both randomness and uncertainty.
+
+## Parameters and model complexity
+
+### Parameters in dynamical systems
+
+In dynamical models, parameters define the structure and dynamics of the scenario. For example, a detailed epidemiological model may need parameters for transmission rates, contact patterns, or intervention schedules. These inputs allow us to examine "what-if" scenarios, enabling decision-makers to predict and manage changes in real-time.
+
+### Parameters in statistical models
+
+Statistical models, meanwhile, use parameters to define probability distributions over possible outcomes rather than scenario-specific details. Here, parameters describe uncertainties, often quantifying risks, making predictions, or summarising patterns in data.
+
+The `odin` package allows us to build and explore the behaviours of dynamical models over time, while `monty` enables us to quantify uncertainties and make probabilistic inferences without simulating detailed dynamics.
+
+## Bayesian or non-Bayesian modelling
+
+Many statistical modelling tools, such as `stan`, `mcstate`, and `BayesTools`, are inherently Bayesian. However, the `monty` package is **Bayesian-agnostic**: it does not require a Bayesian interpretation. This flexibility allows users to work with probabilistic models either within a Bayesian framework or outside of it, depending on the needs of the analysis.
+
+## Probability densities: normalised vs. unnormalised
+
+A key concept in Bayesian inference and Monte Carlo methods is the distinction between **normalised** and **unnormalised** probability densities.
+
+### Unnormalised density
+
+In Bayesian inference, rather than using the "full" Bayesian as above, we are often working with defining posterior density as being proportional to a product of densities:
+
+$$
+p(\theta | y) \propto p(y|\theta) p(\theta)
+$$
+
+where:
+
+  - $p(y|\theta)$ is the likelihood, and
+  - $p(\theta)$ is the prior distribution, as above.
+
+Note that it does not involve anymore $p(y)$ the normalising constant. The reason is that since calculating $p(y)$ can be very difficult, we often work with the **unnormalised posterior density** $p(y|\theta)p(\theta)$. This unnormalised form is sufficient for many Monte Carlo methods, where only relative densities matter.
+
+### Normalised density
+
+A **normalised density** integrates to 1 over its entire parameter space. This is necessary for direct probability interpretations and for certain Bayesian methods, like model comparison using [Bayes factors](https://en.wikipedia.org/wiki/Bayes_factor).
+
+## Monte Carlo methods with unnormalised densities
+
+Monte Carlo algorithms, such as **Metropolis-Hastings** and **Importance Sampling**, often operate using unnormalised densities, focusing on relative probabilities rather than absolute values. This makes them efficient for high-dimensional problems where calculating a normalising constant is untractable.
+
+## Linking to Physics: partition functions and probability densities
+
+Normalisation in probability densities has parallels in physics. In statistical mechanics, the **partition function** $Z$ normalises probabilities over all possible states, much like normalising densities in Bayesian inference. This connection highlights why algorithms like **Metropolis-Hastings** use unnormalised densities: they mirror physical systems where absolute energies are less relevant than energy differences.

monty.qmd

@@ -1,16 +1,95 @@
-# Monty
+# Getting started with monty
 
-Before showing how to fit odin models to data, we'll focus for a bit on monty itself.  If you are anxious to get on and fit the model from @sec-data-filter, you might skip ahead to @sec-inference, where we resume this.
+```{r}
+#| include: false
+source("common.R")
+```
 
-There are a bunch of things we want to cover in this chapter so likely it will be split into a few:
+The `monty` R package is designed for modular work with statistical distributions, leveraging Monte Carlo methods for sampling. It enables the construction of increasingly complex models through four approaches:
 
-* defining a model in monty, both easily with `monty_function` and less easily with `monty_model`
-* fitting models with MCMC
-* writing priors with the `monty_dsl`
-* using observers to follow along with the model (if we can come up with a good example)
-* working with nested models (once we write the code that would allow this)
-* running MCMC chains in parallel
-* exporting chains to work with other packages
-* different samplers
+1. Simple R code,
+2. A dedicated domain-specific language (DSL),
+3. Integration of `odin` models,
+4. Composition of two `monty` models.
 
-Hmm, there really is quite a bit of ground to cover here!
+The most basic method to define a `monty` model is through a simple R function, as shown in this chapter. More advanced examples are covered in the @sec-monty-models in particular in the context of Bayesian statistics. @sec-monty-dsl introduces the monty DSL for more versatile model building using a small probabilistic DSL similar to BUGs. @sec-monty-combine explains how to compose new monty models by combining two existing ones. The last part of this book starting from @sec-inference demonstrates how to incorporate `odin` models into `monty` workflows.
+
+## A simple example
+
+```{r}
+library(monty)
+```
+
+We can define a simple Gaussian [mixture model](https://en.wikipedia.org/wiki/Mixture_model) of two sub-populations using the `monty_model_function()` from the package. The model represents the distribution of a quantity $l$, sampled with a probability $p$ from a normal distribution of mean $m_{1}$ and with a probability $1-p$ from another normal distribution of mean $m_{2}$. Both subpopulations have variance equal to 1. 
+
+To build our `monty` model, we start by defining an R function returning the log-density of our statistical model. This function has for arguments $l$, the quantity of interest and the three parameters $p$, $m_{1}$ and $m_{1}$ defining the two Normally-distributed subpopulations.
+
+```{r}
+fn <- function(l, p, m1, m2) {
+  log(p*dnorm(l,mean = m1) + (1-p)*dnorm(l, mean = m2))
+}
+```
+
+Assuming that the population is 75% from subpopulation 1 (normally distributed with mean 3) and 25% from subpopulation 2 (also normally distributed but with mean 7), we can build our `monty` model by indicating that the parameters of the subpopulations are fixed at these given values.
+
+```{r}
+l_distribution <- monty_model_function(fn,
+                                       fixed = list(p = 0.75, m1 = 3, m2=7))
+```
+
+We have just created a `monty` model. 
+
+```{r}
+l_distribution
+```
+
+We can plot the density of our model for values of $l$ between 0 and 10 and eye-check that we can see the two subpopulations in the correct proportions and with the 3 and 7 as modes.
+
+```{r}
+l <- seq(from = 0, to = 10, by = 0.1)
+plot(l,
+     exp(Vectorize(l_distribution$density)(l)),
+     type = "l",
+     ylab = "density")
+```
+
+## Sampling from our example distribution
+
+We want now to sample from this model, using the `monty_sample()` function. For this we need to tell `monty` which sampler we want to use to explore our distribution. There are a variety of samplers avalaible and you can learn about them in @sec-samplers. One of the simplest is the random walk [Metropolis-Hastings](https://en.wikipedia.org/wiki/Metropolis%E2%80%93Hastings_algorithm) algorithm that should work almost out of the box (though not necesseraly efficiently) in most cases.
+
+The random walk sampler uses a variance-covariance (VCV) matrix to guide its exploration, determining the 'jump' from the current point to the next in a random walk by drawing from a multivariate normal distribution parameterised by this matrix. For a single-parameter model, we use a 1x1 matrix (`matrix(1)`) as our VCV matrix.
+
+The choice of the VCV matrix is critical for the sampler's efficiency, especially in more complex cases where the tuning of this matrix can significantly affect performance. A well-chosen VCV matrix optimises moving accross the parameter space, making the random walk sampler more effective in exploring the distribution of interest.
+
+```{r}
+sampler <- monty_sampler_random_walk(matrix(1))
+```
+
+We are now ready to sample from our distribution using `monty_sample()` with 2000 samples, starting our MCMC chain at the value 3 and running 4 chains simultaneously.
+
+```{r}
+samples <- monty_sample(l_distribution, sampler, 2000, initial = 3, n_chains = 4)
+```
+
+We can visualise our 4 chains.
+```{r}
+matplot(samples$density, type = "l", lty = 1,
+        xlab = "sample", ylab = "log posterior density", col = "#00000055")
+```
+
+We can also check that our samples are correctly representing our distribution:
+```{r}
+hist(samples$pars["l",,], breaks=100)
+```
+
+## Connection of `monty` with other software
+
+`monty` includes a simple probabilistic domain-specific language (DSL) that is inspired by languages of the BUGS family such as `stan` and [Statistical Rethinking](https://xcelab.net/rm/). It is designed to make some tasks a bit easier, particularly when defining priors for your model.  We expect that this DSL is not sufficiently advanced to represent most interesting models but it may get more clever and flexible in the future. In particular we do not expect the DSL to be useful in writing likelihood functions for comparison to data; we expect that if your model is simple enough for this you would be better off using `stan` or some similarly flexible system.
+
+*mention SSM, dr jacoby, mcstate, BayesTools
+
+## Going further
+
+We've discussed a little bit about the philosophy of `monty` and built a simple model using an R function with `monty_function`. The example introduces the basics of defining and sampling from a `monty` model.
+
+For more advanced applications, refer to @sec-monty-dsl for constructing models using `monty`'s domain-specific language, which enables greater flexibility for complex structures. To combine multiple models, see @sec-monty-combine, which demonstrates how to compose and integrate existing `monty` models. Lastly, @sec-inference illustrates how to incorporate `odin` models into `monty` workflows, expanding the package’s potential for Bayesian inference with more sophisticated, system-based models. Each section provides detailed examples to guide you in leveraging `monty`’s full capabilities.