diff --git a/mlss/Bayesian Machine Learning.jl b/mlss/Bayesian Machine Learning.jl index 9ed7281f..86bd3875 100644 --- a/mlss/Bayesian Machine Learning.jl +++ b/mlss/Bayesian Machine Learning.jl @@ -1,5 +1,5 @@ ### A Pluto.jl notebook ### -# v0.20.19 +# v0.20.21 #> [frontmatter] #> image = "https://github.com/bmlip/course/blob/v2/assets/figures/scientific-inquiry-loop-w-BML-eqs.png?raw=true" @@ -45,29 +45,6 @@ title("Bayesian Machine Learning") # ╔═╡ 6be2e966-4048-44d0-a37e-95060e3fe30b PlutoUI.TableOfContents() -# ╔═╡ 6a23df9e-d294-11ef-3ddf-a51d4cea00fc -md""" -## Preliminaries - -##### Goals - - * Introduction to Bayesian (i.e., probabilistic) modeling - -##### Materials - - * Mandatory - - * These lecture notes - - * Optional - - * Bishop pp. 68-74 (on the coin toss example) - * [Ariel Caticha - 2012 - Entropic Inference and the Foundations of Physics](https://github.com/bmlip/course/blob/main/assets/files/Caticha-2012-Entropic-Inference-and-the-Foundations-of-Physics.pdf), pp.35-44 (section 2.9, on deriving Bayes rule for updating probabilities) - - - -""" - # ╔═╡ eca027f8-40c9-4e53-85b5-d08b8fe9dd97 challenge_statement("Predicting a Coin Toss",header_level=1) @@ -86,7 +63,7 @@ What is the probability that heads comes up next? ##### Solution -Later in this lecture. +[Later in this lecture](#Challenge-Revisited:-Predicting-a-Coin-Toss). """ # ╔═╡ 6a24b9e4-d294-11ef-3ead-9d272fbf89be @@ -120,6 +97,14 @@ In principle, based on the model evaluation results, you may want to re-specify """ +# ╔═╡ 3cbe6b9b-4289-4724-ac3a-1b0313b3041b +md""" +Crucially, **all information processing in Bayesian machine learning is governed strictly by probability theory**. Consequently, parameter estimation, model evaluation, and model application are all carried out through probabilistic inference. +""" + +# ╔═╡ 55dec435-aa78-41ba-aad5-9d79ce292f42 +keyconcept("", "Bayesian machine learning is a subfield of machine learning that commits entirely to probability theory as the framework for all information-processing tasks. This is well justified, because probability theory is the optimal calculus for representing and manipulating states of knowledge.") + # ╔═╡ 6a24cee0-d294-11ef-35cb-71ab9ef935e5 md""" Next, we discuss these four stages in a bit more detail. @@ -166,14 +151,11 @@ md""" md""" ## 2. Parameter Estimation -You must now specify a likelihood function for the parameters from the data-generating distribution. Note that, for a given (i.e., *observed*) data set ``D=\{x_1,x_2,\dots,x_N\}`` with *independent* observations ``x_n``, the likelihood factorizes as +You must now specify a likelihood function for the parameters from the data-generating distribution. Note that, for a given (i.e., *observed*) data set ``D=\{x_1,x_2,\dots,x_N\}`` with *independent* observations ``x_n`` (drawn from the same model ``\theta``), the likelihood factorizes as ```math p(D|\theta) = \prod_{n=1}^N p(x_n|\theta)\,. ``` - -So, usually you select the data-generating distribution for one observation ``x_n`` and then use (in-)dependence assumptions to combine these models into a likelihood function for the model parameters. - """ # ╔═╡ 6a25379a-d294-11ef-3e07-87819f6d75cb @@ -194,40 +176,10 @@ Note that there's **no need for you to design some clever parameter estimation a """ -# ╔═╡ 6a2552ac-d294-11ef-08d6-179e068bc297 -md""" -This parameter estimation "recipe" works if the right-hand side (RHS) factors can be evaluated; the computational details can be quite challenging and this is what machine learning is about. - - -""" - -# ╔═╡ ce75e785-868f-4361-93f8-c582ac1b891b -keyconcept(" ", - md""" - - Bayesian Machine learning is EASY, apart from computational details :) - - """ -) - # ╔═╡ 6a2561c0-d294-11ef-124d-373846e3120c md""" ## 3. Model Evaluation -In the framework above, parameter estimation was executed by "perfect" Bayesian reasoning. So is everything settled now? - -""" - -# ╔═╡ 6a257020-d294-11ef-0490-e151934b2f42 -md""" -No, there appears to be one remaining problem: how good really were our assumptions ``p(D|\theta)`` and ``p(\theta)`` in the model specification phase? We want to "score" the model performance. - -""" - -# ╔═╡ 6a257f34-d294-11ef-2928-fbb800e81124 -md""" -Note that this question is only interesting in practice if we have alternative models to choose from. After all, if you don't have an alternative model, any value for the model evidence would still not lead you to switch to another model. - """ # ╔═╡ 6a25a11e-d294-11ef-1c51-09482dad86f2 @@ -275,24 +227,6 @@ p(m_k|D)&= \frac{p(m_k,D) }{p(D)} \\ ``` """) -# ╔═╡ 6a261278-d294-11ef-25a0-5572de58ad06 -md""" -You *can* evaluate the RHS of this equation since *you* selected the model priors ``p(m_k)``, the parameter priors ``p(\theta|m_k)``, and the likelihoods ``p(D|\theta,m_k)``. - -""" - -# ╔═╡ 6a26549a-d294-11ef-1f10-15c4d14ae41f -md""" -Note that, to evaluate the model posterior, you must calculate the **model evidence** ``p(D|m_k)``, which can be interpreted as a likelihood function for model ``m_k``. - -""" - -# ╔═╡ 6a262182-d294-11ef-23e9-ed45e1da9f46 -md""" -You can now compare posterior distributions ``p(m_k|D)`` for a set of models ``\{m_k\}`` and decide on the merits of each model relative to alternative models. This procedure is called **Bayesian model comparison**. - -""" - # ╔═╡ 6a2672d6-d294-11ef-1886-3195c9c7cfa9 md""" Again, **no need to invent a special algorithm for estimating the performance of your model**. Straightforward application of probability theory takes care of all that. @@ -309,89 +243,10 @@ keyconcept(" ", ) -# ╔═╡ 6a2664c6-d294-11ef-0a49-5192e17fb9ea -md""" - -Compare the calculations between parameter estimation and model evaluation -```math -\begin{align*} -p(\theta|D) &\propto p(D|\theta) p(\theta) \; &&\text{(parameter estimation)} \\ -p(m_k|D) &\propto p(D|m_k) p(m_k) \; &&\text{(model evaluation)} -\end{align*} -``` - -""" - -# ╔═╡ 6a26a31e-d294-11ef-2c2f-b349d0859a27 -md""" -With the (relative) performance evaluation scores of your model in hand, you could now re-specify your model (hopefully an improved model) and *repeat* the design process until the model performance score is acceptable (see the 4-step [Bayesian modeling process](#Bayesian-modeling-recipe) above). - -""" - -# ╔═╡ 6a269568-d294-11ef-02e3-13402d296391 -md""" -In principle, you could proceed with asking how good your choice for the candidate model set ``\mathcal{M}`` was. You would have to provide a set of alternative model sets ``\{\mathcal{M}_1,\mathcal{M}_2,\ldots,\mathcal{M}_M\}`` with priors ``p(\mathcal{M}_m)`` for each set and compute posteriors ``p(\mathcal{M}_m|D)``. And so forth ... - -""" - -# ╔═╡ 6a26b7bc-d294-11ef-03e7-2715b6f8dcc7 -md""" -### Bayes Factors - -""" - -# ╔═╡ 6a26f244-d294-11ef-0488-c1e4ec6e739d -md""" -As an aside, in the (statistics and machine learning) literature, performance comparison between two models is often reported by the [Bayes Factor](https://en.wikipedia.org/wiki/Bayes_factor), which is defined as the ratio of model evidences: - -```math -\begin{align*} -\mathrm{BF_{12}} \triangleq \frac{p(D|m_1)}{p(D|m_2)} -= \underbrace{\frac{p(m_1|D)}{p(m_2|D)}}_{\substack{\text{posterior} \\ \text{ratio}}} \cdot \underbrace{\frac{p(m_2)}{p(m_1)}}_{\substack{\text{prior} \\ \text{ratio}}} -\end{align*} -``` -""" - -# ╔═╡ 99db44c9-185c-4f39-ae5e-1a4cd751d980 -details("Prove this yourself, and click for solution", -md""" -```math -\begin{align*} -\mathrm{BF_{12}} &= \frac{p(D|m_1)}{p(D|m_2)} \\ -&= \frac{p(D,m_1)}{p(m_1)} \bigg/ \frac{p(D,m_2)}{p(m_2)} \\ -&= \frac{p(m_1|D) p(D)}{p(m_1)} \cdot \frac{p(m_2)}{p(m_2|D) p(D)} \\ -&= \frac{p(m_1|D)}{p(m_2|D)} \cdot \frac{p(m_2)}{p(m_1)} -\end{align*} -``` -""") - -# ╔═╡ d22f58ac-9f68-41cb-8e61-cf74d3692c44 -md""" -Hence, for equal model priors (``p(m_1)=p(m_2)=0.5``), the Bayes Factor reports the posterior probability ratio for the two models. - -In principle, any hard decision on which is the better model has to accept some *ad hoc* arguments. [Jeffreys (1961)](https://www.amazon.com/Theory-Probability-Classic-Physical-Sciences/dp/0198503687/ref=sr_1_1?qid=1663516628&refinements=p_27%3Athe+late+Harold+Jeffreys&s=books&sr=1-1&text=the+late+Harold+Jeffreys) advises to use the **log-Bayes factor**, - -```math -\mathrm{logBF}_{12} := ^{10}\log\frac{p(D|m_1)}{p(D|m_2)} \,, -``` - -to quantify evidence for preferring model ``m_1`` over ``m_2`` by the following interpretation: - -| ``\mathrm{logBF}_{12}`` | Evidence for ``m_1`` | -|:---------------------|:----------------------------| -| 0 to 0.5 | not worth mentioning | -| 0.5 to 1 | substantial | -| 1 to 2 | strong | -| >2 | decisive | - -""" - # ╔═╡ 6a2707e6-d294-11ef-02ad-31bf84662c70 md""" ## 4. Apply Model (Prediction) -Once we are satisfied with the evidence for a (trained) model, we can apply the model to our prediction/classification/etc task. - """ # ╔═╡ 6a271a56-d294-11ef-0046-add807cc0b4f @@ -421,80 +276,18 @@ In the last equation, the simplification ``p(x|\theta,D) = p(x|\theta)`` follows # ╔═╡ 6a273ae0-d294-11ef-2c00-9b3eaed93f6d md""" -Again, **no need to invent a special prediction algorithm**. Probability theory takes care of all that. The complexity of prediction is just computational, namely, how to carry out the marginalization over ``\theta``. - -""" - -# ╔═╡ 6a274948-d294-11ef-0563-1796b8883306 -md""" -Note that the application of the learned posterior ``p(\theta|D)`` not necessarily has to be a prediction task. We use it here as an example, but other applications (e.g., classification, regression etc.) are of course also possible. - -""" - -# ╔═╡ 6a275a52-d294-11ef-1323-9d83972f611a -md""" -### Prediction with multiple models - -When you have a posterior ``p(m_k|D)`` for the models, you don't *need* to choose one model for the prediction task. You can do prediction by **Bayesian model averaging**, which combines the predictive power from all models: - -```math -\begin{align*} -p(x|D) &= \sum_k \int p(x,\theta,m_k|D)\,\mathrm{d}\theta \\ - &= \sum_k \int p(x|\theta,m_k) \,p(\theta|m_k,D)\, p(m_k|D) \,\mathrm{d}\theta \\ - &= \sum_k \underbrace{p(m_k|D)}_{\substack{\text{model}\\\text{posterior}}} \cdot \int \underbrace{p(\theta|m_k,D)}_{\substack{\text{parameter}\\\text{posterior}}} \, \underbrace{p(x|\theta,m_k)}_{\substack{\text{data generating}\\\text{distribution}}} \,\mathrm{d}\theta -\end{align*} -``` - -""" - -# ╔═╡ 6a27684e-d294-11ef-040e-c302cdad714a -md""" -Alternatively, if you do need to work with one model (e.g. due to computational resource constraints), you can for instance select the model with largest posterior ``p(m_k|D)`` and use that model for prediction. This is called **Bayesian model selection**. - -""" - -# ╔═╡ 6a2777d0-d294-11ef-1ac3-add102c097d6 -md""" -Bayesian model averaging is the principal way to apply PT to machine learning. You don't throw away information by discarding lesser performant models, but rather use PT (marginalization of models) to compute - -```math -p(\text{what-I-am-interested-in} \,|\, \text{all available information})\,. -``` - -""" - -# ╔═╡ 6a278784-d294-11ef-11ae-65bd398910d5 -md""" -## We're Done! - -In principle, you now have the recipe in your hands to solve all your prediction/classification/regression (etc.) problems by the same Bayesian modeling method: - -""" - -# ╔═╡ c03229ef-3e0f-4612-909b-97f488a1e4c9 -md""" - REPEAT - 1. Model Specification - 2. Parameter Estimation - 3. Model Evaluation - UNTIL model performance is satisfactory - 4. Apply Model -""" - -# ╔═╡ 6a27951c-d294-11ef-2e1a-b5a4ce84aceb -md""" -Crucially, there is no need to invent clever machine learning algorithms, and there is no need to invent a clever prediction algorithm nor a need to invent a model performance criterion. Instead, you propose a model and, from there on, you let PT reason about everything that you care about. - -""" - -# ╔═╡ 6a27a28a-d294-11ef-1f33-41b444761429 -md""" -Your problems are only of computational nature. Perhaps the integral to compute the evidence may not be analytically tractable, etc. +**We're DONE!** Again, there is no need to invent a special prediction algorithm. Probability theory takes care of all that. Your problems are only of computational nature. Perhaps the integral to compute the evidence may not be analytically tractable, how to carry out the marginalization over ``\theta``, etc. """ -# ╔═╡ 55dec435-aa78-41ba-aad5-9d79ce292f42 -keyconcept("", "Bayesian machine learning is a subfield of machine learning that commits entirely to probability theory as the framework for all information-processing tasks. This is well justified, because probability theory is the optimal calculus for representing and manipulating states of knowledge.") +# ╔═╡ ce75e785-868f-4361-93f8-c582ac1b891b +keyconcept(" ", + md""" + + Bayesian Machine learning is EASY, apart from computational details :) + + """ +) # ╔═╡ 6a27b114-d294-11ef-099d-1d55968934a6 md""" @@ -506,13 +299,8 @@ I'd like to convince you that $(HTML("Bayesia # ╔═╡ 6a27beca-d294-11ef-1895-d57b11b827c1 md""" -Consider a model ``p(x,\theta|m)`` and a data set ``D = \{x_1,x_2, \ldots,x_N\}``. - -""" +Consider a model ``p(x,\theta|m)`` and a data set ``D = \{x_1,x_2, \ldots,x_N\}``. Given the data set ``D``, the log-evidence for model ``m`` decomposes as -# ╔═╡ cc8af69e-6d00-4327-aaa2-0b1023052b8a -md""" -Given the data set ``D``, the log-evidence for model ``m`` decomposes as """ # ╔═╡ c454be00-05e7-42f6-a243-bf559ed6eff7 @@ -540,19 +328,15 @@ md""" # ╔═╡ 6a27efc6-d294-11ef-2dc2-3b2ef95e72f5 md""" +#### accuracy The "accuracy" term (also known as data fit) measures how well the model predicts the data set ``D``. We want this term to be high because good models should predict the data ``D`` well. Indeed, higher accuracy leads to higher model evidence. To achieve high accuracy, applying Bayes' rule will shift the posterior ``p(\theta|D)`` away from the prior towards the likelihood function ``p(D|\theta)``. """ # ╔═╡ 6a280132-d294-11ef-10ac-f3890cb3f78b md""" -The second term ("complexity", also known as "information gain") is technically a [Kullback-Leibler divergence](https://en.wikipedia.org/wiki/Kullback%E2%80%93Leibler_divergence) (KLD) between the posterior and prior distributions, see [OPTIONAL SLIDE](#KLD) below. The KLD is an information-theoretic quantity that can be interpreted as a "distance" measure between two distributions. In other words, the complexity term measures how much the beliefs about ``\theta`` changed, due to learning from the data ``D``. Generally, we like the complexity term to be low, because moving away means forgetting previously acquired information represented by the prior. Indeed, lower complexity leads to higher model evidence. - -""" - -# ╔═╡ 6a2814b0-d294-11ef-3a76-9b93c1fcd4d5 -md""" -Models with high evidence ``p(D|m)`` prefer both high accuracy and low complexity. Therefore, models with high evidence tend to predict the training data ``D`` well (high accuracy), yet also try to preserve the information encoded by the prior (low complexity). These types of models are said to *generalize* well, since they can be applied to different data sets without specific adaptations for each data set. +#### complexity +The second term ("complexity", also known as "information gain") is technically a [Kullback-Leibler divergence](https://en.wikipedia.org/wiki/Kullback%E2%80%93Leibler_divergence) (KLD) between the posterior and prior distributions, see [this mini lecture](https://bmlip.github.io/course/minis/KL%20Divergence.html). The KLD is an information-theoretic quantity that can be interpreted as a "distance" measure between two distributions. In other words, the complexity term measures how much the beliefs about ``\theta`` changed, due to learning from the data ``D``. Generally, we like the complexity term to be low, because moving away means forgetting previously acquired information represented by the prior. Indeed, lower complexity leads to higher model evidence. """ @@ -570,21 +354,8 @@ keyconcept(" ", """ ) -# ╔═╡ 6a2879e6-d294-11ef-37db-df7babe24d25 -md""" -Put provocatively, this highlights that the common machine learning practice of splitting a dataset into training, validation, and test sets is, in essence, an ad hoc workaround, a substitute for formulating the learning task properly as a Bayesian inference problem. - -""" - -# ╔═╡ 6a2889ae-d294-11ef-2439-e1a541a5ccd7 -md""" -## Bayesian Modeling and the Scientific Method Revisited - -The Bayesian modeling approach provides a unified framework for the Scientific Inquiry method. We can now add equations to the design loop. (Trial design to be discussed in [Intelligent Agent lesson](https://bmlip.github.io/course/lectures/Intelligent%20Agents%20and%20Active%20Inference.html).) - -![](https://github.com/bmlip/course/blob/v2/assets/figures/scientific-inquiry-loop-w-BML-eqs.png?raw=true) - -""" +# ╔═╡ f5d8d021-3157-464f-93a2-b3054779e55f +NotebookCard("https://bmlip.github.io/course/minis/KL%20Divergence.html") # ╔═╡ c050f468-7eec-403f-9304-552bd0d9b222 html""" @@ -602,13 +373,24 @@ challenge_solution("Predicting a Coin Toss",header_level=1) # ╔═╡ 6a2898ea-d294-11ef-39ec-31e4bac1e048 md""" -At the beginning of this lesson, we posed the following challenge: +## Data Generation -We observe a the following sequence of heads (outcome = ``1``) and tails (outcome = ``0``) when tossing the same coin repeatedly +Let's generate a sequence of ``N`` coin tosses ``D=\{x_1,\ldots,x_N\}``, where each throw is drawn from a Bernoulli distribution +```math +p(x_n|\mu=0.4)=0.4^{x_n} \cdot 0.6^{1-x_n} \,, +``` +and where ``x_n`` denotes outcomes by ```math -D=\{1011001\}\,. +x_n = \begin{cases} 1 & \text{if heads comes up} \\ + 0 & \text{otherwise (tails)} \end{cases} ``` +So, this coin is biased! +""" + + +# ╔═╡ e5d9be85-3dfc-4145-a6d1-c872ccae4bc8 +md""" What is the probability that heads comes up next? We solve this in the next slides ... @@ -618,19 +400,7 @@ What is the probability that heads comes up next? We solve this in the next slid md""" ## 1. Model Specification for Coin Toss -We observe a sequence of ``N`` coin tosses ``D=\{x_1,\ldots,x_N\}`` with ``n`` heads. - -""" - -# ╔═╡ 6a28b44c-d294-11ef-15da-81be8753d311 -md""" -Let us denote outcomes by -```math -x_k = \begin{cases} 1 & \text{if heads comes up} \\ - 0 & \text{otherwise (tails)} \end{cases} - -``` """ @@ -686,7 +456,7 @@ so we get a closed-form posterior. # ╔═╡ 6a28f466-d294-11ef-3af9-e34de9736c71 md""" -``\alpha`` and ``\beta`` are called **hyperparameters**, since they parameterize the distribution for another parameter (``\mu``). E.g., ``\alpha=\beta=1`` leads to a uniform prior for ``\mu``. We use Julia below to visualize some priors ``\mathrm{Beta}(\mu|\alpha,\beta)`` for different values of ``\alpha, \beta``. +``\alpha`` and ``\beta`` are called **hyperparameters**, since they parameterize the distribution for another parameter (``\mu``). E.g., ``\alpha=\beta=1`` leads to a uniform prior for ``\mu``. Below, we visualize some priors ``\mathrm{Beta}(\mu|\alpha,\beta)`` for different values of ``\alpha, \beta``. """ @@ -721,38 +491,35 @@ let ) end -# ╔═╡ ee3da94c-5e87-4a0b-8373-c01e339d28aa +# ╔═╡ 6aa75f5b-1701-437e-ba1c-6583698416d3 md""" -##### Pick Your Own Parameters - -α = $(@bind beta_pdf_a NumberField(.1:.1:1000; default=6.0)), -β = $(@bind beta_pdf_b NumberField(.1:.1:1000; default=2.0)) -""" +Concretely, let's compare two models ``m_1``, and ``m_2`` with the same likelihood function, but with different priors: -# ╔═╡ 261620b0-9580-4d9e-b7de-d7972ea549cd -let - α = isnan(beta_pdf_a) ? 0.1 : beta_pdf_a - β = isnan(beta_pdf_b) ? 0.1 : beta_pdf_b - x = 0:0.01:1 +```math +\begin{aligned} +p(\mu|m_1) &= \mathrm{Beta}(\mu|\alpha=100,\beta=500) \\ +p(\mu|m_2) &= \mathrm{Beta}(\mu|\alpha=8,\beta=13). \\ +\end{aligned} +``` - y = pdf.(Beta(α, β), x) - plot(x, y; - # label="α=$α, β=$β", - ylim=(0, clamp(maximum(y)*1.05, 4, 4)), - label=nothing, - xlabel="μ", - ylabel="Density", - size=(600,250) - ) - -end +""" -# ╔═╡ 6a294790-d294-11ef-270b-5b2152431426 +# ╔═╡ e55126ef-e956-464d-8ae0-32b077649f21 md""" -#### Choosing a Prior - -Before observing any data, you can express your state-of-knowledge about the coin by choosing values for ``\alpha`` and ``\beta`` that reflect your beliefs. Stronger yet, you *must* choose values for ``\alpha`` and ``\beta``, because the Bayesian framework does not allow you to walk away from your responsibility to explicitly state your beliefs before the experiment. +> #### We can already guess which one is better! +> +> You can verify that model ``m_2`` has the best prior, since +> +> ```math +> \begin{align*} +> p(x_n=1|m_1) &= \left.\frac{\alpha}{\alpha+\beta}\right|_{m_1} = 100/600 \approx 0.17 \\ +> p(x_n=1|m_2) &= \left.\frac{\alpha}{\alpha+\beta}\right|_{m_2} = 8/21 \approx 0.38 \,, +> \end{align*} +> ``` +> +> (but you are not supposed to know that the real coin has a probability ``0.4`` for heads.) +> """ @@ -775,7 +542,7 @@ where ``B(\alpha,\beta) \triangleq \frac{\Gamma(\alpha) \Gamma(\beta)}{\Gamma(\a """ # ╔═╡ b426df32-5629-4773-b862-101cfbd82d42 -details("Prove this, and click for solution", +details("Click for proof", md""" ```math \begin{flalign*} @@ -799,11 +566,17 @@ p(\mu|D) = \mathrm{Beta}(\mu|\,n+\alpha, N-n+\beta) """ +# ╔═╡ 6a2af90a-d294-11ef-07bd-018326577791 +md""" + +For each model, we plot the parameter **posteriors** ``p(\mu|D_n,m_\bullet)`` computed after ``n`` iterations. +""" + # ╔═╡ 6a29d548-d294-11ef-1361-ad2230cad02b md""" ## 3. Model Evaluation for Coin Toss -It follow from the above calculation that the evidence for model ``m`` can be analytically expressed as +It follows from the above calculation that the evidence for model ``m`` can be analytically expressed as ```math \begin{align} @@ -812,8 +585,13 @@ p(D|m) &= \frac{B(n+\alpha,N-n+\beta)}{B(\alpha,\beta)} \\ \end{align} ``` -The model evidence is a scalar. In the absence of an alternative model, the value of the model evidence is not very useful. However, you may want to compare the model evidence of this model to the evidence for another model on the same data set. +The model evidence is a scalar. + +Let's check this by plotting over time the relative Bayesian evidences for each model: +```math +\frac{p(D_n|m_i)}{\sum_{i=1}^2 p(D_n|m_i)} +``` """ # ╔═╡ 6a29e25e-d294-11ef-15ce-5bf3d8cdb64c @@ -846,6 +624,9 @@ The above integral computes the mean of a beta distribution, which is given by ` """ +# ╔═╡ 5483148f-b385-4afa-ad85-70efe08ba299 +TODO("below should incorporate the slider again and execute the predictions for both models.") + # ╔═╡ 6a2a0f18-d294-11ef-02c2-ef117377ca66 md""" Finally, we're ready to solve our challenge: for ``D=\{1011001\}`` and uniform prior (``\alpha=\beta=1``), we get @@ -856,758 +637,249 @@ Finally, we're ready to solve our challenge: for ``D=\{1011001\}`` and uniform p In other words, given the model assumptions (the Bernoulli data-generating distribution and Beta prior as specified above), and the observations ``D=\{1011001\}``, the probability for observing heads (outcome=``1``) on the next toss is ``\frac{5}{9}``. -""" - -# ╔═╡ 6a2a1daa-d294-11ef-2a67-9f2ac60a14c5 -md""" -Be aware that there is no such thing as an "objective" or "correct" prediction. Every prediction is conditional on the selected model and the used data set. - -""" - -# ╔═╡ 6a2a2af2-d294-11ef-0072-bdc3c6f95bb3 -md""" -## What did we learn from the data? - -What did we learn from the data? Before seeing any data, we think that the probability for throwing heads is - -```math -\left. p(x_\bullet=1|D) \right|_{n=N=0} = \left.\frac{n+\alpha}{N+\alpha+\beta}\right|_{n=N=0} = \frac{\alpha}{\alpha + \beta}\,. -``` - -""" - -# ╔═╡ 6a2a389e-d294-11ef-1b8c-b55de794b65c -md""" -Hence, ``\alpha`` and ``\beta`` can be interpreted as prior pseudo-counts for heads and tails, respectively. - -""" - -# ╔═╡ 6a2a465e-d294-11ef-2aa0-43c954a6439e -md""" -If we were to assume zero pseudo-counts, i.e. ``\alpha=\beta \rightarrow 0``, then our prediction for throwing heads after ``N`` coin tosses is completely based on the data, given by - -```math -\left. p(x_\bullet=1|D) \right|_{\alpha=\beta \rightarrow 0} = \left.\frac{n+\alpha}{N+\alpha+\beta}\right|_{\alpha=\beta \rightarrow 0} = \frac{n}{N}\,. -``` - -""" - -# ╔═╡ 48fd2dff-796d-48bc-b5a8-bee270d119fd -md""" -Note the following decomposition -""" - -# ╔═╡ e3f9e571-2248-403c-8ab8-f6b99597f595 -md""" -```math -\begin{flalign*} - p(x_\bullet=1|\,D) &= \frac{n+\alpha}{N+\alpha+\beta} \\ - &= \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} + \underbrace{\underbrace{\frac{N}{N+\alpha+\beta}}_{\text{gain}}\cdot \underbrace{\biggl( \underbrace{\frac{n}{N}}_{\substack{\text{data-based}\\\text{prediction}}} - \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} \biggr)}_{\text{prediction error}}}_{\text{correction}} -\end{flalign*} -``` - -""" - -# ╔═╡ 90f691ad-046c-4595-99b0-19a1d6cb599e -details("Prove this yourself, and click for solution", -md""" -```math -\begin{align*} - p(x_\bullet=1|\,D) &= \frac{n+\alpha}{N+\alpha+\beta} \\ - &= \frac{\alpha}{N+\alpha+\beta} + \frac{n}{N+\alpha+\beta} \\ - &= \frac{\alpha}{N+\alpha+\beta}\cdot \frac{\alpha+\beta}{\alpha+\beta} + \frac{n}{N+\alpha+\beta}\cdot \frac{N}{N} \\ - &= \frac{\alpha}{\alpha+\beta}\cdot \frac{\alpha+\beta}{N+\alpha+\beta} + \frac{N}{N+\alpha+\beta}\cdot \frac{n}{N} \\ - &= \frac{\alpha}{\alpha+\beta}\cdot \biggl(1-\frac{N}{N+\alpha+\beta} \biggr) + \frac{N}{N+\alpha+\beta}\cdot \frac{n}{N} \\ - &= \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} + \underbrace{\underbrace{\frac{N}{N+\alpha+\beta}}_{\text{gain}}\cdot \underbrace{\biggl( \underbrace{\frac{n}{N}}_{\substack{\text{data-based}\\\text{prediction}}} - \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} \biggr)}_{\text{prediction error}}}_{\text{correction}} -\end{align*} -``` - """) - -# ╔═╡ 6a2a9faa-d294-11ef-1284-cfccb1da444e -md""" -Let's interpret this decomposition of the posterior prediction. Before the data ``D`` was observed, our model generated a *prior prediction* ``p(x_\bullet=1) = \frac{\alpha}{\alpha+\beta}``. Next, the degree to which the actually observed data matches this prediction is represented by the *prediction error* ``\frac{n}{N} - \frac{\alpha}{\alpha-\beta}``. The prior prediction is then updated to a *posterior prediction* ``p(x_\bullet=1|D)`` by adding a fraction of the prediction error to the prior prediction. Hence, the data plays the role of "correcting" the prior prediction. - -Note that, since ``0\leq \underbrace{\frac{N}{N+\alpha+\beta}}_{\text{gain}} \lt 1``, the Bayesian prediction lies between (fuses) the prior and data-based predictions. - -""" - -# ╔═╡ 6a2aad42-d294-11ef-3129-3be5be8c82d6 -md""" -For large ``N``, the gain goes to ``1`` and ``\left. p(x_\bullet=1|D)\right|_{N\rightarrow \infty} \rightarrow \frac{n}{N}`` goes to the data-based prediction (the observed relative frequency). - -""" - -# ╔═╡ 6a2abb16-d294-11ef-0243-d376e8a39bb0 -code_example("Bayesian Evolution for the Coin Toss"; big=true) - -# ╔═╡ 54cb380d-8864-4158-ba68-55027ae68971 -md""" -Let's code an example for a sequence of coin tosses, where we assume that the true coin generates data ``x_n \in \{0,1\}`` by a Bernoulli distribution: - -```math -p(x_n|\mu=0.4)=0.4^{x_n} \cdot 0.6^{1-x_n} -``` - -So, this coin is biased! - -""" - -# ╔═╡ 9da43d0f-e605-41b7-9bc6-db5be95bc87f -secret_distribution = Bernoulli(0.4); - -# ╔═╡ b791e819-f5a0-4c44-983b-07d8497516fb -@mdx """ - -```math -D=\\{$(Int.(rand(MersenneTwister(234), secret_distribution, intro_N)))\\}\\,. -``` - -""" - -# ╔═╡ 280e0819-674e-4a58-854d-5a66e0777074 -md""" -### Two Models -To predict the outcomes of future coin tosses, we'll compare **two models**: ``m_1`` and ``m_2``. Both models have the same data-generating distribution (also Bernoulli): - -```math -p(x_n|\mu,m_k) = \mu^{x_n} (1-\mu)^{1-x_n} \quad \text{for }k=1,2 \,, -``` - - -but they have different priors: - -```math -\begin{aligned} -p(\mu|m_1) &= \mathrm{Beta}(\mu|\alpha=100,\beta=500) \\ -p(\mu|m_2) &= \mathrm{Beta}(\mu|\alpha=8,\beta=13). \\ -\end{aligned} -``` - -""" - -# ╔═╡ e47b6eb6-2bb3-4c2d-bda6-f1535f2f94c4 -priors = [ - Beta(100., 500.), - Beta(8., 13.) -]; - -# ╔═╡ e55126ef-e956-464d-8ae0-32b077649f21 -md""" -> #### We can already guess which one is better! -> -> You can verify that model ``m_2`` has the best prior, since -> -> ```math -> \begin{align*} -> p(x_n=1|m_1) &= \left.\frac{\alpha}{\alpha+\beta}\right|_{m_1} = 100/600 \approx 0.17 \\ -> p(x_n=1|m_2) &= \left.\frac{\alpha}{\alpha+\beta}\right|_{m_2} = 8/21 \approx 0.38 \,, -> \end{align*} -> ``` -> -> (but you are not supposed to know that the real coin has a probability ``0.4`` for heads.) -> - -""" - -# ╔═╡ f67136ff-f33c-436e-823b-9c530d257ab0 -md""" -### Simulation - -Let's run ``500`` tosses: -""" - -# ╔═╡ d1d2bb84-7083-435a-9c19-4c02074143e3 -n_tosses = 500; - -# ╔═╡ 9c751f8e-f7ed-464f-b63c-41e318bbff2d -samples = rand(secret_distribution, n_tosses) - -# ╔═╡ 6922b899-499e-4f73-a6be-ca427fdc14ea -md""" -### Bayesian Machine Learning - -Now, let's update our posteriors iteratively. For every toss, we do a **Bayesian update step** to compute the new posterior for the distribution of ``\mu``. This posterior then becomes the **prior** for the next step. -""" - -# ╔═╡ 6a2af90a-d294-11ef-07bd-018326577791 -md""" - -For each model, we plot the parameter **posteriors** ``p(\mu|D_n,m_\bullet)`` computed after ``n`` iterations. -""" - -# ╔═╡ d484c41d-9834-4528-bf47-93ab4e35ebaa -md""" -Select iteration: $(@bind toss_index_1 Slider(0:n_tosses; show_value=true)) -""" - -# ╔═╡ 5ea6cefa-621a-4afc-bf9e-02d42b1d53f8 - - -# ╔═╡ 6a2b2d44-d294-11ef-33ba-15db357708b1 -md""" - -#### What happens with model 1? - -Note that both posteriors move toward the "correct" value (``\mu=0.4``). However, the posterior for ``m_1`` (blue) moves much slower because we assumed far more pseudo-observations for ``m_1`` than for ``m_2``. - - -""" - -# ╔═╡ e9b32823-efd2-4a27-b529-4f49752c00bb -keyconcept( - "", - "As additional observations are acquired, the influence of the prior progressively diminishes, while the likelihood increasingly dominates the posterior." -) - -# ╔═╡ 56b9aba3-6ead-498c-8670-ad93a1953b2a - - -# ╔═╡ c28b7130-f7fb-41ee-852e-9964b091d7fb -md""" -### Implementation: Iterative Bayesian Updating - -""" - -# ╔═╡ e99e7650-bb72-4576-8f2a-c3994533b644 -function handle_coin_toss(prior::Beta, observation::Bool) - posterior = Beta(prior.α + observation, prior.β + (1 - observation)) - return posterior -end; - -# ╔═╡ 51829800-1781-49ae-8ee7-ac15c0bfcb88 -# computes log10 of Gamma function -function log10gamma(num::Real)::Real - num = convert(BigInt, num) - return log10(gamma(num)) -end - -# ╔═╡ 7a624d2f-812a-47a0-a609-9fe299de94f5 -function log_evidence_prior(prior::Beta, N::Int64, n::Int64)::Real - log10gamma(prior.α + prior.β) - - log10gamma(prior.α) - - log10gamma(prior.β) + - log10gamma(n+prior.α) + - log10gamma((N-n)+prior.β) - - log10gamma(N+prior.α+prior.β) -end - -# ╔═╡ 3a903a4d-1fb0-4566-8151-9c86dfc40ceb -begin - # save a sequence of posterior distributions for every prior, starting with the prior itself - prior_distributions = [d for d in priors] - posterior_distributions = [[d] for d in priors] - log_evidences = [[] for _ in priors] - - - # for every sample we want to update our posterior - for (N, sample) in enumerate(samples) - # at every sample we want to update all distributions - for (i, prior) in enumerate(prior_distributions) - - # do bayesian updating - posterior = handle_coin_toss(prior, sample) - - # add posterior to vector of posterior distributions - push!(posterior_distributions[i], posterior) - - # compute log evidence and add to vector - log_evidence = log_evidence_prior(posterior_distributions[i][N], N, sum(samples[1:N])) - push!(log_evidences[i], log_evidence) - - # the prior for the next sample is the posterior from the current sample - prior_distributions[i] = posterior - end - end -end; - -# ╔═╡ 6a2b1106-d294-11ef-0d64-dbc26ba3eb44 -# Animate posterior distributions over time in a gif - -let i = toss_index_1 - p = plot(title="n = $i$(i == 0 ? " (priors)" : "")") - for (j,post) in enumerate(posterior_distributions) - plot!(post[i+1], xlims = (0, 1), fill=(0, .2,), label="Posterior model $j", linewidth=2, ylims=(0,28), xlabel="μ", legend=:topright) - end - vline!([mean(secret_distribution)]; style=:dash, color="purple", label="True parameter") -end - -# ╔═╡ f956e217-3dce-446a-8660-25f2c9cb05e2 -md""" -We now have a sequence of **posterior distributions** and **log evidences** for each model. Notice that the first "posterior" in each sequence is the prior, and the last posterior is the posterior that takes all data into account. This is our updated model, using Bayesian reasoning. - -_Click on the vectors below to see their values._ -""" - -# ╔═╡ 2c90eee1-b5d9-434d-bccc-64de8b458a48 -posterior_distributions - -# ╔═╡ 9d82af33-8e91-48e9-8c34-fa6ea31492c2 -log_evidences - -# ╔═╡ 69eaf045-e766-4b7f-a9e8-8eac674ca2ae - - -# ╔═╡ 6a2b3ba4-d294-11ef-3c28-176be260cb15 -md""" - -### Evidence Visualised - -We have an intuition that ``m_2`` is superior over ``m_1``. Let's check this by plotting over time the relative Bayesian evidences for each model: - -```math -\frac{p(D_n|m_i)}{\sum_{i=1}^2 p(D_n|m_i)} -``` - -""" - -# ╔═╡ ebcfcd1b-7fc8-42b7-a35e-4530f798cfdf -md""" -Select iteration: $(@bind toss_index_2 Slider(1:n_tosses; show_value=true)) -""" - -# ╔═╡ 188b5bea-6765-4dcf-9369-3b1fdbe94494 -let i = toss_index_2 - evidences = map(model -> exp.(model), log_evidences) - - plot(title=string(L"\frac{p_i(\mathbf{x}_{1:n})}{\sum_i p_i(\mathbf{x}_{1:n})}"," n = ", i), ylims=(0, 1), legend=:topleft) - total = sum(e[i] for e in evidences) - bar!([(e[i] / total) for e in evidences], group=["Model $i" for i in eachindex(priors)]) -end - -# ╔═╡ ee4006aa-e54b-4501-93ee-60b34bdf5c7b -exercise_statement("Convergence to 0"; header_level=4) - -# ╔═╡ 6a2b9676-d294-11ef-241a-89ff7aa676f9 -md""" -Over time, the relative evidence of model ``m_1`` converges to ``0``. Can you explain this behavior? - -""" - -# ╔═╡ 9c5d7c89-f65c-4f52-9e49-14692bed2452 -md""" -# Maximum Likelihood Estimation -""" - -# ╔═╡ 6a2bb18a-d294-11ef-23bb-99082caf6e01 -md""" -## From Posterior to Point-Estimate - -In the example above, Bayesian parameter estimation and prediction were tractable in closed form. This is often not the case. In that case, we will need to approximate some of the computations. - -""" - -# ╔═╡ 6a2bd3ac-d294-11ef-0543-6fe202ca35b6 -md""" -Recall Bayesian prediction - -```math -p(x|D) = \int p(x|\theta)p(\theta|D)\,\mathrm{d}{\theta} -``` - -""" - -# ╔═╡ 6a2bf332-d294-11ef-1ff1-cdbfb7732cf1 -md""" -If we approximate the posterior by a delta function, i.e., ``p(\theta|D) = \delta(\theta-\hat\theta)`` for one "best" value ``\hat\theta``, then the predictive distribution collapses to - -```math -p(x|D)= \int p(x|\theta)\,\delta(\theta-\hat\theta)\,\mathrm{d}{\theta} = p(x|\hat\theta) -``` - -""" - -# ╔═╡ 6a2c008e-d294-11ef-2f07-11cdfb2bddca -md""" -This is just the data-generating distribution ``p(x|\theta)`` evaluated at ``\theta=\hat\theta``, which is easy to evaluate. - -""" - -# ╔═╡ 6a2c11e6-d294-11ef-173b-23fc6dbfefca -md""" -The next question is how to get the parameter estimate ``\hat{\theta}``? (See next slide). - -""" - -# ╔═╡ 6a2c229e-d294-11ef-2f24-ebe43cbfbfa4 -md""" -## Some Well-known Point-Estimates - -- **Bayes estimate** (the mean of the posterior) - -```math -\hat \theta_{\text{Bayes}} = \int \theta \, p\left( \theta |D \right) -\,\mathrm{d}{\theta} -``` - -""" - -# ╔═╡ 6a2c3036-d294-11ef-23cb-c3b36c475e8f -md""" -- **Maximum A Posteriori** (MAP) estimate - -```math -\hat \theta_{\text{map}}= \arg\max _{\theta} p\left( \theta |D \right) = -\arg \max_{\theta} p\left(D |\theta \right) \, p\left(\theta \right) -``` - -""" - -# ╔═╡ 6a2c4058-d294-11ef-2312-d9c672d49701 -md""" -- **Maximum Likelihood** (ML) estimate - -```math -\hat \theta_{ml} = \arg \max_{\theta} p\left(D |\theta\right) -``` - -Note that Maximum Likelihood (ML) is MAP with a uniform prior. MAP is sometimes called a 'penalized' ML procedure: - -```math -\hat \theta_{map} = \arg \max _\theta \{ \underbrace{\log -p\left( D|\theta \right)}_{\text{log-likelihood}} + \underbrace{\log -p\left( \theta \right)}_{\text{penalty}} \} -``` - -ML is the most common approximation to the full Bayesian posterior. - -""" - -# ╔═╡ 6a2c505c-d294-11ef-1c92-c1b0e9d50da5 -md""" -## Bayesian vs Maximum Likelihood Learning - -Consider the task: predict a future observation ``x`` from an observed data set ``D``. Let us compare full Bayesian modeling with the maximum likelihood approach. - -""" - -# ╔═╡ 7c8b1add-085a-41ba-9d6c-b26d3eef22e4 -md""" - -| | **Bayesian** | **Maximum Likelihood** | -|:----|:---------|:-----| -| 1. **Model Specification** | Choose a model ``m`` with data-generating distribution ``p(x\|\theta, m)`` and parameter prior ``p(\theta\|m)``. | Choose a model ``m`` with same data generating distribution ``p(x\|\theta, m)``. No need for priors. | -| 2. **Learning** | Use Bayes rule to find the parameter posterior: $(HTML("
"))``p(\theta\|D) \propto p(D\|\theta) p(\theta)`` | By Maximum Likelihood (ML) optimization: $(HTML("
")) ``\hat \theta = \arg \max_{\theta} p(D\|\theta)`` | -| 3. **Prediction** | ``p(x\|D) = \int p(x\|\theta) p(\theta\|D) \,\mathrm{d}\theta`` | ``p(x\|D) = p(x\|\hat\theta)`` | - - -""" - -# ╔═╡ 6a2c5e08-d294-11ef-213d-97bcfa16eb5a -md""" -## Report Card on Maximum Likelihood Estimation - - - -""" - -# ╔═╡ 6a2c7230-d294-11ef-05a2-3ff2f65d10e0 -md""" -(good!). ML works rather well if we have a lot of data because the influence of the prior diminishes with more data. - -""" - -# ╔═╡ 6a2c7f5a-d294-11ef-2e17-9108a39df280 -md""" -(good!). Computationally often do-able. Useful fact that makes the optimization easier (since ``\log`` is monotonously increasing): - -```math -\arg\max_\theta \log p(D|\theta) = \arg\max_\theta p(D|\theta) -``` - -""" - -# ╔═╡ 6a2c8f4a-d294-11ef-213c-dfa929a403bc -md""" -(bad). ML cannot be used for model comparison! In ML estimation, the Bayesian model evidence is undefined because no prior distribution is specified. Even if we attempt to simulate ML as a special case of Bayesian inference by using a uniform prior, the evidence still collapses: a uniform prior over the entire real line is not a proper probability distribution, since its integral does not evaluate to 1. Consequently, when performing ML estimation, Bayesian model evidence cannot be used to evaluate model performance: - -```math -\begin{align*} -\underbrace{p(D|m)}_{\substack{\text{Bayesian}\\ \text{evidence}}} &= \int p(D|\theta) \cdot p(\theta|m)\,\mathrm{d}\theta \\ - &= \lim_{(b-a)\rightarrow \infty} \int p(D|\theta)\cdot \underbrace{\text{Uniform}(\theta|a,b)}_{\text{"ML prior"}}\,\mathrm{d}\theta \\ - &= \lim_{(b-a)\rightarrow \infty} \frac{1}{b-a}\underbrace{\int_a^b p(D|\theta)\,\mathrm{d}\theta}_{<\infty} \\ - &= 0 -\end{align*} -``` - -In fact, this is a serious disadvantage because Bayesian evidence is a principled performance assessment criterion that follows from straightforward PT. In practice, when estimating parameters by maximum likelihood, we often evaluate model performance by an *ad hoc* performance measure such as mean-squared-error on a testing data set. - -""" - -# ╔═╡ 6a2ca496-d294-11ef-0043-1f350b36773e -keyconcept(" ", - md""" - Maximum likelihood estimation is, at best, an approximation to Bayesian learning. Still, it is a very popular method, and with good reason: when plenty of data are available, it often provides a practical and effective solution. - """ -) - - -# ╔═╡ 47842de0-d17e-460e-b3b7-b2e642569e25 -md""" -# Summary -""" - -# ╔═╡ b273c8bc-3819-4f63-801a-acf0ee78ef1d -keyconceptsummary() - -# ╔═╡ 4bfd141f-fe2f-46a0-aa35-872cab45ea00 -exercises(header_level=1) - -# ╔═╡ f2969d91-4a5b-4665-9fa5-521db750302f -md""" - -##### Bayes estimate (**) - -(##) The Bayes estimate is a summary of a posterior distribution by a delta distribution on its mean, i.e., - -```math -\hat \theta_{bayes} = \int \theta \, p\left( \theta |D \right) -\,\mathrm{d}{\theta} -``` - -Prove that the Bayes estimate minimizes the mean-squared error, i.e., Prove that - -```math -\hat \theta_{bayes} = \arg\min_{\hat \theta} \int_\theta (\hat \theta -\theta)^2 p \left( \theta |D \right) \,\mathrm{d}{\theta} -``` -""" - -# ╔═╡ 7dd9a456-9dca-47c8-98c5-51f87f28e6a4 -hide_solution( -md""" -To minimize the expected mean-squared error we will look for ``\hat{\theta}`` that makes the gradient of the integral with respect to ``\hat{\theta}`` vanish. - -```math -\begin{align*} - \nabla_{\hat{\theta}} \int_\theta (\hat \theta -\theta)^2 p \left( \theta |D \right) \,\mathrm{d}{\theta} &= 0 \\ - \int_\theta \nabla_{\hat{\theta}} (\hat \theta -\theta)^2 p \left( \theta |D \right) \,\mathrm{d}{\theta} &= 0 \\ - \int_\theta 2(\hat \theta -\theta) p \left( \theta |D \right) \,\mathrm{d}{\theta} &= 0 \\ - \int_\theta \hat \theta p \left( \theta |D \right) \,\mathrm{d}{\theta} &= \int_\theta \theta p \left( \theta |D \right) \,\mathrm{d}{\theta} \\ - \hat \theta \underbrace{\int_\theta p \left( \theta |D \right) \,\mathrm{d}{\theta}}_{1} &= \int_\theta \theta p \left( \theta |D \right) \,\mathrm{d}{\theta} \\ - \Rightarrow \hat \theta &= \int_\theta \theta p \left( \theta |D \right) \,\mathrm{d}{\theta} -\end{align*} -``` -""" -) - -# ╔═╡ b2820dfd-b3ca-477b-8cb7-c430e0fe18dd -md""" - -##### Coin Toss MAP and ML (**) - -Consider the coin toss example with model -```math -\begin{align} -p(x_k|\mu) &= \mu^{x_k} (1-\mu)^{1-x_k} \\ -p(\mu) &= \mathrm{Beta}(\mu|\alpha,\beta) \,. -\end{align} -``` -and a given data set ``D=\{x_1, x_2,\ldots,x_N\}``. - -- (a) Derive the Maximum Likelihood estimate for ``\mu``. -- (b) Derive the MAP estimate for ``\mu``. -- (c) Do these two estimates ever coincide (if so, under what circumstances)? +""" +# ╔═╡ 6a2a1daa-d294-11ef-2a67-9f2ac60a14c5 +md""" +Be aware that there is no such thing as an "objective" or "correct" prediction. Every prediction is conditional on the selected model and the used data set. """ -# ╔═╡ 664d4183-edb6-4818-a44b-bf4c0a22a33c -hide_solution( +# ╔═╡ 6a2a2af2-d294-11ef-0072-bdc3c6f95bb3 md""" -- (a) The likelihood is given by ``p(D|\mu) = \mu^n\cdot (1-\mu)^{(N-n)}``. It follows that +## All Learning is Interpretable as Correcting Prediction Errors +What did we learn from the data? Before seeing any data, we think that the probability for throwing heads is ```math -\begin{align*} - \nabla \log p(D|\mu) &= 0 \\ - \nabla \left( n\log \mu + (N-n)\log(1-\mu)\right) &= 0\\ - \frac{n}{\mu} - \frac{N-n}{1-\mu} &= 0 \\ - \rightarrow \hat{\mu}_{\text{ML}} &= \frac{n}{N} - \end{align*} +\left. p(x_\bullet=1|D) \right|_{n=N=0} = \left.\frac{n+\alpha}{N+\alpha+\beta}\right|_{n=N=0} = \frac{\alpha}{\alpha + \beta}\,. ``` -- (b) We can write the posterior as as - +""" -```math -\begin{align*} - p(\mu|D) &\propto p(D|\mu)p(\mu) \\ - &\propto \mu^n (1-\mu)^{N-n} \mu^{\alpha-1} (1-\mu)^{\beta-1} \\ - &\propto \mathcal{B}(\mu|n+\alpha,N-n+\beta) - \end{align*} -``` +# ╔═╡ 6a2a389e-d294-11ef-1b8c-b55de794b65c +md""" +Hence, ``\alpha`` and ``\beta`` can be interpreted as prior pseudo-counts for heads and tails, respectively. -The MAP estimate for a beta distribution ``\mathcal{B}(a,b)`` is located at ``\frac{a - 1}{a+b-2}``, see [wikipedia](https://en.wikipedia.org/wiki/Beta_distribution). Hence, +""" +# ╔═╡ 6a2a465e-d294-11ef-2aa0-43c954a6439e +md""" +If we were to assume zero pseudo-counts, i.e. ``\alpha=\beta \rightarrow 0``, then our prediction for throwing heads after ``N`` coin tosses is completely based on the data, given by ```math -\begin{align*} -\hat{\mu}_{\text{MAP}} &= \frac{(n+\alpha)-1}{(n+\alpha) + (N-n+\beta) -2} \\ - &= \frac{n+\alpha-1}{N + \alpha +\beta -2} -\end{align*} +\left. p(x_\bullet=1|D) \right|_{\alpha=\beta \rightarrow 0} = \left.\frac{n+\alpha}{N+\alpha+\beta}\right|_{\alpha=\beta \rightarrow 0} = \frac{n}{N}\,. ``` -- (c) As ``N`` gets larger, the MAP estimate approaches the ML estimate. In the limit the MAP solution converges to the ML solution. - - """ -) -# ╔═╡ ecb036da-a0a2-4919-b1aa-bc33b6ba7e73 +# ╔═╡ 48fd2dff-796d-48bc-b5a8-bee270d119fd md""" +Note the following decomposition +""" -##### Model Comparison (**) - -A model ``m_1`` is described by a single parameter ``\theta``, with ``0 \leq \theta \leq 1``. The system can produce data ``x \in \{0,1\}``. The sampling distribution and prior are given by - +# ╔═╡ e3f9e571-2248-403c-8ab8-f6b99597f595 +md""" ```math -\begin{align*} -p(x|\theta,m_1) &= \theta^x (1-\theta)^{(1-x)} \\ -p(\theta|m_1) &= 6\theta(1-\theta) -\end{align*} +\begin{flalign*} + p(x_\bullet=1|\,D) &= \frac{n+\alpha}{N+\alpha+\beta} \\ + &= \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} + \underbrace{\underbrace{\frac{N}{N+\alpha+\beta}}_{\text{gain}}\cdot \underbrace{\biggl( \underbrace{\frac{n}{N}}_{\substack{\text{data-based}\\\text{prediction}}} - \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} \biggr)}_{\text{prediction error}}}_{\text{correction}} +\end{flalign*} ``` -- (a) Work out the probability ``p(x=1|m_1)``. - -- (b) Determine the posterior ``p(\theta|x=1,m_1)``. - -Now consider a second model ``m_2`` with the following sampling distribution and prior on ``0 \leq \theta \leq 1``: +""" +# ╔═╡ 90f691ad-046c-4595-99b0-19a1d6cb599e +details("Prove this yourself, and click for solution", +md""" ```math \begin{align*} -p(x|\theta,m_2) &= (1-\theta)^x \theta^{(1-x)} \\ -p(\theta|m_2) &= 2\theta + p(x_\bullet=1|\,D) &= \frac{n+\alpha}{N+\alpha+\beta} \\ + &= \frac{\alpha}{N+\alpha+\beta} + \frac{n}{N+\alpha+\beta} \\ + &= \frac{\alpha}{N+\alpha+\beta}\cdot \frac{\alpha+\beta}{\alpha+\beta} + \frac{n}{N+\alpha+\beta}\cdot \frac{N}{N} \\ + &= \frac{\alpha}{\alpha+\beta}\cdot \frac{\alpha+\beta}{N+\alpha+\beta} + \frac{N}{N+\alpha+\beta}\cdot \frac{n}{N} \\ + &= \frac{\alpha}{\alpha+\beta}\cdot \biggl(1-\frac{N}{N+\alpha+\beta} \biggr) + \frac{N}{N+\alpha+\beta}\cdot \frac{n}{N} \\ + &= \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} + \underbrace{\underbrace{\frac{N}{N+\alpha+\beta}}_{\text{gain}}\cdot \underbrace{\biggl( \underbrace{\frac{n}{N}}_{\substack{\text{data-based}\\\text{prediction}}} - \underbrace{\frac{\alpha}{\alpha+\beta}}_{\substack{\text{prior}\\\text{prediction}}} \biggr)}_{\text{prediction error}}}_{\text{correction}} \end{align*} ``` + """) -- (c) Determine the probability ``p(x=1|m_2)``. - -Now assume that the model priors are given by - -```math -\begin{align*} - p(m_1) &= 1/3 \\ - p(m_2) &= 2/3 - \end{align*} -``` +# ╔═╡ 6a2a9faa-d294-11ef-1284-cfccb1da444e +md""" +Let's interpret this decomposition of the posterior prediction. Before the data ``D`` was observed, our model generated a *prior prediction* ``p(x_\bullet=1) = \frac{\alpha}{\alpha+\beta}``. Next, the degree to which the actually observed data matches this prediction is represented by the *prediction error* ``\frac{n}{N} - \frac{\alpha}{\alpha-\beta}``. The prior prediction is then updated to a *posterior prediction* ``p(x_\bullet=1|D)`` by adding a fraction ``0\leq \frac{N}{N+\alpha+\beta} \lt 1`` of the prediction error to the prior prediction. Hence, the **prediction error** plays the role of "correcting" the prior prediction. -- (d) Compute the probability ``p(x=1)`` by "Bayesian model averaging", i.e., by weighing the predictions of both models appropriately. +""" +# ╔═╡ 496a4329-bee2-4265-8854-e6019b358d01 +md""" +# Example: Bayesian (Linear) Regression +""" -- (e) Compute the fraction of posterior model probabilities ``\frac{p(m_1|x=1)}{p(m_2|x=1)}``. +# ╔═╡ 679ef9d1-cc1c-4fc1-bf82-caa967c196c2 +md""" +# Example: Generative Classification +""" +# ╔═╡ 47842de0-d17e-460e-b3b7-b2e642569e25 +md""" +# Summary +""" -- (f) Which model do you prefer after observation ``x=1``? +# ╔═╡ b273c8bc-3819-4f63-801a-acf0ee78ef1d +keyconceptsummary() +# ╔═╡ 6a2cb25e-d294-11ef-1d88-1fc784b33df0 +md""" +# Optional Slides """ -# ╔═╡ de08c2a1-c5e3-4add-8b22-2c633247da48 -hide_solution( +# ╔═╡ ee3da94c-5e87-4a0b-8373-c01e339d28aa md""" -- (a) Work out the probability ``p(x=1|m_1)``. - -```math -\begin{align*} - p(x=1|m_1) &= \int_0^1 p(x=1|\theta,m_1) p(\theta|m_1) \mathrm{d}\theta \\ - &= \int \theta \cdot 6\theta (1-\theta) \mathrm{d}\theta \\ - &= 6 \cdot \left(\frac{1}{3}\theta^3 - \frac{1}{4}\theta^4\right) \bigg|_0^1 \\ - &= 6 \cdot (\frac{1}{3} - \frac{1}{4}) = \frac{1}{2} -\end{align*} -``` - -- (b) Determine the posterior ``p(\theta|x=1,m_1)``. +##### Pick Your Own Parameters for the Beta distribution -```math -\begin{align*} - p(\theta|x=1,m_1) &= \frac{p(x=1|\theta) p(\theta|m_1)}{p(x=1|m_1)} \\ - &= 2\cdot \theta \cdot 6\theta (1-\theta) \\ - &= \begin{cases} 12 \theta^2 (1-\theta) & \text{if }0 \leq \theta \leq 1 \\ - 0 & \text{otherwise} \end{cases} - \end{align*} -``` -Now consider a second model ``m_2`` with the following sampling distribution and prior on ``0 \leq \theta \leq 1``: +α = $(@bind beta_pdf_a NumberField(.1:.1:1000; default=6.0)), +β = $(@bind beta_pdf_b NumberField(.1:.1:1000; default=2.0)) +""" -```math -\begin{align*} -p(x|\theta,m_2) &= (1-\theta)^x \theta^{(1-x)} \\ -p(\theta|m_2) &= 2\theta -\end{align*} -``` +# ╔═╡ 261620b0-9580-4d9e-b7de-d7972ea549cd +let + α = isnan(beta_pdf_a) ? 0.1 : beta_pdf_a + β = isnan(beta_pdf_b) ? 0.1 : beta_pdf_b + x = 0:0.01:1 -- (c) Determine the probability ``p(x=1|m_2)``. + y = pdf.(Beta(α, β), x) + plot(x, y; + # label="α=$α, β=$β", + ylim=(0, clamp(maximum(y)*1.05, 4, 4)), + label=nothing, + xlabel="μ", + ylabel="Density", + size=(600,250) + ) + +end -```math -\begin{align*} - p(x=1|m_2) &= \int_0^1 p(x=1|\theta,m_2) p(\theta|m_2) \mathrm{d}\theta \\ - &= \int (1-\theta) \cdot 2\theta \mathrm{d}\theta \\ - &= 2 \cdot \left( \frac{1}{2}\theta^2 - \frac{1}{3}\theta^3 \right) \bigg|_0^1 \\ - &= 2 \cdot (\frac{1}{2} - \frac{1}{3}) = \frac{1}{3} - \end{align*} -``` +# ╔═╡ 1f92c406-6792-4af6-9132-35efd8223bc5 +md""" +# Code +""" -Now assume that the model priors are given by +# ╔═╡ 7a764a14-a5df-4f76-8836-f0a571fc3519 +wideq(x) = PlutoUI.ExperimentalLayout.Div([x]; style="min-width: max-content;") |> WideCell -```math -\begin{align*} - p(m_1) &= 1/3 \\ - p(m_2) &= 2/3 - \end{align*} -``` +# ╔═╡ c28b7130-f7fb-41ee-852e-9964b091d7fb +md""" +### Coin Toss Implementation -- (d) Compute the probability ``p(x=1)`` by "Bayesian model averaging", i.e., by weighing the predictions of both models appropriately. +""" -```math -\begin{align*} - p(x=1) &= \sum_{k=1}^2 p(x=1|m_k) p(m_k) \\ - &= \frac{1}{2} \cdot \frac{1}{3} + \frac{1}{3} \cdot \frac{2}{3} = \frac{7}{18} - \end{align*} -``` +# ╔═╡ 9da43d0f-e605-41b7-9bc6-db5be95bc87f +secret_distribution = Bernoulli(0.4); -- (e) Compute the fraction of posterior model probabilities ``\frac{p(m_1|x=1)}{p(m_2|x=1)}``. +# ╔═╡ b791e819-f5a0-4c44-983b-07d8497516fb +@mdx """ ```math -\frac{p(m_1|x=1)}{p(m_2|x=1)} = \frac{p(x=1|m_1) p(m_1)}{p(x=1|m_2) p(m_2)} = \frac{\frac{1}{2} \cdot \frac{1}{3}}{\frac{1}{3} \cdot \frac{2}{3}} =\frac{3}{4} +D=\\{$(Int.(rand(MersenneTwister(234), secret_distribution, intro_N)))\\}\\,. ``` -- (f) Which model do you prefer after observation ``x=1``? +""" -In principle, the observation ``x=1`` favors model ``m_2``, since ``p(m_2|x=1) = \frac{4}{3} \times p(m_1|x=1)``. However, note that ``\log_{10} \frac{3}{4} \approx -0.125``, so the extra evidence for ``m_2`` relative to ``m_1`` is very low. At this point, after 1 observation, we have no preference for a model yet. +# ╔═╡ e47b6eb6-2bb3-4c2d-bda6-f1535f2f94c4 +priors = [ + Beta(100., 500.), + Beta(8., 13.) +]; -""") +# ╔═╡ d1d2bb84-7083-435a-9c19-4c02074143e3 +n_tosses = 500; -# ╔═╡ 6a2cb25e-d294-11ef-1d88-1fc784b33df0 +# ╔═╡ d484c41d-9834-4528-bf47-93ab4e35ebaa md""" -# Optional Slides - +Select iteration: $(@bind toss_index_1 Slider(0:n_tosses; show_value=true)) """ -# ╔═╡ 1edae118-dcc7-4169-95cf-f36025f2c336 +# ╔═╡ ebcfcd1b-7fc8-42b7-a35e-4530f798cfdf md""" -## Working with Distributions in code - -Take a look at this mini lecture to see some simple examples of using distributions in Julia: +Select iteration: $(@bind toss_index_2 Slider(1:n_tosses; show_value=true)) """ -# ╔═╡ 275a9a69-3135-4cbd-8a35-b1abee4af83f -NotebookCard("https://bmlip.github.io/course/minis/Distributions%20in%20Julia.html") +# ╔═╡ 9c751f8e-f7ed-464f-b63c-41e318bbff2d +samples = rand(secret_distribution, n_tosses) -# ╔═╡ 6a2ccd16-d294-11ef-22ee-a5cff62ccd9c -md""" -## The Kullback-Leibler Divergence$(HTML("")) +# ╔═╡ e99e7650-bb72-4576-8f2a-c3994533b644 +function handle_coin_toss(prior::Beta, observation::Bool) + posterior = Beta(prior.α + observation, prior.β + (1 - observation)) + return posterior +end; -""" +# ╔═╡ 758a13d7-5c4b-498c-8c37-c3ba9b4daf48 +# computes log10 of Gamma function +function log10gamma(num::Real)::Real + num = convert(BigInt, num) + return log10(gamma(num)) +end -# ╔═╡ f5d8d021-3157-464f-93a2-b3054779e55f -NotebookCard("https://bmlip.github.io/course/minis/KL%20Divergence.html") +# ╔═╡ 7a624d2f-812a-47a0-a609-9fe299de94f5 +function log_evidence_prior(prior::Beta, N::Int64, n::Int64)::Real + log10gamma(prior.α + prior.β) - + log10gamma(prior.α) - + log10gamma(prior.β) + + log10gamma(n+prior.α) + + log10gamma((N-n)+prior.β) - + log10gamma(N+prior.α+prior.β) +end; -# ╔═╡ 1f92c406-6792-4af6-9132-35efd8223bc5 -md""" -# Code -""" +# ╔═╡ 3a903a4d-1fb0-4566-8151-9c86dfc40ceb +begin + # save a sequence of posterior distributions for every prior, starting with the prior itself + prior_distributions = [d for d in priors] + posterior_distributions = [[d] for d in priors] + log_evidences = [[] for _ in priors] -# ╔═╡ 7a764a14-a5df-4f76-8836-f0a571fc3519 -wideq(x) = PlutoUI.ExperimentalLayout.Div([x]; style="min-width: max-content;") |> WideCell + + # for every sample we want to update our posterior + for (N, sample) in enumerate(samples) + # at every sample we want to update all distributions + for (i, prior) in enumerate(prior_distributions) + + # do bayesian updating + posterior = handle_coin_toss(prior, sample) + + # add posterior to vector of posterior distributions + push!(posterior_distributions[i], posterior) + + # compute log evidence and add to vector + log_evidence = log_evidence_prior(posterior_distributions[i][N], N, sum(samples[1:N])) + push!(log_evidences[i], log_evidence) + + # the prior for the next sample is the posterior from the current sample + prior_distributions[i] = posterior + end + end +end; + +# ╔═╡ 6a2b1106-d294-11ef-0d64-dbc26ba3eb44 +# Animate posterior distributions over time in a gif + +let i = toss_index_1 + p = plot() + for (j,post) in enumerate(posterior_distributions) + plot!(post[i+1], xlims = (0, 1), fill=(0, .2,), label="Posterior model $j", linewidth=2, ylims=(0,28), xlabel="μ", legend=:topright) + end + vline!([mean(secret_distribution)]; style=:dash, color="purple", label="True parameter") +end + +# ╔═╡ 188b5bea-6765-4dcf-9369-3b1fdbe94494 +let i = toss_index_2 + evidences = map(model -> exp.(model), log_evidences) + + plot(ylims=(0, 1), legend=:topleft) + total = sum(e[i] for e in evidences) + bar!([(e[i] / total) for e in evidences], group=["Model $i" for i in eachindex(priors)]) +end # ╔═╡ 00000000-0000-0000-0000-000000000001 PLUTO_PROJECT_TOML_CONTENTS = """ @@ -3142,7 +2414,6 @@ version = "1.9.2+0" # ╔═╡ Cell order: # ╟─6a23b828-d294-11ef-371a-05d061144a43 # ╟─6be2e966-4048-44d0-a37e-95060e3fe30b -# ╟─6a23df9e-d294-11ef-3ddf-a51d4cea00fc # ╟─eca027f8-40c9-4e53-85b5-d08b8fe9dd97 # ╟─4f6a2d4f-bd89-4b0c-b544-397de2e34e72 # ╟─b791e819-f5a0-4c44-983b-07d8497516fb @@ -3151,6 +2422,8 @@ version = "1.9.2+0" # ╟─6a24c3e6-d294-11ef-3581-2755a9ba15ba # ╟─e2de9415-7bd8-4e95-abeb-53fc068ee950 # ╟─6a24c9f4-d294-11ef-20cc-172ea50da901 +# ╟─3cbe6b9b-4289-4724-ac3a-1b0313b3041b +# ╟─55dec435-aa78-41ba-aad5-9d79ce292f42 # ╟─6a24cee0-d294-11ef-35cb-71ab9ef935e5 # ╟─6a24d478-d294-11ef-2a75-9d03a5ba7ff8 # ╟─6a24fde8-d294-11ef-29bf-ad3e20a53c29 @@ -3160,75 +2433,55 @@ version = "1.9.2+0" # ╟─6a25307e-d294-11ef-0662-3db678b32e99 # ╟─6a25379a-d294-11ef-3e07-87819f6d75cb # ╟─6a254460-d294-11ef-1890-230b75b6b9ee -# ╟─6a2552ac-d294-11ef-08d6-179e068bc297 -# ╟─ce75e785-868f-4361-93f8-c582ac1b891b # ╟─6a2561c0-d294-11ef-124d-373846e3120c -# ╟─6a257020-d294-11ef-0490-e151934b2f42 -# ╟─6a257f34-d294-11ef-2928-fbb800e81124 # ╟─6a25a11e-d294-11ef-1c51-09482dad86f2 # ╟─6a25edfc-d294-11ef-3411-6f74c376461e # ╟─53de7edd-6c28-49a7-9f54-cf7b8ca42aeb # ╟─288fbee6-0783-4447-b5d0-f5c2b29b39c7 # ╟─74fa1925-0d9f-47f6-a6bd-b822948a4fbc -# ╟─6a261278-d294-11ef-25a0-5572de58ad06 -# ╟─6a26549a-d294-11ef-1f10-15c4d14ae41f -# ╟─6a262182-d294-11ef-23e9-ed45e1da9f46 # ╟─6a2672d6-d294-11ef-1886-3195c9c7cfa9 # ╟─6aa2399d-a949-40f9-8ee6-b0c2be1dc478 -# ╟─6a2664c6-d294-11ef-0a49-5192e17fb9ea -# ╟─6a26a31e-d294-11ef-2c2f-b349d0859a27 -# ╟─6a269568-d294-11ef-02e3-13402d296391 -# ╟─6a26b7bc-d294-11ef-03e7-2715b6f8dcc7 -# ╟─6a26f244-d294-11ef-0488-c1e4ec6e739d -# ╟─99db44c9-185c-4f39-ae5e-1a4cd751d980 -# ╟─d22f58ac-9f68-41cb-8e61-cf74d3692c44 # ╟─6a2707e6-d294-11ef-02ad-31bf84662c70 # ╟─6a271a56-d294-11ef-0046-add807cc0b4f # ╟─f6ee5570-9b92-42b6-baf3-3eed5352a060 # ╟─6a273ae0-d294-11ef-2c00-9b3eaed93f6d -# ╟─6a274948-d294-11ef-0563-1796b8883306 -# ╟─6a275a52-d294-11ef-1323-9d83972f611a -# ╟─6a27684e-d294-11ef-040e-c302cdad714a -# ╟─6a2777d0-d294-11ef-1ac3-add102c097d6 -# ╟─6a278784-d294-11ef-11ae-65bd398910d5 -# ╟─c03229ef-3e0f-4612-909b-97f488a1e4c9 -# ╟─6a27951c-d294-11ef-2e1a-b5a4ce84aceb -# ╟─6a27a28a-d294-11ef-1f33-41b444761429 -# ╟─55dec435-aa78-41ba-aad5-9d79ce292f42 +# ╟─ce75e785-868f-4361-93f8-c582ac1b891b # ╟─6a27b114-d294-11ef-099d-1d55968934a6 # ╟─6a27beca-d294-11ef-1895-d57b11b827c1 -# ╟─cc8af69e-6d00-4327-aaa2-0b1023052b8a # ╟─c454be00-05e7-42f6-a243-bf559ed6eff7 # ╟─6a9ad1c4-dfb2-4987-9ddc-da6131605083 # ╟─6a27efc6-d294-11ef-2dc2-3b2ef95e72f5 # ╟─6a280132-d294-11ef-10ac-f3890cb3f78b -# ╟─6a2814b0-d294-11ef-3a76-9b93c1fcd4d5 # ╟─6a282892-d294-11ef-2c12-4b1c7374617c # ╟─6a286b04-d294-11ef-1b34-8b7a85c0048c -# ╟─6a2879e6-d294-11ef-37db-df7babe24d25 -# ╟─6a2889ae-d294-11ef-2439-e1a541a5ccd7 +# ╟─f5d8d021-3157-464f-93a2-b3054779e55f # ╟─c050f468-7eec-403f-9304-552bd0d9b222 # ╟─1dbc69a3-b3ec-44de-af7c-944ebc01f523 # ╟─6a2898ea-d294-11ef-39ec-31e4bac1e048 +# ╟─e5d9be85-3dfc-4145-a6d1-c872ccae4bc8 # ╟─6a28a704-d294-11ef-1bf2-efbdb0cb4cbc -# ╟─6a28b44c-d294-11ef-15da-81be8753d311 # ╟─6a28c9b4-d294-11ef-222b-97bf0912efe7 # ╟─6a28d81e-d294-11ef-2a9f-d32daa5556ae # ╟─6a28e674-d294-11ef-391b-0d33fd609fb8 # ╟─6a28f466-d294-11ef-3af9-e34de9736c71 # ╟─51bed1cc-c960-46fe-bc09-2b684df3b0cc # ╟─513414c7-0a54-4767-a583-7d779f8fbc55 -# ╟─ee3da94c-5e87-4a0b-8373-c01e339d28aa -# ╟─261620b0-9580-4d9e-b7de-d7972ea549cd -# ╟─6a294790-d294-11ef-270b-5b2152431426 +# ╟─6aa75f5b-1701-437e-ba1c-6583698416d3 +# ╟─e55126ef-e956-464d-8ae0-32b077649f21 # ╟─b872cd69-d534-4b04-bb76-d85bb7ef0ea9 # ╟─1ba1939d-9986-4b97-9273-4f2434f1d385 # ╟─b426df32-5629-4773-b862-101cfbd82d42 # ╟─181ade96-8e1e-4186-9227-c1561352529d +# ╟─6a2af90a-d294-11ef-07bd-018326577791 +# ╟─d484c41d-9834-4528-bf47-93ab4e35ebaa +# ╟─6a2b1106-d294-11ef-0d64-dbc26ba3eb44 # ╟─6a29d548-d294-11ef-1361-ad2230cad02b +# ╟─ebcfcd1b-7fc8-42b7-a35e-4530f798cfdf +# ╟─188b5bea-6765-4dcf-9369-3b1fdbe94494 # ╟─6a29e25e-d294-11ef-15ce-5bf3d8cdb64c # ╟─6a29f1c2-d294-11ef-147f-877f99e5b57c # ╟─6a2a000e-d294-11ef-17d6-bdcddeedc65d +# ╠═5483148f-b385-4afa-ad85-70efe08ba299 # ╟─6a2a0f18-d294-11ef-02c2-ef117377ca66 # ╟─6a2a1daa-d294-11ef-2a67-9f2ac60a14c5 # ╟─6a2a2af2-d294-11ef-0072-bdc3c6f95bb3 @@ -3238,68 +2491,13 @@ version = "1.9.2+0" # ╟─e3f9e571-2248-403c-8ab8-f6b99597f595 # ╟─90f691ad-046c-4595-99b0-19a1d6cb599e # ╟─6a2a9faa-d294-11ef-1284-cfccb1da444e -# ╟─6a2aad42-d294-11ef-3129-3be5be8c82d6 -# ╟─6a2abb16-d294-11ef-0243-d376e8a39bb0 -# ╟─54cb380d-8864-4158-ba68-55027ae68971 -# ╠═9da43d0f-e605-41b7-9bc6-db5be95bc87f -# ╟─280e0819-674e-4a58-854d-5a66e0777074 -# ╠═e47b6eb6-2bb3-4c2d-bda6-f1535f2f94c4 -# ╟─e55126ef-e956-464d-8ae0-32b077649f21 -# ╟─f67136ff-f33c-436e-823b-9c530d257ab0 -# ╠═9c751f8e-f7ed-464f-b63c-41e318bbff2d -# ╟─d1d2bb84-7083-435a-9c19-4c02074143e3 -# ╟─6922b899-499e-4f73-a6be-ca427fdc14ea -# ╟─6a2af90a-d294-11ef-07bd-018326577791 -# ╟─6a2b1106-d294-11ef-0d64-dbc26ba3eb44 -# ╟─d484c41d-9834-4528-bf47-93ab4e35ebaa -# ╟─5ea6cefa-621a-4afc-bf9e-02d42b1d53f8 -# ╟─6a2b2d44-d294-11ef-33ba-15db357708b1 -# ╟─e9b32823-efd2-4a27-b529-4f49752c00bb -# ╟─56b9aba3-6ead-498c-8670-ad93a1953b2a -# ╟─c28b7130-f7fb-41ee-852e-9964b091d7fb -# ╠═3a903a4d-1fb0-4566-8151-9c86dfc40ceb -# ╠═e99e7650-bb72-4576-8f2a-c3994533b644 -# ╟─7a624d2f-812a-47a0-a609-9fe299de94f5 -# ╟─51829800-1781-49ae-8ee7-ac15c0bfcb88 -# ╟─f956e217-3dce-446a-8660-25f2c9cb05e2 -# ╠═2c90eee1-b5d9-434d-bccc-64de8b458a48 -# ╠═9d82af33-8e91-48e9-8c34-fa6ea31492c2 -# ╟─69eaf045-e766-4b7f-a9e8-8eac674ca2ae -# ╟─6a2b3ba4-d294-11ef-3c28-176be260cb15 -# ╟─188b5bea-6765-4dcf-9369-3b1fdbe94494 -# ╟─ebcfcd1b-7fc8-42b7-a35e-4530f798cfdf -# ╟─ee4006aa-e54b-4501-93ee-60b34bdf5c7b -# ╟─6a2b9676-d294-11ef-241a-89ff7aa676f9 -# ╟─9c5d7c89-f65c-4f52-9e49-14692bed2452 -# ╟─6a2bb18a-d294-11ef-23bb-99082caf6e01 -# ╟─6a2bd3ac-d294-11ef-0543-6fe202ca35b6 -# ╟─6a2bf332-d294-11ef-1ff1-cdbfb7732cf1 -# ╟─6a2c008e-d294-11ef-2f07-11cdfb2bddca -# ╟─6a2c11e6-d294-11ef-173b-23fc6dbfefca -# ╟─6a2c229e-d294-11ef-2f24-ebe43cbfbfa4 -# ╟─6a2c3036-d294-11ef-23cb-c3b36c475e8f -# ╟─6a2c4058-d294-11ef-2312-d9c672d49701 -# ╟─6a2c505c-d294-11ef-1c92-c1b0e9d50da5 -# ╟─7c8b1add-085a-41ba-9d6c-b26d3eef22e4 -# ╟─6a2c5e08-d294-11ef-213d-97bcfa16eb5a -# ╟─6a2c7230-d294-11ef-05a2-3ff2f65d10e0 -# ╟─6a2c7f5a-d294-11ef-2e17-9108a39df280 -# ╟─6a2c8f4a-d294-11ef-213c-dfa929a403bc -# ╟─6a2ca496-d294-11ef-0043-1f350b36773e +# ╠═496a4329-bee2-4265-8854-e6019b358d01 +# ╠═679ef9d1-cc1c-4fc1-bf82-caa967c196c2 # ╟─47842de0-d17e-460e-b3b7-b2e642569e25 # ╟─b273c8bc-3819-4f63-801a-acf0ee78ef1d -# ╟─4bfd141f-fe2f-46a0-aa35-872cab45ea00 -# ╟─f2969d91-4a5b-4665-9fa5-521db750302f -# ╟─7dd9a456-9dca-47c8-98c5-51f87f28e6a4 -# ╟─b2820dfd-b3ca-477b-8cb7-c430e0fe18dd -# ╟─664d4183-edb6-4818-a44b-bf4c0a22a33c -# ╟─ecb036da-a0a2-4919-b1aa-bc33b6ba7e73 -# ╟─de08c2a1-c5e3-4add-8b22-2c633247da48 # ╟─6a2cb25e-d294-11ef-1d88-1fc784b33df0 -# ╟─1edae118-dcc7-4169-95cf-f36025f2c336 -# ╟─275a9a69-3135-4cbd-8a35-b1abee4af83f -# ╟─6a2ccd16-d294-11ef-22ee-a5cff62ccd9c -# ╟─f5d8d021-3157-464f-93a2-b3054779e55f +# ╟─ee3da94c-5e87-4a0b-8373-c01e339d28aa +# ╟─261620b0-9580-4d9e-b7de-d7972ea549cd # ╟─1f92c406-6792-4af6-9132-35efd8223bc5 # ╠═b32d2d36-f3a5-406f-adf6-c7b8ebe6cc77 # ╠═17182feb-2c58-4485-aadc-114003376607 @@ -3307,5 +2505,14 @@ version = "1.9.2+0" # ╠═caba8eee-dfea-45bc-a8a7-1dd20a1fa994 # ╠═3987d441-b9c8-4bb1-8b2d-0cc78d78819e # ╟─7a764a14-a5df-4f76-8836-f0a571fc3519 +# ╟─c28b7130-f7fb-41ee-852e-9964b091d7fb +# ╠═9da43d0f-e605-41b7-9bc6-db5be95bc87f +# ╠═e47b6eb6-2bb3-4c2d-bda6-f1535f2f94c4 +# ╠═d1d2bb84-7083-435a-9c19-4c02074143e3 +# ╠═9c751f8e-f7ed-464f-b63c-41e318bbff2d +# ╠═3a903a4d-1fb0-4566-8151-9c86dfc40ceb +# ╠═e99e7650-bb72-4576-8f2a-c3994533b644 +# ╠═7a624d2f-812a-47a0-a609-9fe299de94f5 +# ╠═758a13d7-5c4b-498c-8c37-c3ba9b4daf48 # ╟─00000000-0000-0000-0000-000000000001 # ╟─00000000-0000-0000-0000-000000000002 diff --git a/mlss/Factor Graphs.jl b/mlss/Factor Graphs.jl index f2eca43e..b8a0f9e0 100644 --- a/mlss/Factor Graphs.jl +++ b/mlss/Factor Graphs.jl @@ -1,5 +1,5 @@ ### A Pluto.jl notebook ### -# v0.20.19 +# v0.20.21 #> [frontmatter] #> image = "https://github.com/bmlip/course/blob/v2/assets/figures/ffg-example-1.png?raw=true" @@ -12,104 +12,21 @@ using Markdown using InteractiveUtils -# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error). -macro bind(def, element) - #! format: off - return quote - local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end - local el = $(esc(element)) - global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el) - el - end - #! format: on -end - # ╔═╡ 5a8dcadb-f0c2-4fb0-b8cd-db8cf49cc292 using BmlipTeachingTools # ╔═╡ 965a08f4-d294-11ef-0604-1586ff37c0d4 using Plots, LinearAlgebra, LaTeXStrings -# ╔═╡ 2cb7d369-e7fd-4d66-8321-66a9197a26bd -using RxInfer, Random - -# ╔═╡ 387f55ba-0fa4-4171-a405-3b1fb4e6b586 -using StableRNGs - # ╔═╡ 96547560-d294-11ef-0fa7-6b6489f7baba title("Factor Graphs") # ╔═╡ af24aa27-b0a1-4c9b-aee0-0e5143d2f47e PlutoUI.TableOfContents() -# ╔═╡ 9654ea3e-d294-11ef-335c-657af1ceaf19 -md""" -## Preliminaries - -##### Goal - - * Introduction to Forney-style factor graphs and message passing-based inference - -##### Materials - - * Mandatory - - * These lecture notes - * Loeliger (2007), [The factor graph approach to model based signal processing](https://github.com/bmlip/course/blob/main/assets/files/Loeliger-2007-The-factor-graph-approach-to-model-based-signal-processing.pdf), pp. 1295-1302 (until section V) - * Optional - - * Frederico Wadehn (2015), [Probabilistic graphical models: Factor graphs and more](https://www.youtube.com/watch?v=Fv2YbVg9Frc&t=31) video lecture (**recommended**) - * References - - * Forney (2001), [Codes on graphs: normal realizations](https://github.com/bmlip/course/blob/main/assets/files/Forney-2001-Codes-on-graphs-normal-realizations.pdf) - -""" - -# ╔═╡ 96552348-d294-11ef-16d8-b53563054687 -md""" -## Why Factor Graphs? - - -A probabilistic inference task derives most of its computational complexity from the need to perform marginalization, i.e., integrating (or summing) over latent or nuisance variables. - -For example, for a model ``p(x_1,x_2,x_3,x_4,x_5)``, the inference task ``p(x_2|x_3)`` involves computing - -```math -p(x_2|x_3) = \frac{p(x_2,x_3)}{p(x_3)} = \frac{\int \cdots \int p(x_1,x_2,x_3,x_4,x_5) \, \mathrm{d}x_1 \mathrm{d}x_4 \mathrm{d}x_5}{\int \cdots \int p(x_1,x_2,x_3,x_4,x_5) \, \mathrm{d}x_1 \mathrm{d}x_2 \mathrm{d}x_4 \mathrm{d}x_5} -``` - -Because marginalization operations (sums or integrals) scale poorly with the number of variables, a phenomenon known as the [curse of dimensionality](https://en.wikipedia.org/wiki/Curse_of_dimensionality), we must often leverage the model’s conditional independence structure to reduce computational complexity and obtain feasible inference procedures. - -Factor graphs provide a computationally efficient approach to solving inference problems **if the probabilistic model can be factorized**. - -""" - - -# ╔═╡ 9655959e-d294-11ef-0ca6-5f20aa579e91 -md""" - $(HTML("Factorization helps.")) For instance, if ``p(x_1,x_2,x_3,x_4,x_5) = p(x_1)p(x_2,x_3)p(x_4)p(x_5|x_4)``, then - -```math -\begin{align} -p(x_2|x_3) &= \frac{\int \cdots \int p(x_1)p(x_2,x_3)p(x_4)p(x_5|x_4) \, \mathrm{d}x_1 \mathrm{d}x_4 \mathrm{d}x_5}{\int \cdots \int p(x_1)p(x_2,x_3)p(x_4)p(x_5|x_4) \, \mathrm{d}x_1 \mathrm{d}x_2 \mathrm{d}x_4 \mathrm{d}x_5} \\ - &= \frac{p(x_2,x_3)}{\int p(x_2,x_3) \mathrm{d}x_2} -\end{align} -``` - -which is computationally much cheaper than the general case above. - -In this lesson, we discuss how computationally efficient inference in *factorized* probability distributions can be automated by message passing-based inference in factor graphs. - -""" - -# ╔═╡ 05db1eab-5b63-4ab9-8b4a-ab2cde554295 -md""" -# Forney-style Factor Graphs -""" - # ╔═╡ 9655b2c2-d294-11ef-057f-9b3984064411 md""" -## Factor Graph Construction Rules +## Forney-style Factor Graph Construction Rules """ @@ -150,7 +67,6 @@ A **configuration** is an assignment of values to all variables. A configuration md""" ## Equality Nodes for Branching Points -Note that a variable can appear in maximally two factors in an FFG (since an edge has only two end points). """ @@ -185,7 +101,7 @@ where ```math f_=(x_2,x_2^\prime,x_2^{\prime\prime}) \triangleq \delta(x_2-x_2^\prime)\, \delta(x_2-x_2^{\prime\prime}) ``` -is a so-called **equality** (or branching) node. +is a so-called **equality** (or branching) node. The constraint ``f_=(x,x^\prime,x^{\prime\prime})`` enforces that ``x=x^\prime=x^{\prime\prime}`` **for every valid configuration**. """, md""" ![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-wEquality-node.png?raw=true) @@ -195,93 +111,6 @@ md""" md""" Note that through introduction of auxiliary variables ``x_2^{\prime}`` and ``x_2^{\prime\prime}`` and a factor ``f_=(x_2,x_2^\prime,x_2^{\prime\prime})``, each variable in ``g`` appears in maximally two factors. -The constraint ``f_=(x,x^\prime,x^{\prime\prime})`` enforces that ``x=x^\prime=x^{\prime\prime}`` **for every valid configuration**. - -Since ``f`` is a marginal of ``g``, i.e., - -```math -f(x_1,x_2,x_3,x_4) = \iint g(x_1,x_2,x_2^\prime,x_2^{\prime\prime},x_3,x_4)\, \mathrm{d}x_2^\prime \mathrm{d}x_2^{\prime\prime} -``` - -it follows that any inference problem on ``f`` can be executed by a corresponding inference problem on ``g``, e.g., - -```math -\begin{align*} -f(x_1 \mid x_2) &\triangleq \frac{\iint f(x_1,x_2,x_3,x_4) \,\mathrm{d}x_3 \mathrm{d}x_4 }{ \int\cdots\int f(x_1,x_2,x_3,x_4) \,\mathrm{d}x_1 \mathrm{d}x_3 \mathrm{d}x_4} \\ - &= \frac{\int\cdots\int g(x_1,x_2,x_2^\prime,x_2^{\prime\prime},x_3,x_4) \,\mathrm{d}x_2^\prime \mathrm{d}x_2^{\prime\prime} \mathrm{d}x_3 \mathrm{d}x_4 }{ \int\cdots\int g(x_1,x_2,x_2^\prime,x_2^{\prime\prime},x_3,x_4) \,\mathrm{d}x_1 \mathrm{d}x_2^\prime \mathrm{d}x_2^{\prime\prime} \mathrm{d}x_3 \mathrm{d}x_4} \\ - &= g(x_1 \mid x_2) -\end{align*} -``` -""" - -# ╔═╡ 9656cf72-d294-11ef-03aa-b715dd686c09 -md""" -## Probabilistic Models as Factor Graphs - -FFGs can be used to express conditional independence (factorization) in probabilistic models. - - -For example, the (previously shown) graph for - -```math -f_a(x_1,x_2,x_3) \cdot f_b(x_3,x_4,x_5) \cdot f_c(x_4) -``` - -could represent the probabilistic model - -```math -p(x_1,x_2,x_3,x_4,x_5) = p(x_1,x_2|x_3) \cdot p(x_3,x_5|x_4) \cdot p(x_4) -``` - -""" - -# ╔═╡ 27b9e811-4542-4fe8-86a5-f8ba09074761 -TwoColumn( -md""" -where we identify -```math -\begin{align*} -f_a(x_1,x_2,x_3) &= p(x_1,x_2|x_3) \\ -f_b(x_3,x_4,x_5) &= p(x_3,x_5|x_4) \\ -f_c(x_4) &= p(x_4) -\end{align*} -``` -""", -@htl """ - - - -""") - -# ╔═╡ 9656d850-d294-11ef-21a1-474b07ea7729 -md""" -This factorized probability distribution is represented by the above FFG. - -""" - -# ╔═╡ 9658329c-d294-11ef-0d03-45e6872c4985 -md""" -## Terminating an FFG - -Consider a model - -```math -f(x_1,x_2,y) = f_a(x_1) f_b(x_1,x_2,y) \,. -``` - -In this model, the variables ``x_2`` and ``y`` appear in only one factor. In the corresponding FFG, edges that only connect to one factor are called a **half-edges**. Half-edges typically represent inputs or outputs of the graph, such as observed variables and external control signals. - -In general, any half-edge can be terminated by a node ``f(\cdot) = 1``, since the model - -```math -f_a(x_1) f_b(x_1,x_2,y) \underbrace{f_c(y)}_{=1} \underbrace{f_d(x_2)}_{=1}\,, -``` -is the same model as ``f_a(x_1) f_b(x_1,x_2,y)``. - - -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-terminal-node.png?raw=true) - -An FFG without half-edges is called a Terminated FFG (TFFG). """ @@ -289,9 +118,7 @@ An FFG without half-edges is called a Terminated FFG (TFFG). md""" ## Representing Observations -An observation, say ``y=3``, can be represented by a **delta node** ``f(y)=\delta(y−3)`` to terminate the half-edge for variable ``y``. - -In an FFG, we visualize a delta node by a small black box, +An observation, say ``y=3``, can be represented by a **delta node** ``f(y)=\delta(y−3)`` to terminate the half-edge for variable ``y``. In an FFG, we visualize a delta node by a small black box, """ @@ -345,12 +172,6 @@ We draw here the FFG for the factorized distribution: # ╔═╡ f46473f8-e749-431d-8350-751115f0aaf0 Resource("https://github.com/bmlip/course/blob/v4/assets/figures/ffg-message-passing.png?raw=true", :style => "background: #eee; padding: 1em; border-radius: 1em;") -# ╔═╡ cb3df230-6c7e-41b9-ba13-3c5f8a7fbb62 -md""" -Note that we drew *directed edges* to distinguish between intermediate results ("messages") ``\overrightarrow{\mu}_\bullet(\cdot)`` that flow in the same direction as the arrow of the edge (later to be called: forward messages) from intermediate results ``\overleftarrow{\mu}_\bullet(\cdot)`` that flow in opposite direction (later to be called: backward messages). This is just a notational convenience since an FFG is computationally an undirected graph. For now, only consider the nodes and edges. The messages ``\overrightarrow{\mu}_\bullet(\cdot)`` will be discussed later. - -""" - # ╔═╡ 96570d3e-d294-11ef-0178-c34dda717495 md""" Due to the factorization of ``f(x_1,x_2,\ldots,x_7)`` and the [Generalized Distributive Law](https://en.wikipedia.org/wiki/Generalized_distributive_law), we can decompose the marginalization operation to the following product-of-sums: @@ -363,882 +184,143 @@ Due to the factorization of ``f(x_1,x_2,\ldots,x_7)`` and the [Generalized Distr \end{align*} ``` -which, in case ``x_i`` has ``10`` values, requires a few hundred additions and is therefore computationally (much!) lighter than executing the whole sum ``\sum_{x_1,\ldots,x_7}f(x_1,x_2,\ldots,x_7)`` - - -""" - -# ╔═╡ 9657b088-d294-11ef-3017-e95c4c69b62b -md""" -As an afterthought, note that applying the distributive law in an FFG for inference is analogous to replacing the sum-of-products - -```math -ac + ad + bc + bd -``` - -by the following product-of-sums: - -```math -(a + b)(c + d) \,. -``` - -Which of these two computations is cheaper to execute? - -""" - -# ╔═╡ 0afe3cdc-15ed-4d9a-848a-d1977d051866 -md""" -## Closing-the-Box and Message Passing Interpretations -""" - -# ╔═╡ 96571c34-d294-11ef-11ef-29beeb1f96c2 -md""" -Note that the intermediate result ``\overrightarrow{\mu}_{X_3}(x_3)`` is obtained by multiplying all enclosed factors (``f_a``, ``f_b, f_c``) by the red dashed box, followed by marginalization (summing) over all enclosed variables (``x_1``, ``x_2``), - -```math -\overrightarrow{\mu}_{X_3}(x_3) = \underbrace{\sum_{x_1}\sum_{x_2}}_{\text{enclosed variables}}\underbrace{f_a(x_1)f_b(x_2) f_c(x_1,x_2,x_3) }_{\text{enclosed factors}} -``` - -This operation is known as **Closing-the-Box**. The result is a new **composite node** that holds the factor ``\overrightarrow{\mu}_{X_3}(x_3)``, and is visually represented by the red dashed box in the factor graph. The composite node ``\overrightarrow{\mu}_{X_3}(x_3)`` depends only on the variable(s) that cross the boundary of the box (in this case ``x_3``) and effectively replaces the internal subgraph contained within the red box. -""" - -# ╔═╡ 253d4703-03d6-4961-8c3b-b70d2cbc0710 -md""" -When closing the box around a terminal node, the result is simply the factor associated with that node, since there are no internal variables to marginalize. -""" - -# ╔═╡ a7b1f559-3c34-491e-83e7-ba95c8c22c80 -md""" - -The Closing-the-box operation can alternatively be interpreted as **passing a message** from the newly created composite node to the rest of the graph. For instance, ``\overrightarrow{\mu}_{X_3}(x_3)`` can be understood in two equivalent ways: - * as a factor associated with the composite node that encloses the subgraph inside the red box. - * as a message sent from this composite node to the variable ``x_3``. - -In both interpretations, the internal details of the subgraph are abstracted away, and the composite node effectively summarizes its contribution to the overall inference process. - - -""" - -# ╔═╡ 70736e62-2b6c-4b3a-ab59-7e51522d620b -md""" - -The complete inference process for computing ``\bar{f}(x_3)`` can be interpreted as a **message passing process**. It begins by sending messages from the terminal nodes and then propagates them through the internal nodes of the factor graph. This continues until both the forward and backward messages for ``x_3`` have been computed. The final result, ``\bar{f}(x_3)``, is obtained by multiplying the forward and backward messages, - -```math -\bar{f}(x_3) = \overrightarrow{\mu}_{X_3}(x_3) \cdot \overleftarrow{\mu}_{X_3}(x_3) -``` - -This message-based interpretation enables modular, local inference that scales efficiently with the structure of the factor graph. - - -""" - - -# ╔═╡ 96575dd4-d294-11ef-31d6-b39b4c4bdea1 -md""" -## Sum-Product Messages - -Let's continue with the message passing interpretation of inference in an FFG. Closing the red box around ``f_a``, ``f_b`` and ``f_c`` leads to an outgoing message ``\overrightarrow{\mu}_{X_3}(x_3)`` for node ``f_c``, given by - -```math -\begin{align} -\underbrace{\overrightarrow{\mu}_{X_3}(x_3)}_{\substack{ \text{outgoing} \\ \text{message} }} -&= \sum_{x_1}\sum_{x_2} \underbrace{\overrightarrow{\mu}_{X_1}(x_1) \overrightarrow{\mu}_{X_2}(x_2)}_{\substack{\text{incoming} \\ \text{messages}}} \underbrace{f_c(x_1,x_2,x_3)}_{\text{factor}} -\end{align} -``` - -This recipe holds generally. For a node ``f(y,x_1,\ldots,x_n)`` with incoming messages ``\overrightarrow{\mu}_{X_1}(x_1)``, ``\overrightarrow{\mu}_{X_2}(x_2)``, ``\ldots``,``\overrightarrow{\mu}_{X_n}(x_n)``, the outgoing message is given by ([Loeliger (2007), pg.1299](https://github.com/bmlip/course/blob/main/assets/files/Loeliger-2007-The-factor-graph-approach-to-model-based-signal-processing.pdf)): - -```math -\underbrace{\overrightarrow{\mu}_{Y}(y)}_{\substack{ \text{outgoing}\\ \text{message}}} = \sum_{x_1,\ldots,x_n} \underbrace{\overrightarrow{\mu}_{X_1}(x_1)\cdots \overrightarrow{\mu}_{X_n}(x_n)}_{\substack{\text{incoming} \\ \text{messages}}} \cdot \underbrace{f(y,x_1,\ldots,x_n)}_{\substack{\text{node}\\ \text{function}}} \tag{SP} -``` - -""" - -# ╔═╡ 5cc2016e-0383-448c-bd33-5b3a687b7436 -TwoColumn( -md""" -Equation (SP) is called a **Sum-Product** message, so named because the computation involves evaluating a sum-of-products. Note that all SP messages in an FFG can be computed from information that is **locally available** at each node. -""", -@htl """ - - - -""") - -# ╔═╡ f65f5d0e-2583-4b88-b9f2-5fee15257c05 -md""" - - -If the factor graph for the whole model has no cycles, i.e., the FFG is a tree, then passing SP messages from the terminal nodes to the internal (latent) variables yields exact Bayesian marginals for all hidden variables. This inference method is known as the **Sum-Product** (SP) algorithm. - -However, if the graph contains cycles, one can view it conceptually as an infinite tree by “unrolling” the cycles. In this loopy setting, SP-based inference is not guaranteed to yield exact marginals. Nevertheless, in practice, if we run the SP algorithm for a limited number of iterations (i.e., a finite unrolling), we often obtain high-quality approximate marginals that are sufficient for many inference tasks. -""" - -# ╔═╡ 7009cdc8-892c-499e-b932-b828fa300b6c -keyconcept("", -md""" -For a node ``f(y,x_1,\ldots,x_n)`` with incoming messages ``\overrightarrow{\mu}_{X_1}(x_1)``, ``\overrightarrow{\mu}_{X_2}(x_2)``, ``\ldots``,``\overrightarrow{\mu}_{X_n}(x_n)``, the outgoing message is given by the **sum-product rule**: - -```math -\overrightarrow{\mu}_{Y}(y)= \sum_{x_1,\ldots,x_n} \overrightarrow{\mu}_{X_1}(x_1)\cdots \overrightarrow{\mu}_{X_n}(x_n) \cdot f(y,x_1,\ldots,x_n) -``` -""") - -# ╔═╡ 91f81188-727c-4754-9a07-e754eef8bbe0 -md""" -## Example: Sum-Product Messages for the Equality Node -""" - -# ╔═╡ 0633afea-5e92-4bad-8402-d159c534af81 -TwoColumn( -md""" -As an example, let´s evaluate the SP messages for the **equality node** - -```math -f_=(x,y,z) = \delta(z-x)\delta(z-y) \,. -``` -""", -@htl """ - - - -""") - -# ╔═╡ f11564db-aafc-4df9-b494-4e5ced9bfcfe -md""" - -Given incoming messages ``\overrightarrow{\mu}_{X}(x)`` and ``\overrightarrow{\mu}_{Y}(y)``, the outgoing SP message ``\overrightarrow{\mu}_{Z}(z)`` to edge ``z`` is given by -```math -\begin{align*} -\overrightarrow{\mu}_{Z}(z) &= \iint \overrightarrow{\mu}_{X}(x) \overrightarrow{\mu}_{Y}(y) \,\delta(z-x)\delta(z-y) \,\mathrm{d}x \mathrm{d}y \\ - &= \overrightarrow{\mu}_{X}(z) \int \overrightarrow{\mu}_{Y}(y) \,\delta(z-y) \,\mathrm{d}y \\ - &= \overrightarrow{\mu}_{X}(z) \overrightarrow{\mu}_{Y}(z) -\end{align*} -``` - -By symmetry, this also implies (for the same equality node) that - -```math -\begin{align*} -\overleftarrow{\mu}_{X}(x) &= \overrightarrow{\mu}_{Y}(x) \overleftarrow{\mu}_{Z}(x) \\ -\overleftarrow{\mu}_{Y}(y) &= \overrightarrow{\mu}_{X}(y) \overleftarrow{\mu}_{Z}(y)\,. -\end{align*} -``` - -It follows that message passing through an equality node is analogous to applying Bayes’ rule: two information sources are fused by multiplying their corresponding probability distributions. - -""" - -# ╔═╡ 9651f976-b834-4b81-8810-649f0290969d -md""" -# RxInfer: A Toolbox for Automated Bayesian inference -""" - -# ╔═╡ 96587a66-d294-11ef-2c7a-9fd7bea76582 -md""" -## Automating Bayesian Inference by Message Passing - -The foregoing message update rules can be worked out in closed form and put into tables (e.g., see Tables 1 through 6 in [Loeliger (2007)](https://github.com/bmlip/course/blob/main/assets/files/Loeliger-2007-The-factor-graph-approach-to-model-based-signal-processing.pdf) for many standard factors such as essential probability distributions and operations such as additions, fixed-gain multiplications, and branching (equality nodes). - -In the optional slides below, we have worked out a few more update rules for the [addition node](#sp-for-addition-node) and the [multiplication node](#sp-for-multiplication-node). - -If the update rules for all node types in a graph have been tabulated, then inference by message passing comes down to executing a set of table-lookup operations, thus creating a completely **automatable Bayesian inference framework**. - -In our research lab [BIASlab](http://biaslab.org) (FLUX 7.060), we are developing [RxInfer](http://rxinfer.com), which is a (Julia) toolbox for automating Bayesian inference by message passing in a factor graph. - -In general, a code package that automates Bayesian inference is called a [Probabilistic Programming](https://en.wikipedia.org/wiki/Probabilistic_programming) Language (PPL). RxInfer is a PPL that automates inference through message passing-based inference in a factor graph. - -""" - -# ╔═╡ 89e2757e-a09f-40c6-8dd7-9b4b4d232e17 -md""" -![](https://github.com/bmlip/course/blob/main/assets/figures/RxInfer-gif.gif?raw=true) -""" - -# ╔═╡ c4b5b124-e52a-41fc-b27e-a58181622e5c -md""" -The figure above (a screen recording from the [RxInfer webpage](http://rxinfer.com)) is an animated GIF illustrating how RxInfer operates. The model is represented as a graph in which each node passes messages to its neighbors. When messages meet on an edge, the belief about the variable associated with that edge is updated. -""" - -# ╔═╡ a1c957c1-69b7-4178-ab59-c0b2439bb01a -code_example("Bayesian Linear Regression by Message Passing"; big=true) - -# ╔═╡ 9658c106-d294-11ef-01db-cfcff611ed81 -md""" -Assume we want to estimate some function ``f: \mathbb{R}^D \rightarrow \mathbb{R}`` from a given data set ``D = \{(x_1,y_1), \ldots, (x_N,y_N)\}``, with ``x_i \in \mathbb{R}^D``, ``y_i \in \mathbb{R}``. - -""" - -# ╔═╡ 96594d44-d294-11ef-22b8-95165fb08ce4 -md""" -### Model Specification - -We will assume a linear model with white Gaussian noise and a Gaussian prior on the coefficients ``w``: - -```math -\begin{align*} - y_i &= w^T x_i + \epsilon_i \\ - \epsilon_i &\sim \mathcal{N}(0, \sigma^2) \\ - w &\sim \mathcal{N}(0,\Sigma) -\end{align*} -``` - -or equivalently - -```math -\begin{align*} -p(w,\epsilon,D) &= \overbrace{p(w)}^{\text{weight prior}} \prod_{i=1}^N \overbrace{p(y_i\,|\,x_i,w,\epsilon_i)}^{\text{data-generating}} \overbrace{p(\epsilon_i)}^{\text{noise prior}} \\ - &= \mathcal{N}(w\,|\,0,\Sigma) \prod_{i=1}^N \delta(y_i - w^T x_i - \epsilon_i) \mathcal{N}(\epsilon_i\,|\,0,\sigma^2) -\end{align*} -``` - -""" - -# ╔═╡ 96597ce0-d294-11ef-3478-25c6bbef601e -md""" -### Inference (parameter estimation) - -We are interested in inferring the posterior ``p(w|D)``. We will execute inference by message passing on the FFG for the model. - -""" - -# ╔═╡ 965998a8-d294-11ef-1d18-85876e3656c5 -md""" -The left figure shows the factor graph for this model for one observation ``(x,y)``. The figure on the right shows the message passing scheme. - -""" - -# ╔═╡ 284a9dd5-1e26-4fd3-bb58-6e7ac0a0872f -Resource("https://github.com/bmlip/course/blob/v2/assets/figures/ffg-bayesian-linear-regression.png?raw=true", :width=>500, :style => "background: white; padding: 1em; border-radius: 1em;") - -# ╔═╡ 42dd67e6-eb0f-4368-9947-47de229f7be1 -md""" -### Modeling a Polynomial -The section above is about a **linear model** for a function on ``\mathbb{R}^D \rightarrow \mathbb{R}``. - -In the example below, we will model a function ``f: \mathbb{R} \rightarrow \mathbb{R}`` with a **polynomial** of degree ``K``. - - -We can use a linear model for a polynomial! The trick is to turn each training data point ``(x,y) \in \mathbb{R} \times \mathbb{R}`` into a feature vector: - -```math -(\,[1, x, x^2, x^3, \cdots, x^K], \ y\,) \quad \in \, \mathbb{R}^{K+1} \times \mathbb{R} -``` - -And we have now transformed our polynomial model on ``\mathbb{R}`` to a linear model on ``\mathbb{R}^{K+1}``. - - -""" - -# ╔═╡ 480165f9-33d9-4db1-bf05-8d99f0d9fb3e -md""" -### Generate the Data Set -For details, see [the Appendix](#Code). - -""" - -# ╔═╡ 1c9c7994-672c-42a3-8ae7-8ce092ada9f0 -begin - N_bond = @bindname Nsamples Slider(1:30; default=20, show_value=true) -end - -# ╔═╡ 8a2019af-9500-42c5-8408-ff93104a2d79 -md""" -Now, we want to find the parameter ``w`` for our polynomial ``f_w`` to model this data. -""" - -# ╔═╡ 1a40ef8d-d677-4bf0-9186-18c5aa43a849 - - -# ╔═╡ 965a1df0-d294-11ef-323c-3da765f1104a -md""" -### Infer Solution with RxInfer - -Now build the factor graph in RxInfer, and perform sum-product message passing to generate a posterior for the weights. - -""" - -# ╔═╡ fd338a30-9622-405a-96fa-caca6bd4ccfb -@model function linear_regression(y,x, Nsamples, Σ, σ²) - - w ~ MvNormalMeanCovariance(zeros(3),Σ) - - for i in 1:Nsamples - y[i] ~ NormalMeanVariance(dot(w, x[i]), σ²) - end -end - -# ╔═╡ 1070063a-ef85-4527-ae82-1f01c1a506ff -prior_Σ = 1e5 * Diagonal(I,3); # prior for the weights - -# ╔═╡ 9431bc9a-bd83-4e4d-b64d-0571c1d01c87 -md""" -It worked! Now we have a **posterior distribution** for ``w``: -""" - -# ╔═╡ b3262127-69e0-4efb-875b-074d1d70437c - - -# ╔═╡ fb61c774-34a3-493a-b149-c870993b6d46 -md""" -#### Plot the Results - -Let's sample ``10`` typical values for the weights ``w`` from this posterior distribution, and plot the corresponding curves ``f_w: x \mapsto w^Tx`` in the scatter plot again. - - -""" - -# ╔═╡ 5bcefd5f-4cd2-4cfe-8c1f-1129e5020d9a -N_bond - -# ╔═╡ 1832bffd-2729-4d3f-86f4-0e2d9ab26ba3 -md""" -Notice how the samples of the functions ``f_w`` lie closer together as we get more observations! -""" - -# ╔═╡ 4a10044c-e044-43e1-bd44-847f56019061 -keyconcept( - "", - md""" - “In `RxInfer`, once the model is specified and the observations are provided, Bayesian inference is carried out automatically via the `infer()` function.” - """ -) - -# ╔═╡ 965a6c20-d294-11ef-1c91-4bd237afbd20 -md""" -## Final thoughts: Modularity and Abstraction - -The great Michael Jordan (no, not [this one](https://youtu.be/cuLprHh_BRg), but [this one](https://people.eecs.berkeley.edu/~jordan/)), wrote: - -> "I basically know of two principles for treating complicated systems in simple ways: the first is the principle of **modularity** and the second is the principle of **abstraction**. I am an apologist for computational probability in machine learning because I believe that probability theory implements these two principles in deep and intriguing ways — namely through factorization and through averaging. Exploiting these two mechanisms as fully as possible seems to me to be the way forward in machine learning." — Michael Jordan, 1997 (quoted in [Fre98](https://mitpress.mit.edu/9780262062022/)). - -Factor graphs capture these ideas elegantly—both visually and computationally. - -**Visually**, the graph structure displays the modularity of conditional independencies in the model. Each node encapsulates internal complexity, and by closing the box, we can hierarchically move to higher levels of abstraction. - -**Computationally**, message-passing inference exploits the distributive law to avoid unnecessary computations. - -Although RxInfer is still under active development, my prediction is that within 5–10 years, RxInfer—or a comparable toolbox—will be able to automate Bayesian inference for virtually any interesting probabilistic model you can conceive. In principle, you will then have all the tools needed to implement the four-step Bayesian ML recipe—model specification, parameter learning, model evaluation, and application—for any (Bayesian) information processing problem. - - -""" - -# ╔═╡ fa5bdb1c-4412-48cc-950c-9ed92b4c9f76 -md""" -# Summary -""" - -# ╔═╡ be670693-2036-46cb-8452-a2d0e1bf1172 -keyconceptsummary() - -# ╔═╡ 25492eea-e649-43f9-b71f-ac6d1a80d0ee -exercises(header_level=1) - -# ╔═╡ a5cd774f-57ad-4cb5-86c0-35987aa6e221 -md""" -##### Message Passing in a State Space Model (*) -""" - -# ╔═╡ b6de3f00-d3b8-44d8-b72a-48cd5628b607 -TwoColumn(md""" Consider the following state-space model: - -```math -\begin{align*} -z_k &= A z_{k-1} + w_k \\ -x_k &= C z_k + v_k -\end{align*} -``` - -where ``k=1,2,\ldots,n`` is the time step counter; ``z_k`` is an *unobserved* state sequence; ``x_k`` is an *observed* sequence; ``w_k \sim \mathcal{N}(0,\Sigma_w)`` and ``v_k \sim \mathcal{N}(0,\Sigma_v)`` are (unobserved) state and observation noise sequences respectively; ``z_0 \sim \mathcal{N}(0,\Sigma_0)`` is the initial state and ``A``, ``C``, ``\Sigma_v``,``\Sigma_w`` and ``\Sigma_0`` are known parameters. """, -@htl """ - -""") - - -# ╔═╡ 05375a01-4d1b-44cc-b1c4-a5eb4b6c5c5b -md""" -- (a) Rewrite the state-space equations as a set of conditional probability distributions. - -```math -\begin{align*} - p(z_k|z_{k-1},A,\Sigma_w) &= \ldots \\ - p(x_k|z_k,C,\Sigma_v) &= \ldots \\ - p(z_0|\Sigma_0) &= \ldots -\end{align*} -``` - -- (b) Define ``z^n \triangleq (z_0,z_1,\ldots,z_n)``, ``x^n \triangleq (x_1,\ldots,x_n)`` and ``\theta=\{A,C,\Sigma_w,\Sigma_v\}``. Now write out the generative model ``p(x^n,z^n|\theta)`` as a product of factors. - -- (c) We are interested in estimating ``z_k`` from a given estimate for ``z_{k-1}`` and the current observation ``x_k``, i.e., we are interested in computing ``p(z_k|z_{k-1},x_k,\theta)``. Can ``p(z_k|z_{k-1},x_k,\theta)`` be expressed as a Gaussian distribution? Explain why or why not in one sentence. - -- (d) Copy the graph onto your exam paper and draw the message passing schedule for computing ``p(z_k|z_{k-1},x_k,\theta)`` by drawing arrows in the factor graph. Indicate the order of the messages by assigning numbers to the arrows. - -- (e) Now assume that our belief about parameter ``\Sigma_v`` is instead given by a distribution ``p(\Sigma_v)`` (rather than a known value). Adapt the factor graph drawing of the previous answer to reflect our belief about ``\Sigma_v``. -""" - -# ╔═╡ 45251c19-6eae-41e7-b0ed-8bd70a67d4e0 -ex_d_sol = TwoColumn( - md""" - - (d) Copy the graph onto your exam paper and draw the message passing schedule for computing ``p(z_k|z_{k-1},x_k,\theta)`` by drawing arrows in the factor graph. Indicate the order of the messages by assigning numbers to the arrows. - - Some permutations of this order are also possible. The most important thing here is that you recognize the tree with ``Z_k`` as a root of the tree and pass messages from the terminals (e.g., ``Z_{k-1}``, ``X_k``, etc.) towards the root. - """, - @htl """ - - """); - -# ╔═╡ 206c34b3-1873-460b-911e-f2cd4f8886af -hide_solution( -md""" - -- (a) Rewrite the state-space equations as a set of conditional probability distributions. - -```math -\begin{align*} - p(z_k|z_{k-1},A,\Sigma_w) &= \ldots \\ - p(x_k|z_k,C,\Sigma_v) &= \ldots \\ - p(z_0|\Sigma_0) &= \ldots -\end{align*} -``` - -This is a linear system with only Gaussian source signals (``w_k`` and ``v_k``), hence the distributions for ``z_k`` and ``x_k`` will also be Gaussian. As a result, we only need to compute the mean and covariance matrix. We begin with the mean for ``p(z_k|z_{k-1},A,\Sigma_w)``: - - -```math -\begin{align*} - E[z_k|z_{k-1},A,\Sigma_w] &= E[A z_{k-1} + w_k|z_{k-1},A,\Sigma_w] \\ - &= E[A z_{k-1}|z_{k-1},A] + E[w_k|\Sigma_w] \\ - &= A z_{k-1} + 0 - \end{align*} -``` - -And now the variance: - - -```math -\begin{align*} - V[z_k|z_{k-1},A,\Sigma_w] &= E[(z_k - E[z_k])(z_k-E[z_k])^T \,|\,z_{k-1},A,\Sigma_w ] \\ &= E[(\overbrace{A z_{k-1} + w_k}^{z_k} - \overbrace{A z_{k-1}}^{E[z_k]})(A z_{k-1} + w_k-A z_{k-1})^T|z_{k-1},A,\Sigma_w] \\ - &= E[w_k w_k^T|\Sigma_w] \\ - &= \Sigma_w - \end{align*} -``` - -You can execute similar computations for the other distributions, leading to - - -```math -\begin{align*} - p(z_k|z_{k-1},A,\Sigma_w) &= \mathcal{N}(z_k|A z_{k-1},\Sigma_w) \\ - p(x_k|z_k,C,\Sigma_v) &= \mathcal{N}(x_k|C z_k,\Sigma_v) \\ - p(z_0|\Sigma_0) &= \mathcal{N}(z_0|0,\Sigma_0) -\end{align*} -``` - -- (b) Define ``z^n \triangleq (z_0,z_1,\ldots,z_n)``, ``x^n \triangleq (x_1,\ldots,x_n)`` and ``\theta=\{A,C,\Sigma_w,\Sigma_v\}``. Now write out the generative model ``p(x^n,z^n|\theta)`` as a product of factors. - -```math -\begin{align*} -p(x^n,z^n|\theta) &= p(z_0|\Sigma_0) \prod_{k=1}^n p(x_k|z_k,C,\Sigma_v) \,p(z_k|z_{k-1},A,\Sigma_w) \\ - &= \mathcal{N}(z_0|0,\Sigma_0) \prod_{k=1}^n \mathcal{N}(x_k|C z_k,\Sigma_v) \,\mathcal{N}(z_k|A z_{k-1},\Sigma_w) -\end{align*} -``` - -- (c) We are interested in estimating ``z_k`` from a given estimate for ``z_{k-1}`` and the current observation ``x_k``, i.e., we are interested in computing ``p(z_k|z_{k-1},x_k,\theta)``. Can ``p(z_k|z_{k-1},x_k,\theta)`` be expressed as a Gaussian distribution? Explain why or why not in one sentence. - -Yes, since the generative model ``p(x^n,z^n|\theta)`` is (one big) Gaussian. - - $ex_d_sol - -- (e) Now assume that our belief about parameter ``\Sigma_v`` is instead given by a distribution ``p(\Sigma_v)`` (rather than a known value). Adapt the factor graph drawing of the previous answer to reflects our belief about ``\Sigma_v``. - -For answer, see drawing for answer (d). - -""") - -# ╔═╡ a6e155eb-7376-4e57-8e63-628934e14e78 -md""" -##### Messages for the Addition Node (*) - -""" - -# ╔═╡ 9dc870d7-a5f3-447c-96ee-ad23199bc253 -TwoColumn( -md""" -Consider an addition node - -```math -f_+(x,y,z) = \delta(z-x-y) -``` -- Derive an expression for the outgoing message ``\overrightarrow{\mu}_{Z}(z)`` in terms of the incoming messages ``\overrightarrow{\mu}_{X}(\cdot)`` and ``\overrightarrow{\mu}_{Y}(\cdot)``. - -""", - -@htl """ - - - -""") - -# ╔═╡ e8a35c28-6d6d-4066-8251-f091f28622a9 -hide_solution( -md""" - -We use the sum-product rule to compute - - -```math -\begin{align*} - \overrightarrow{\mu}_{Z}(z) &= \iint \overrightarrow{\mu}_{X}(x) \overrightarrow{\mu}_{Y}(y) \,\delta(z-x-y) \,\mathrm{d}x \mathrm{d}y \\ - &= \int \overrightarrow{\mu}_{X}(x) \overrightarrow{\mu}_{Y}(z-x) \,\mathrm{d}x \,, - \end{align*} -``` - -i.e., ``\overrightarrow{\mu}_{Z}`` is the convolution of the messages ``\overrightarrow{\mu}_{X}`` and ``\overrightarrow{\mu}_{Y}``. - - """) - -# ╔═╡ 965a8a1a-d294-11ef-1d2f-65abf76665e8 -md""" -# Optional Slides - -""" - -# ╔═╡ 965aa14c-d294-11ef-226f-65d587fefa64 -md""" -## $(HTML("Sum-Product Messages for Multiplication Nodes")) -""" - -# ╔═╡ 56e8a1bd-ef80-4265-b926-e5e9e085b72f -TwoColumn( -md""" -Next, let us consider a **multiplication** by a fixed (invertible matrix) gain ``f_A(x,y) = \delta(y-Ax)`` -""", -md""" -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-gain-node.png?raw=true) -""" -) - -# ╔═╡ 965ab77c-d294-11ef-2510-95b1a998589f -md""" -```math -\begin{align*} -\overrightarrow{\mu}_{Y}(y) &= \int \overrightarrow{\mu}_{X}(x) \,\delta(y-Ax) \,\mathrm{d}x \\ -&= \int \overrightarrow{\mu}_{X}(x) \,|A|^{-1}\delta(x-A^{-1}y) \,\mathrm{d}x \\ -&= |A|^{-1}\overrightarrow{\mu}_{X}(A^{-1}y) \,. -\end{align*} -``` - -""" - -# ╔═╡ 965af708-d294-11ef-112c-f5470031dbbe -md""" -For a Gaussian message input message ``\overrightarrow{\mu}_{X}(x) = \mathcal{N}(x|\overrightarrow{m}_{X},\overrightarrow{V}_{X})``, the output message is also Gaussian with - -```math -\begin{align*} -\overrightarrow{m}_{Y} = A\overrightarrow{m}_{X} \,,\,\text{and}\,\, -\overrightarrow{V}_{Y} = A\overrightarrow{V}_{X}A^T -\end{align*} -``` - -since - -```math -\begin{align*} -\overrightarrow{\mu}_{Y}(y) &= |A|^{-1}\overrightarrow{\mu}_{X}(A^{-1}y) \\ - &\propto \exp \left( -\frac{1}{2} \left( A^{-1}y - \overrightarrow{m}_{X}\right)^T \overrightarrow{V}_{X}^{-1} \left( A^{-1}y - \overrightarrow{m}_{X}\right)\right) \\ - &= \exp \big( -\frac{1}{2} \left( y - A\overrightarrow{m}_{X}\right)^T \underbrace{A^{-T}\overrightarrow{V}_{X}^{-1} A^{-1}}_{(A \overrightarrow{V}_{X} A^T)^{-1}} \left( y - A\overrightarrow{m}_{X}\right)\big) \\ - &\propto \mathcal{N}(y| A\overrightarrow{m}_{X},A\overrightarrow{V}_{X}A^T) \,. -\end{align*} -``` - -""" - -# ╔═╡ 965b11a4-d294-11ef-1d04-dbdf39ce91a3 -md""" -**Exercise**: Prove that, for the same factor ``\delta(y-Ax)`` and Gaussian messages, the (backward) sum-product message ``\overleftarrow{\mu}_{X}`` is given by - -```math -\begin{align*} -\overleftarrow{\xi}_{X} &= A^T\overleftarrow{\xi}_{Y} \\ -\overleftarrow{W}_{X} &= A^T\overleftarrow{W}_{Y}A -\end{align*} -``` - -where ``\overleftarrow{\xi}_X \triangleq \overleftarrow{W}_X \overleftarrow{m}_X`` and ``\overleftarrow{W}_{X} \triangleq \overleftarrow{V}_{X}^{-1}`` (and similarly for ``Y``). - -""" - -# ╔═╡ 965b25ac-d294-11ef-0b9a-9d5a50a76069 -md""" -## $(HTML("Code example: Gaussian forward and backward messages for the Addition node")) - -""" - -# ╔═╡ 8dd6874c-12b7-47b6-b589-009849198024 -TwoColumn( -md""" -Let's calculate the Gaussian forward and backward messages for the addition node in RxInfer. -""", -md""" -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-addition-node.png?raw=true) -""" -) - -# ╔═╡ bfbf3d09-23f5-4f54-96f6-bfe536cfc228 -md"Forward message on ``Z``:" - -# ╔═╡ e7e4b6d0-bdf0-4a93-9a73-7971e6e33065 -@call_rule typeof(+)(:out, Marginalisation) (m_in1 = NormalMeanVariance(1.0, 1.0), m_in2 = NormalMeanVariance(2.0, 1.0)) - -# ╔═╡ 2f5415e5-70b1-47ea-9790-7ac953bca538 -md"Backward message on ``X``:" - -# ╔═╡ 1b76ab6c-ffa2-40eb-a6c6-55d7097a5108 -@call_rule typeof(+)(:in1, Marginalisation) (m_out = NormalMeanVariance(3.0, 1.0), m_in2 = NormalMeanVariance(2.0, 1.0)) - -# ╔═╡ 965b886e-d294-11ef-1b10-0319896874cf -md""" -## Code Example: forward and backward messages for the Matrix Multiplication node - -""" - -# ╔═╡ a3e11d46-5a22-4eb6-ba91-7258ba3c667e -TwoColumn( -md""" -In the same way we can also investigate the forward and backward messages for the matrix multiplication ("gain") node -""", -md""" -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-gain-node.png?raw=true) -""" -) - -# ╔═╡ 0efe10d8-1d0e-4a8f-8005-25ee261322b8 -md"Forward message on ``Y``:" - -# ╔═╡ 1be3121d-be18-46a1-9af9-f108a2257c22 -@call_rule typeof(*)(:out, Marginalisation) (m_A = PointMass(4.0), m_in = NormalMeanVariance(1.0, 1.0)) - -# ╔═╡ e5658c95-6cd0-426f-b819-31f9f2c7eaf4 -md"Backward message on ``X``:" - -# ╔═╡ 94ca674e-1a01-424c-8657-6510be7097c3 -@call_rule typeof(*)(:in, Marginalisation) (m_out = NormalMeanVariance(2.0, 1.0), m_A = PointMass(4.0)) - -# ╔═╡ 965c18f8-d294-11ef-2456-b945a46241f4 -md""" -## Example: Sum-Product Algorithm to infer a posterior - -""" - -# ╔═╡ e0add49a-94ac-4247-8554-5a50d4abbebb -TwoColumn( -md""" -Consider a generative model - -```math -p(x,y_1,y_2) = p(x)\,p(y_1|x)\,p(y_2|x) . -``` - -This model expresses the assumption that ``Y_1`` and ``Y_2`` are independent measurements of ``X``. +which, in case ``x_i`` has ``10`` values, **requires a few hundred additions and is therefore computationally (much!) lighter** than executing the whole sum ``\sum_{x_1,\ldots,x_7}f(x_1,x_2,\ldots,x_7)`` + -""", -md""" -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-observations.png?raw=true) """ -) -# ╔═╡ d05277c1-fb9e-4b2b-bcbc-d8be5e63cab5 -TwoColumn( +# ╔═╡ 679d05b5-9fee-46e7-9627-c886c34a612d +keyconcept("", md""" -Assume that we are interested in the posterior for ``X`` after observing ``Y_1= \hat y_1`` and ``Y_2= \hat y_2``. The posterior for ``X`` can be inferred by applying the sum-product algorithm to the following graph: +By naturally exploiting the distributive law, message-passing inference in factor graphs can drastically reduce the computational complexity of Bayesian inference in sparse models. +""") -""", +# ╔═╡ 0afe3cdc-15ed-4d9a-848a-d1977d051866 md""" -![](https://github.com/bmlip/course/blob/v2/assets/figures/ffg-observations-2.png?raw=true) +## Closing-the-Box and Message Passing Interpretation """ -) -# ╔═╡ 965c5f28-d294-11ef-324e-4df3e38b5045 +# ╔═╡ 96571c34-d294-11ef-11ef-29beeb1f96c2 md""" -## Code for Sum-Product Algorithm to infer a posterior - -We'll use RxInfer to build the above graph, and perform sum-product message passing to infer the posterior ``p(x|y_1,y_2)``. We assume ``p(y_1|x)`` and ``p(y_2|x)`` to be Gaussian likelihoods with known variances: +Note that the intermediate result ``\overrightarrow{\mu}_{X_3}(x_3)`` is obtained by multiplying all enclosed factors (``f_a``, ``f_b, f_c``) by the red dashed box, followed by marginalization (summing) over all enclosed variables (``x_1``, ``x_2``), ```math -\begin{align*} - p(y_1\,|\,x) &= \mathcal{N}(y_1\,|\,x, v_{y1}) \\ - p(y_2\,|\,x) &= \mathcal{N}(y_2\,|\,x, v_{y2}) -\end{align*} +\overrightarrow{\mu}_{X_3}(x_3) = \underbrace{\sum_{x_1}\sum_{x_2}}_{\text{enclosed variables}}\underbrace{f_a(x_1)f_b(x_2) f_c(x_1,x_2,x_3) }_{\text{enclosed factors}} ``` -Under this model, the posterior is given by: +This operation is known as **Closing-the-Box**. The result is a new **composite node** that holds the factor ``\overrightarrow{\mu}_{X_3}(x_3)``, and is visually represented by the red dashed box in the factor graph. The composite node ``\overrightarrow{\mu}_{X_3}(x_3)`` depends only on the variable(s) that cross the boundary of the box (in this case ``x_3``) and effectively replaces the internal subgraph contained within the red box. +""" + +# ╔═╡ a7b1f559-3c34-491e-83e7-ba95c8c22c80 +md""" + +The Closing-the-box operation can alternatively be interpreted as **passing a message** from the newly created composite node to the rest of the graph.Closing the red box around ``f_a``, ``f_b`` and ``f_c`` leads to an outgoing message ``\overrightarrow{\mu}_{X_3}(x_3)`` for node ``f_c``, given by ```math -\begin{align*} - p(x\,|\,y_1,y_2) &\propto \overbrace{p(y_1\,|\,x)\,p(y_2\,|\,x)}^{\text{likelihood}}\,\overbrace{p(x)}^{\text{prior}} \\ - &=\mathcal{N}(x\,|\,\hat{y}_1, v_{y1})\, \mathcal{N}(x\,|\,\hat{y}_2, v_{y2}) \, \mathcal{N}(x\,|\,m_x, v_x) -\end{align*} +\begin{align} +\underbrace{\overrightarrow{\mu}_{X_3}(x_3)}_{\substack{ \text{outgoing} \\ \text{message} }} +&= \sum_{x_1}\sum_{x_2} \underbrace{\overrightarrow{\mu}_{X_1}(x_1) \overrightarrow{\mu}_{X_2}(x_2)}_{\substack{\text{incoming} \\ \text{messages}}} \underbrace{f_c(x_1,x_2,x_3)}_{\text{factor}} +\end{align} ``` -so we can validate the answer by solving the Gaussian multiplication manually. - """ -# ╔═╡ d27f7af6-e094-44fa-8ba4-4ad2fa38f8bc -y1_hat = 1.0; y2_hat = 2.0; - -# ╔═╡ 90d62ba0-ca97-43f6-8f5a-0c1086a13f3d +# ╔═╡ 96575dd4-d294-11ef-31d6-b39b4c4bdea1 md""" -Construct the factor graph - -""" - -# ╔═╡ 053e9dde-c088-4f15-9ca6-98b8185a8a11 -@model function my_model(y1,y2) +## Sum-Product Messages - # `x` is the hidden states - x ~ NormalMeanVariance(0.0, 4.0) +This recipe holds generally. For a node ``f(y,x_1,\ldots,x_n)`` with incoming messages ``\overrightarrow{\mu}_{X_1}(x_1)``, ``\overrightarrow{\mu}_{X_2}(x_2)``, ``\ldots``,``\overrightarrow{\mu}_{X_n}(x_n)``, the outgoing message is given by ([Loeliger (2007), pg.1299](https://github.com/bmlip/course/blob/main/assets/files/Loeliger-2007-The-factor-graph-approach-to-model-based-signal-processing.pdf)): - # `y1` and `y2` are "clamped" observations - y1 ~ NormalMeanVariance(x, 1.0) - y2 ~ NormalMeanVariance(x, 2.0) - - return x -end +```math +\underbrace{\overrightarrow{\mu}_{Y}(y)}_{\substack{ \text{outgoing}\\ \text{message}}} = \sum_{x_1,\ldots,x_n} \underbrace{\overrightarrow{\mu}_{X_1}(x_1)\cdots \overrightarrow{\mu}_{X_n}(x_n)}_{\substack{\text{incoming} \\ \text{messages}}} \cdot \underbrace{f(y,x_1,\ldots,x_n)}_{\substack{\text{node}\\ \text{function}}} \tag{SP} +``` -# ╔═╡ 07b09ac1-7fa7-4b62-b130-97315adb6fa7 -result = infer(model=my_model(), data=(y1=y1_hat, y2 = y2_hat,)) +""" -# ╔═╡ defb2149-294b-47a8-99ed-1b3746b275f1 -Text("Sum-product message passing result: p(x|y1,y2) = \n\t𝒩($( - round(mean(result.posteriors[:x]); digits=3) -),$( - round(var(result.posteriors[:x]); digits=3) -))") +# ╔═╡ 5cc2016e-0383-448c-bd33-5b3a687b7436 +TwoColumn( +md""" +Equation (SP) is called a **Sum-Product** message, so named because the computation involves evaluating a sum-of-products. Note that all SP messages in an FFG can be computed from information that is **locally available** at each node. +""", +@htl """ -# ╔═╡ b3656d6c-4717-4fcd-90c6-ae4f4aa5e1be + +""") -# ╔═╡ b15f28ce-c8c1-439b-aeca-74a58d2557e2 +# ╔═╡ f65f5d0e-2583-4b88-b9f2-5fee15257c05 md""" -We calculate mean and variance of `p(x|y1,y2)` manually by multiplying 3 Gaussians (see lesson 4 for details) -""" -# ╔═╡ 86e67c05-068d-4de4-80f3-1a20cc8a43ea -v = 1 / (1/4 + 1/1 + 1/2) -# ╔═╡ fffa27d5-eb68-4dd3-9995-4a53fba6c1e4 -m = v * (0/4 + y1_hat/1.0 + y2_hat/2.0) - -# ╔═╡ 578ec319-337d-4396-bb75-eaf99d95a38d -Text("Manual result: p(x|y1,y2) = \n\t𝒩($(m), $(v))") +If the factor graph for the whole model has no cycles, i.e., the FFG is a tree, then passing SP messages from the terminal nodes to the internal (latent) variables yields exact Bayesian marginals for all hidden variables. This inference method is known as the **Sum-Product** (SP) algorithm. -# ╔═╡ 89da2fc0-a7c8-4a9d-82d9-622a311d010d -md""" -# Code """ -# ╔═╡ 981b08cc-7fb4-4880-8e8a-0b60a5dd72a2 -stable_rand(args...; seed=nothing) = rand(StableRNG(543432 + hash(seed)), args...) - -# ╔═╡ 997235c1-08bd-4dbc-b1bc-cb10a3b83da4 +# ╔═╡ 7009cdc8-892c-499e-b932-b828fa300b6c +keyconcept("", md""" -## Data generation details +For a node ``f(y,x_1,\ldots,x_n)`` with incoming messages ``\overrightarrow{\mu}_{X_1}(x_1)``, ``\overrightarrow{\mu}_{X_2}(x_2)``, ``\ldots``,``\overrightarrow{\mu}_{X_n}(x_n)``, the outgoing message is given by the **sum-product rule**: -We first generate data by a "secret" function ``f`` that is parameterized by weights ``w^*``: +```math +\overrightarrow{\mu}_{Y}(y)= \sum_{x_1,\ldots,x_n} \overrightarrow{\mu}_{X_1}(x_1)\cdots \overrightarrow{\mu}_{X_n}(x_n) \cdot f(y,x_1,\ldots,x_n) +``` +""") +# ╔═╡ 9651f976-b834-4b81-8810-649f0290969d +md""" +# RxInfer: A Toolbox for Automated Bayesian inference """ -# ╔═╡ aec4726a-954e-4e76-aae5-2dd6c979b12d -secret_true_w = [1.0; 2.0; 0.25]; - -# ╔═╡ 96ef3cfb-ca18-46d6-bcac-0122c2c85fba -f(x::Vector)::Real = secret_true_w' * x; - -# ╔═╡ 79a0d02b-368f-4371-854c-cf2cea9328e5 -f([3.0^0, 3.0^1, 3.0^2]) +# ╔═╡ 96587a66-d294-11ef-2c7a-9fd7bea76582 +md""" +## Automating Bayesian Inference by Message Passing -# ╔═╡ cb4427f9-0cb1-4393-b4de-14a4d64cc29c +If the update rules for common node types in a graph have been tabulated, then inference by message passing comes down to executing a set of table-lookup operations, thus creating a completely **automatable Bayesian inference framework**. +In our research lab [BIASlab](http://biaslab.org) at TU Eindhoven, we are developing [RxInfer](http://rxinfer.com), which is a (Julia) toolbox for automating Bayesian inference by message passing in a factor graph. -# ╔═╡ aca1f927-bc3b-48f6-af5c-12ee2ea4a49b -N_bond +In general, a code package that automates Bayesian inference is called a [Probabilistic Programming](https://en.wikipedia.org/wiki/Probabilistic_programming) Language (PPL). RxInfer is a PPL that automates inference through message passing-based inference in a factor graph. -# ╔═╡ 99265e22-e8dc-40fe-989f-0d2a6c72faac -z = stable_rand(Uniform(0, 10), Nsamples; seed=1234) +""" -# ╔═╡ f6fc4fad-70fb-432f-b77d-8e6ad42eef6c +# ╔═╡ 89e2757e-a09f-40c6-8dd7-9b4b4d232e17 md""" -Create the feature vector ``x = [1.0, z, z^2]``: +![](https://github.com/bmlip/course/blob/main/assets/figures/RxInfer-gif.gif?raw=true) """ -# ╔═╡ e20e9048-1271-41c7-97d3-635f320aa365 -x_train = [[1.0, z, z^2] for z in z] +# ╔═╡ c4b5b124-e52a-41fc-b27e-a58181622e5c +md""" +The figure above (a screen recording from the [RxInfer webpage](http://rxinfer.com)) is an animated GIF illustrating how RxInfer operates. The model is represented as a graph in which each node passes messages to its neighbors. When messages meet on an edge, the belief about the variable associated with that edge is updated. +""" -# ╔═╡ 2f86a4cf-2075-45b5-bf2d-2d4d6888461a -x_train +# ╔═╡ 4a10044c-e044-43e1-bd44-847f56019061 +keyconcept( + "", + md""" + “In `RxInfer`, once the model is specified and the observations are provided, Bayesian inference is carried out automatically via the `infer()` function.” + """ +) -# ╔═╡ 3a045b5c-9d87-46a6-a404-85c4bd77dd61 +# ╔═╡ 1f16c8db-f2da-4a72-8a69-743ec857ac3b md""" -Now we can generate the observed ``y`` coordinates in the data set: +Although RxInfer is still under active development, my prediction is that within 5–10 years, RxInfer, or a comparable toolbox, will be able to automate Bayesian inference for virtually any interesting probabilistic model you can conceive. In principle, you will then have all the tools needed to implement the four-step Bayesian ML recipe—model specification, parameter learning, model evaluation, and application, for any (Bayesian) information processing problem. """ -# ╔═╡ ba7a2dbd-f068-4249-bc29-77f2d0804676 -data_noise_σ² = 2.0; +# ╔═╡ fa5bdb1c-4412-48cc-950c-9ed92b4c9f76 +md""" +# Summary +""" -# ╔═╡ 34ebbbe1-2a6b-422b-aeb1-cd2953acddca -# y[i] = w' * x[i] + ϵ -y_train = f.(x_train) + stable_rand(Normal(0, sqrt(data_noise_σ²)), Nsamples; seed=4566) +# ╔═╡ be670693-2036-46cb-8452-a2d0e1bf1172 +keyconceptsummary() -# ╔═╡ 7764541a-c11e-4e12-bbac-f8906cbc5dc6 -scatter(z, y_train; - xlim=(-.2,10.2), - ylim=(-1,51), - label="data", - xlabel=L"z", - ylabel=L"f([1.0, z, z^2]) + \epsilon" -) +# ╔═╡ 89da2fc0-a7c8-4a9d-82d9-622a311d010d +md""" +# Code +""" -# ╔═╡ 22f22f59-c320-4654-b472-b64cc3a001ff -y_train - -# ╔═╡ c03b1140-adce-467a-b953-50ad1bf3bc34 -results = infer( - model = linear_regression( - Nsamples=length(x_train), - Σ=prior_Σ, σ²=data_noise_σ² - ), - data = (y = y_train, x = x_train), - returnvars = (w = KeepLast(),), - iterations = 20, -) +# ╔═╡ cb4427f9-0cb1-4393-b4de-14a4d64cc29c -# ╔═╡ 83a70a4b-b114-4351-8fa2-dd565ebc9916 -convert(MvNormal, results.posteriors[:w]) - -# ╔═╡ 92f7bcfd-00a4-4cb7-a3eb-c1e101fdbcf6 -w_samples = rand(results.posteriors[:w], 10) |> eachcol .|> collect - -# ╔═╡ 965a37e8-d294-11ef-340f-0930b229dd32 -let - plt = scatter( - z, y_train; - xlim=(-.2,10.2), - ylim=(-1,51), - label="data", - xlabel=L"z", - ylabel=L"f([1.0, z, z^2]) + \epsilon" - ) - z_test = collect(0:0.2:12) - x_test = [[1.0; z; z^2] for z in z_test] - for w in w_samples - f_est(x) = w'*x - plot!(plt, z_test, map(f_est, x_test), alpha=0.5, label=nothing); - end - plt -end # ╔═╡ 00000000-0000-0000-0000-000000000001 PLUTO_PROJECT_TOML_CONTENTS = """ @@ -1247,16 +329,11 @@ BmlipTeachingTools = "656a7065-6f73-6c65-7465-6e646e617262" LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" -Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" -RxInfer = "86711068-29c9-4ff7-b620-ae75d7495b3d" -StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3" [compat] BmlipTeachingTools = "~1.3.1" LaTeXStrings = "~1.4.0" Plots = "~1.40.17" -RxInfer = "~4.6.2" -StableRNGs = "~1.0.3" """ # ╔═╡ 00000000-0000-0000-0000-000000000002 @@ -1265,22 +342,7 @@ PLUTO_MANIFEST_TOML_CONTENTS = """ julia_version = "1.12.1" manifest_format = "2.0" -project_hash = "46a2a325723fda574bb6c8846f74c34db6ab166e" - -[[deps.ADTypes]] -git-tree-sha1 = "27cecae79e5cc9935255f90c53bb831cc3c870d7" -uuid = "47edcb42-4c32-4615-8424-f2b9edc5f35b" -version = "1.18.0" - - [deps.ADTypes.extensions] - ADTypesChainRulesCoreExt = "ChainRulesCore" - ADTypesConstructionBaseExt = "ConstructionBase" - ADTypesEnzymeCoreExt = "EnzymeCore" - - [deps.ADTypes.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9" - EnzymeCore = "f151be2c-9106-41f4-ab19-57ee4f262869" +project_hash = "1d1644a5dd8f949f47ef0cb736ed9c8e8eeaf6cc" [[deps.AbstractPlutoDingetjes]] deps = ["Pkg"] @@ -1288,17 +350,6 @@ git-tree-sha1 = "6e1d2a35f2f90a4bc7c2ed98079b2ba09c35b83a" uuid = "6e696c72-6542-2067-7265-42206c756150" version = "1.3.2" -[[deps.Adapt]] -deps = ["LinearAlgebra", "Requires"] -git-tree-sha1 = "7e35fca2bdfba44d797c53dfe63a51fabf39bfc0" -uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" -version = "4.4.0" -weakdeps = ["SparseArrays", "StaticArrays"] - - [deps.Adapt.extensions] - AdaptSparseArraysExt = "SparseArrays" - AdaptStaticArraysExt = "StaticArrays" - [[deps.AliasTables]] deps = ["PtrArrays", "Random"] git-tree-sha1 = "9876e1e164b144ca45e9e3198d0b689cadfed9ff" @@ -1309,56 +360,6 @@ version = "1.1.3" uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" version = "1.1.2" -[[deps.ArnoldiMethod]] -deps = ["LinearAlgebra", "Random", "StaticArrays"] -git-tree-sha1 = "d57bd3762d308bded22c3b82d033bff85f6195c6" -uuid = "ec485272-7323-5ecc-a04f-4719b315124d" -version = "0.4.0" - -[[deps.ArrayInterface]] -deps = ["Adapt", "LinearAlgebra"] -git-tree-sha1 = "d2cd034553ee6ca084edaaf8ed6c9d50fd01555d" -uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" -version = "7.21.0" - - [deps.ArrayInterface.extensions] - ArrayInterfaceBandedMatricesExt = "BandedMatrices" - ArrayInterfaceBlockBandedMatricesExt = "BlockBandedMatrices" - ArrayInterfaceCUDAExt = "CUDA" - ArrayInterfaceCUDSSExt = ["CUDSS", "CUDA"] - ArrayInterfaceChainRulesCoreExt = "ChainRulesCore" - ArrayInterfaceChainRulesExt = "ChainRules" - ArrayInterfaceGPUArraysCoreExt = "GPUArraysCore" - ArrayInterfaceMetalExt = "Metal" - ArrayInterfaceReverseDiffExt = "ReverseDiff" - ArrayInterfaceSparseArraysExt = "SparseArrays" - ArrayInterfaceStaticArraysCoreExt = "StaticArraysCore" - ArrayInterfaceTrackerExt = "Tracker" - - [deps.ArrayInterface.weakdeps] - BandedMatrices = "aae01518-5342-5314-be14-df237901396f" - BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" - CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" - CUDSS = "45b445bb-4962-46a0-9369-b4df9d0f772e" - ChainRules = "082447d4-558c-5d27-93f4-14fc19e9eca2" - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" - Metal = "dde4c033-4e86-420c-a63e-0dd931031962" - ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" - SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" - StaticArraysCore = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" - Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" - -[[deps.ArrayLayouts]] -deps = ["FillArrays", "LinearAlgebra", "StaticArrays"] -git-tree-sha1 = "355ab2d61069927d4247cd69ad0e1f140b31e30d" -uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" -version = "1.12.0" -weakdeps = ["SparseArrays"] - - [deps.ArrayLayouts.extensions] - ArrayLayoutsSparseArraysExt = "SparseArrays" - [[deps.Artifacts]] uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" version = "1.11.0" @@ -1367,41 +368,11 @@ version = "1.11.0" uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" version = "1.11.0" -[[deps.BayesBase]] -deps = ["Distributions", "DomainSets", "LinearAlgebra", "Random", "SpecialFunctions", "StaticArrays", "Statistics", "StatsAPI", "StatsBase", "StatsFuns", "TinyHugeNumbers"] -git-tree-sha1 = "5b723bf6b1081cab4d263e425be097224e0f434f" -uuid = "b4ee3484-f114-42fe-b91c-797d54a0c67e" -version = "1.5.8" -weakdeps = ["FastCholesky"] - - [deps.BayesBase.extensions] - FastCholeskyExt = "FastCholesky" - [[deps.BitFlags]] git-tree-sha1 = "0691e34b3bb8be9307330f88d1a3c3f25466c24d" uuid = "d1d4a3ce-64b1-5f1a-9ba4-7e7e69966f35" version = "0.1.9" -[[deps.BitSetTuples]] -deps = ["TupleTools"] -git-tree-sha1 = "aa19428fb6ad21db22f8568f068de4f443d3bacc" -uuid = "0f2f92aa-23a3-4d05-b791-88071d064721" -version = "1.1.5" - -[[deps.BlockArrays]] -deps = ["ArrayLayouts", "FillArrays", "LinearAlgebra"] -git-tree-sha1 = "79e651aa489a7879107d66e3d1948e9aa1b4055e" -uuid = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" -version = "1.7.2" - - [deps.BlockArrays.extensions] - BlockArraysAdaptExt = "Adapt" - BlockArraysBandedMatricesExt = "BandedMatrices" - - [deps.BlockArrays.weakdeps] - Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" - BandedMatrices = "aae01518-5342-5314-be14-df237901396f" - [[deps.BmlipTeachingTools]] deps = ["HypertextLiteral", "InteractiveUtils", "Markdown", "PlutoTeachingTools", "PlutoUI", "Reexport"] git-tree-sha1 = "806eadb642467b05f9d930f0d127f1e6fa5130f0" @@ -1420,23 +391,6 @@ git-tree-sha1 = "fde3bf89aead2e723284a8ff9cdf5b551ed700e8" uuid = "83423d85-b0ee-5818-9007-b63ccbeb887a" version = "1.18.5+0" -[[deps.ChunkCodecCore]] -git-tree-sha1 = "51f4c10ee01bda57371e977931de39ee0f0cdb3e" -uuid = "0b6fb165-00bc-4d37-ab8b-79f91016dbe1" -version = "1.0.0" - -[[deps.ChunkCodecLibZlib]] -deps = ["ChunkCodecCore", "Zlib_jll"] -git-tree-sha1 = "cee8104904c53d39eb94fd06cbe60cb5acde7177" -uuid = "4c0bbee4-addc-4d73-81a0-b6caacae83c8" -version = "1.0.0" - -[[deps.ChunkCodecLibZstd]] -deps = ["ChunkCodecCore", "Zstd_jll"] -git-tree-sha1 = "34d9873079e4cb3d0c62926a225136824677073f" -uuid = "55437552-ac27-4d47-9aa3-63184e8fd398" -version = "1.0.0" - [[deps.CodecZlib]] deps = ["TranscodingStreams", "Zlib_jll"] git-tree-sha1 = "962834c22b66e32aa10f7611c08c8ca4e20749a9" @@ -1464,33 +418,19 @@ deps = ["ColorTypes", "FixedPointNumbers", "LinearAlgebra", "Requires", "Statist git-tree-sha1 = "8b3b6f87ce8f65a2b4f857528fd8d70086cd72b1" uuid = "c3611d14-8923-5661-9e6a-0046d554d3a4" version = "0.11.0" -weakdeps = ["SpecialFunctions"] [deps.ColorVectorSpace.extensions] SpecialFunctionsExt = "SpecialFunctions" + [deps.ColorVectorSpace.weakdeps] + SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" + [[deps.Colors]] deps = ["ColorTypes", "FixedPointNumbers", "Reexport"] git-tree-sha1 = "37ea44092930b1811e666c3bc38065d7d87fcc74" uuid = "5ae59095-9a9b-59fe-a467-6f913c188581" version = "0.13.1" -[[deps.Combinatorics]] -git-tree-sha1 = "8010b6bb3388abe68d95743dcbea77650bb2eddf" -uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" -version = "1.0.3" - -[[deps.CommonSubexpressions]] -deps = ["MacroTools"] -git-tree-sha1 = "cda2cfaebb4be89c9084adaca7dd7333369715c5" -uuid = "bbf7d656-a473-5ed7-a52c-81e309532950" -version = "0.3.1" - -[[deps.CommonWorldInvalidations]] -git-tree-sha1 = "ae52d1c52048455e85a387fbee9be553ec2b68d0" -uuid = "f70d9fcc-98c5-4d4a-abd7-e4cdeebd8ca8" -version = "1.0.0" - [[deps.Compat]] deps = ["TOML", "UUIDs"] git-tree-sha1 = "9d8a54ce4b17aa5bdce0ea5c34bc5e7c340d16ad" @@ -1506,28 +446,12 @@ deps = ["Artifacts", "Libdl"] uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" version = "1.3.0+1" -[[deps.CompositeTypes]] -git-tree-sha1 = "bce26c3dab336582805503bed209faab1c279768" -uuid = "b152e2b5-7a66-4b01-a709-34e65c35f657" -version = "0.1.4" - [[deps.ConcurrentUtilities]] deps = ["Serialization", "Sockets"] git-tree-sha1 = "d9d26935a0bcffc87d2613ce14c527c99fc543fd" uuid = "f0e56b4a-5159-44fe-b623-3e5288b988bb" version = "2.5.0" -[[deps.ConstructionBase]] -git-tree-sha1 = "b4b092499347b18a015186eae3042f72267106cb" -uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" -version = "1.6.0" -weakdeps = ["IntervalSets", "LinearAlgebra", "StaticArrays"] - - [deps.ConstructionBase.extensions] - ConstructionBaseIntervalSetsExt = "IntervalSets" - ConstructionBaseLinearAlgebraExt = "LinearAlgebra" - ConstructionBaseStaticArraysExt = "StaticArrays" - [[deps.Contour]] git-tree-sha1 = "439e35b0b36e2e5881738abc8857bd92ad6ff9a8" uuid = "d38c429a-6771-53c6-b99e-75d170b6e991" @@ -1561,130 +485,16 @@ git-tree-sha1 = "9e2f36d3c96a820c678f2f1f1782582fcf685bae" uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab" version = "1.9.1" -[[deps.Dictionaries]] -deps = ["Indexing", "Random", "Serialization"] -git-tree-sha1 = "a86af9c4c4f33e16a2b2ff43c2113b2f390081fa" -uuid = "85a47980-9c8c-11e8-2b9f-f7ca1fa99fb4" -version = "0.4.5" - -[[deps.DiffResults]] -deps = ["StaticArraysCore"] -git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" -uuid = "163ba53b-c6d8-5494-b064-1a9d43ac40c5" -version = "1.1.0" - -[[deps.DiffRules]] -deps = ["IrrationalConstants", "LogExpFunctions", "NaNMath", "Random", "SpecialFunctions"] -git-tree-sha1 = "23163d55f885173722d1e4cf0f6110cdbaf7e272" -uuid = "b552c78f-8df3-52c6-915a-8e097449b14b" -version = "1.15.1" - -[[deps.DifferentiationInterface]] -deps = ["ADTypes", "LinearAlgebra"] -git-tree-sha1 = "529bebbc74b36a4cfea09dd2aecb1288cd713a6d" -uuid = "a0c0ee7d-e4b9-4e03-894e-1c5f64a51d63" -version = "0.7.9" - - [deps.DifferentiationInterface.extensions] - DifferentiationInterfaceChainRulesCoreExt = "ChainRulesCore" - DifferentiationInterfaceDiffractorExt = "Diffractor" - DifferentiationInterfaceEnzymeExt = ["EnzymeCore", "Enzyme"] - DifferentiationInterfaceFastDifferentiationExt = "FastDifferentiation" - DifferentiationInterfaceFiniteDiffExt = "FiniteDiff" - DifferentiationInterfaceFiniteDifferencesExt = "FiniteDifferences" - DifferentiationInterfaceForwardDiffExt = ["ForwardDiff", "DiffResults"] - DifferentiationInterfaceGPUArraysCoreExt = "GPUArraysCore" - DifferentiationInterfaceGTPSAExt = "GTPSA" - DifferentiationInterfaceMooncakeExt = "Mooncake" - DifferentiationInterfacePolyesterForwardDiffExt = ["PolyesterForwardDiff", "ForwardDiff", "DiffResults"] - DifferentiationInterfaceReverseDiffExt = ["ReverseDiff", "DiffResults"] - DifferentiationInterfaceSparseArraysExt = "SparseArrays" - DifferentiationInterfaceSparseConnectivityTracerExt = "SparseConnectivityTracer" - DifferentiationInterfaceSparseMatrixColoringsExt = "SparseMatrixColorings" - DifferentiationInterfaceStaticArraysExt = "StaticArrays" - DifferentiationInterfaceSymbolicsExt = "Symbolics" - DifferentiationInterfaceTrackerExt = "Tracker" - DifferentiationInterfaceZygoteExt = ["Zygote", "ForwardDiff"] - - [deps.DifferentiationInterface.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - DiffResults = "163ba53b-c6d8-5494-b064-1a9d43ac40c5" - Diffractor = "9f5e2b26-1114-432f-b630-d3fe2085c51c" - Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" - EnzymeCore = "f151be2c-9106-41f4-ab19-57ee4f262869" - FastDifferentiation = "eb9bf01b-bf85-4b60-bf87-ee5de06c00be" - FiniteDiff = "6a86dc24-6348-571c-b903-95158fe2bd41" - FiniteDifferences = "26cc04aa-876d-5657-8c51-4c34ba976000" - ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" - GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" - GTPSA = "b27dd330-f138-47c5-815b-40db9dd9b6e8" - Mooncake = "da2b9cff-9c12-43a0-ae48-6db2b0edb7d6" - PolyesterForwardDiff = "98d1487c-24ca-40b6-b7ab-df2af84e126b" - ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" - SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" - SparseConnectivityTracer = "9f842d2f-2579-4b1d-911e-f412cf18a3f5" - SparseMatrixColorings = "0a514795-09f3-496d-8182-132a7b665d35" - StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" - Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" - Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" - -[[deps.Distributed]] -deps = ["Random", "Serialization", "Sockets"] -uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" -version = "1.11.0" - -[[deps.Distributions]] -deps = ["AliasTables", "FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SpecialFunctions", "Statistics", "StatsAPI", "StatsBase", "StatsFuns"] -git-tree-sha1 = "3bc002af51045ca3b47d2e1787d6ce02e68b943a" -uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" -version = "0.25.122" - - [deps.Distributions.extensions] - DistributionsChainRulesCoreExt = "ChainRulesCore" - DistributionsDensityInterfaceExt = "DensityInterface" - DistributionsTestExt = "Test" - - [deps.Distributions.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d" - Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" - [[deps.DocStringExtensions]] git-tree-sha1 = "7442a5dfe1ebb773c29cc2962a8980f47221d76c" uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" version = "0.9.5" -[[deps.DomainIntegrals]] -deps = ["CompositeTypes", "DomainSets", "FastGaussQuadrature", "GaussQuadrature", "HCubature", "IntervalSets", "LinearAlgebra", "QuadGK", "SpecialFunctions", "StaticArrays"] -git-tree-sha1 = "934bf806ef2948114243f25e84a3ddf775d0f1a6" -uuid = "cc6bae93-f070-4015-88fd-838f9505a86c" -version = "0.5.2" - -[[deps.DomainSets]] -deps = ["CompositeTypes", "IntervalSets", "LinearAlgebra", "StaticArrays"] -git-tree-sha1 = "c249d86e97a7e8398ce2068dce4c078a1c3464de" -uuid = "5b8099bc-c8ec-5219-889f-1d9e522a28bf" -version = "0.7.16" - - [deps.DomainSets.extensions] - DomainSetsMakieExt = "Makie" - DomainSetsRandomExt = "Random" - - [deps.DomainSets.weakdeps] - Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" - Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" - [[deps.Downloads]] deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" version = "1.6.0" -[[deps.EnumX]] -git-tree-sha1 = "bddad79635af6aec424f53ed8aad5d7555dc6f00" -uuid = "4e289a0a-7415-4d19-859d-a7e5c4648b56" -version = "1.0.5" - [[deps.EpollShim_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] git-tree-sha1 = "8a4be429317c42cfae6a7fc03c31bad1970c310d" @@ -1703,12 +513,6 @@ git-tree-sha1 = "27af30de8b5445644e8ffe3bcb0d72049c089cf1" uuid = "2e619515-83b5-522b-bb60-26c02a35a201" version = "2.7.3+0" -[[deps.ExponentialFamily]] -deps = ["BayesBase", "BlockArrays", "Distributions", "DomainSets", "FastCholesky", "FillArrays", "ForwardDiff", "HCubature", "HypergeometricFunctions", "IntervalSets", "IrrationalConstants", "LinearAlgebra", "LogExpFunctions", "PositiveFactorizations", "Random", "SparseArrays", "SpecialFunctions", "StaticArrays", "StatsBase", "StatsFuns", "TinyHugeNumbers"] -git-tree-sha1 = "8351e116e111c97ad57718851da00ae5a5f92e0c" -uuid = "62312e5e-252a-4322-ace9-a5f4bf9b357b" -version = "2.1.1" - [[deps.FFMPEG]] deps = ["FFMPEG_jll"] git-tree-sha1 = "83dc665d0312b41367b7263e8a4d172eac1897f4" @@ -1721,72 +525,10 @@ git-tree-sha1 = "3a948313e7a41eb1db7a1e733e6335f17b4ab3c4" uuid = "b22a6f82-2f65-5046-a5b2-351ab43fb4e5" version = "7.1.1+0" -[[deps.FastCholesky]] -deps = ["LinearAlgebra", "PositiveFactorizations"] -git-tree-sha1 = "1c0a81e006e40e9fcbd5f6f6cb42ac2700f86889" -uuid = "2d5283b6-8564-42b6-bb00-83ed8e915756" -version = "1.4.3" -weakdeps = ["StaticArraysCore"] - - [deps.FastCholesky.extensions] - StaticArraysCoreExt = "StaticArraysCore" - -[[deps.FastGaussQuadrature]] -deps = ["LinearAlgebra", "SpecialFunctions", "StaticArrays"] -git-tree-sha1 = "0044e9f5e49a57e88205e8f30ab73928b05fe5b6" -uuid = "442a2c76-b920-505d-bb47-c5924d526838" -version = "1.1.0" - -[[deps.FileIO]] -deps = ["Pkg", "Requires", "UUIDs"] -git-tree-sha1 = "d60eb76f37d7e5a40cc2e7c36974d864b82dc802" -uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549" -version = "1.17.1" -weakdeps = ["HTTP"] - - [deps.FileIO.extensions] - HTTPExt = "HTTP" - [[deps.FileWatching]] uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" version = "1.11.0" -[[deps.FillArrays]] -deps = ["LinearAlgebra"] -git-tree-sha1 = "173e4d8f14230a7523ae11b9a3fa9edb3e0efd78" -uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" -version = "1.14.0" -weakdeps = ["PDMats", "SparseArrays", "Statistics"] - - [deps.FillArrays.extensions] - FillArraysPDMatsExt = "PDMats" - FillArraysSparseArraysExt = "SparseArrays" - FillArraysStatisticsExt = "Statistics" - -[[deps.FiniteDiff]] -deps = ["ArrayInterface", "LinearAlgebra", "Setfield"] -git-tree-sha1 = "9340ca07ca27093ff68418b7558ca37b05f8aeb1" -uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" -version = "2.29.0" - - [deps.FiniteDiff.extensions] - FiniteDiffBandedMatricesExt = "BandedMatrices" - FiniteDiffBlockBandedMatricesExt = "BlockBandedMatrices" - FiniteDiffSparseArraysExt = "SparseArrays" - FiniteDiffStaticArraysExt = "StaticArrays" - - [deps.FiniteDiff.weakdeps] - BandedMatrices = "aae01518-5342-5314-be14-df237901396f" - BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" - SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" - StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - -[[deps.FixedArguments]] -deps = ["TupleTools"] -git-tree-sha1 = "befa1ad59c77643dec6fc20d71fd6f5c3afcdadd" -uuid = "4130a065-6d82-41fe-881e-7a5c65156f7d" -version = "0.1.1" - [[deps.FixedPointNumbers]] deps = ["Statistics"] git-tree-sha1 = "05882d6995ae5c12bb5f36dd2ed3f61c98cbb172" @@ -1804,16 +546,6 @@ git-tree-sha1 = "9c68794ef81b08086aeb32eeaf33531668d5f5fc" uuid = "1fa38f19-a742-5d3f-a2b9-30dd87b9d5f8" version = "1.3.7" -[[deps.ForwardDiff]] -deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions"] -git-tree-sha1 = "ba6ce081425d0afb2bedd00d9884464f764a9225" -uuid = "f6369f11-7733-5829-9624-2563aa707210" -version = "1.2.2" -weakdeps = ["StaticArrays"] - - [deps.ForwardDiff.extensions] - ForwardDiffStaticArraysExt = "StaticArrays" - [[deps.FreeType2_jll]] deps = ["Artifacts", "Bzip2_jll", "JLLWrappers", "Libdl", "Zlib_jll"] git-tree-sha1 = "2c5512e11c791d1baed2049c5652441b28fc6a31" @@ -1826,11 +558,6 @@ git-tree-sha1 = "7a214fdac5ed5f59a22c2d9a885a16da1c74bbc7" uuid = "559328eb-81f9-559d-9380-de523a88c83c" version = "1.0.17+0" -[[deps.Future]] -deps = ["Random"] -uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" -version = "1.11.0" - [[deps.GLFW_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Libglvnd_jll", "Xorg_libXcursor_jll", "Xorg_libXi_jll", "Xorg_libXinerama_jll", "Xorg_libXrandr_jll", "libdecor_jll", "xkbcommon_jll"] git-tree-sha1 = "fcb0584ff34e25155876418979d4c8971243bb89" @@ -1849,12 +576,6 @@ git-tree-sha1 = "27299071cc29e409488ada41ec7643e0ab19091f" uuid = "d2c73de3-f751-5644-a686-071e5b155ba9" version = "0.73.17+0" -[[deps.GaussQuadrature]] -deps = ["SpecialFunctions"] -git-tree-sha1 = "eb6f1f48aa994f3018cbd029a17863c6535a266d" -uuid = "d54b0c1a-921d-58e0-8e36-89d8069c0969" -version = "0.5.8" - [[deps.GettextRuntime_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Libiconv_jll"] git-tree-sha1 = "45288942190db7c5f760f59c04495064eedf9340" @@ -1873,46 +594,17 @@ git-tree-sha1 = "50c11ffab2a3d50192a228c313f05b5b5dc5acb2" uuid = "7746bdde-850d-59dc-9ae8-88ece973131d" version = "2.86.0+0" -[[deps.GraphPPL]] -deps = ["BitSetTuples", "DataStructures", "Dictionaries", "MacroTools", "MetaGraphsNext", "NamedTupleTools", "Static", "StaticArrays", "TupleTools", "Unrolled"] -git-tree-sha1 = "db4aece54ddddaa9e8d2880eb7cfc6f29bd1a650" -uuid = "b3f8163a-e979-4e85-b43e-1f63d8c8b42c" -version = "4.6.5" - - [deps.GraphPPL.extensions] - GraphPPLDistributionsExt = "Distributions" - GraphPPLGraphVizExt = "GraphViz" - GraphPPLPlottingExt = ["Cairo", "GraphPlot"] - - [deps.GraphPPL.weakdeps] - Cairo = "159f3aea-2a34-519c-b102-8c37f9878175" - Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" - GraphPlot = "a2cc645c-3eea-5389-862e-a155d0052231" - GraphViz = "f526b714-d49f-11e8-06ff-31ed36ee7ee0" - [[deps.Graphite2_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] git-tree-sha1 = "8a6dbda1fd736d60cc477d99f2e7a042acfa46e8" uuid = "3b182d85-2403-5c21-9c21-1e1f0cc25472" version = "1.3.15+0" -[[deps.Graphs]] -deps = ["ArnoldiMethod", "DataStructures", "Distributed", "Inflate", "LinearAlgebra", "Random", "SharedArrays", "SimpleTraits", "SparseArrays", "Statistics"] -git-tree-sha1 = "7a98c6502f4632dbe9fb1973a4244eaa3324e84d" -uuid = "86223c79-3864-5bf0-83f7-82e725a168b6" -version = "1.13.1" - [[deps.Grisu]] git-tree-sha1 = "53bb909d1151e57e2484c3d1b53e19552b887fb2" uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe" version = "1.0.2" -[[deps.HCubature]] -deps = ["Combinatorics", "DataStructures", "LinearAlgebra", "QuadGK", "StaticArrays"] -git-tree-sha1 = "19ef9f0cb324eed957b7fe7257ac84e8ed8a48ec" -uuid = "19dc6840-f33b-545b-b366-655c7e3ffd49" -version = "1.7.0" - [[deps.HTTP]] deps = ["Base64", "CodecZlib", "ConcurrentUtilities", "Dates", "ExceptionUnwrapping", "Logging", "LoggingExtras", "MbedTLS", "NetworkOptions", "OpenSSL", "PrecompileTools", "Random", "SimpleBufferStream", "Sockets", "URIs", "UUIDs"] git-tree-sha1 = "5e6fe50ae7f23d171f44e311c2960294aaa0beb5" @@ -1925,17 +617,6 @@ git-tree-sha1 = "f923f9a774fcf3f5cb761bfa43aeadd689714813" uuid = "2e76f6c2-a576-52d4-95c1-20adfe4de566" version = "8.5.1+0" -[[deps.HashArrayMappedTries]] -git-tree-sha1 = "2eaa69a7cab70a52b9687c8bf950a5a93ec895ae" -uuid = "076d061b-32b6-4027-95e0-9a2c6f6d7e74" -version = "0.2.0" - -[[deps.HypergeometricFunctions]] -deps = ["LinearAlgebra", "OpenLibm_jll", "SpecialFunctions"] -git-tree-sha1 = "68c173f4f449de5b438ee67ed0c9c748dc31a2ec" -uuid = "34004b35-14d8-5ef3-9330-4cdb6864b03a" -version = "0.3.28" - [[deps.Hyperscript]] deps = ["Test"] git-tree-sha1 = "179267cfa5e712760cd43dcae385d7ea90cc25a4" @@ -1954,52 +635,16 @@ git-tree-sha1 = "b6d6bfdd7ce25b0f9b2f6b3dd56b2673a66c8770" uuid = "b5f81e59-6552-4d32-b1f0-c071b021bf89" version = "0.2.5" -[[deps.IfElse]] -git-tree-sha1 = "debdd00ffef04665ccbb3e150747a77560e8fad1" -uuid = "615f187c-cbe4-4ef1-ba3b-2fcf58d6d173" -version = "0.1.1" - -[[deps.Indexing]] -git-tree-sha1 = "ce1566720fd6b19ff3411404d4b977acd4814f9f" -uuid = "313cdc1a-70c2-5d6a-ae34-0150d3930a38" -version = "1.1.1" - -[[deps.Inflate]] -git-tree-sha1 = "d1b1b796e47d94588b3757fe84fbf65a5ec4a80d" -uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9" -version = "0.1.5" - [[deps.InteractiveUtils]] deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" version = "1.11.0" -[[deps.IntervalSets]] -git-tree-sha1 = "5fbb102dcb8b1a858111ae81d56682376130517d" -uuid = "8197267c-284f-5f27-9208-e0e47529a953" -version = "0.7.11" -weakdeps = ["Random", "RecipesBase", "Statistics"] - - [deps.IntervalSets.extensions] - IntervalSetsRandomExt = "Random" - IntervalSetsRecipesBaseExt = "RecipesBase" - IntervalSetsStatisticsExt = "Statistics" - [[deps.IrrationalConstants]] git-tree-sha1 = "b2d91fe939cae05960e760110b328288867b5758" uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" version = "0.2.6" -[[deps.JLD2]] -deps = ["ChunkCodecLibZlib", "ChunkCodecLibZstd", "FileIO", "MacroTools", "Mmap", "OrderedCollections", "PrecompileTools", "ScopedValues"] -git-tree-sha1 = "da2e9b4d1abbebdcca0aa68afa0aa272102baad7" -uuid = "033835bb-8acc-5ee8-8aae-3f567f8a3819" -version = "0.6.2" -weakdeps = ["UnPack"] - - [deps.JLD2.extensions] - UnPackExt = "UnPack" - [[deps.JLFzf]] deps = ["REPL", "Random", "fzf_jll"] git-tree-sha1 = "82f7acdc599b65e0f8ccd270ffa1467c21cb647b" @@ -2076,24 +721,6 @@ version = "0.16.10" SymEngine = "123dc426-2d89-5057-bbad-38513e3affd8" tectonic_jll = "d7dd28d6-a5e6-559c-9131-7eb760cdacc5" -[[deps.LazyArrays]] -deps = ["ArrayLayouts", "FillArrays", "LinearAlgebra", "MacroTools", "SparseArrays"] -git-tree-sha1 = "79ee64f6ba0a5a49930f51c86f60d7526b5e12c8" -uuid = "5078a376-72f3-5289-bfd5-ec5146d43c02" -version = "2.8.0" - - [deps.LazyArrays.extensions] - LazyArraysBandedMatricesExt = "BandedMatrices" - LazyArraysBlockArraysExt = "BlockArrays" - LazyArraysBlockBandedMatricesExt = "BlockBandedMatrices" - LazyArraysStaticArraysExt = "StaticArrays" - - [deps.LazyArrays.weakdeps] - BandedMatrices = "aae01518-5342-5314-be14-df237901396f" - BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" - BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" - StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - [[deps.LibCURL]] deps = ["LibCURL_jll", "MozillaCACerts_jll"] uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" @@ -2159,12 +786,6 @@ git-tree-sha1 = "2a7a12fc0a4e7fb773450d17975322aa77142106" uuid = "38a345b3-de98-5d2b-a5d3-14cd9215e700" version = "2.41.2+0" -[[deps.LineSearches]] -deps = ["LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "Printf"] -git-tree-sha1 = "4adee99b7262ad2a1a4bbbc59d993d24e55ea96f" -uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" -version = "7.4.0" - [[deps.LinearAlgebra]] deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" @@ -2211,12 +832,6 @@ deps = ["Base64", "JuliaSyntaxHighlighting", "StyledStrings"] uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" version = "1.11.0" -[[deps.MatrixCorrectionTools]] -deps = ["LinearAlgebra"] -git-tree-sha1 = "73f93b21eae5714c282396bfae9d9f13d6ad04b6" -uuid = "41f81499-25de-46de-b591-c3cfc21e9eaf" -version = "1.2.0" - [[deps.MbedTLS]] deps = ["Dates", "MbedTLS_jll", "MozillaCACerts_jll", "NetworkOptions", "Random", "Sockets"] git-tree-sha1 = "c067a280ddc25f196b5e7df3877c6b226d390aaf" @@ -2234,12 +849,6 @@ git-tree-sha1 = "c13304c81eec1ed3af7fc20e75fb6b26092a1102" uuid = "442fdcdd-2543-5da2-b0f3-8c86c306513e" version = "0.3.2" -[[deps.MetaGraphsNext]] -deps = ["Graphs", "JLD2", "SimpleTraits"] -git-tree-sha1 = "c3f7e597f1cf5fe04e68e7907af47f055cad211c" -uuid = "fa8bd995-216d-47f1-8a91-f3b68fbeb377" -version = "0.7.4" - [[deps.Missings]] deps = ["DataAPI"] git-tree-sha1 = "ec4f7fbeab05d7747bdf98eb74d130a2a2ed298d" @@ -2254,23 +863,12 @@ version = "1.11.0" uuid = "14a3606d-f60d-562e-9121-12d972cd8159" version = "2025.5.20" -[[deps.NLSolversBase]] -deps = ["ADTypes", "DifferentiationInterface", "Distributed", "FiniteDiff", "ForwardDiff"] -git-tree-sha1 = "25a6638571a902ecfb1ae2a18fc1575f86b1d4df" -uuid = "d41bc354-129a-5804-8e4c-c37616107c6c" -version = "7.10.0" - [[deps.NaNMath]] deps = ["OpenLibm_jll"] git-tree-sha1 = "9b8215b1ee9e78a293f99797cd31375471b2bcae" uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" version = "1.1.3" -[[deps.NamedTupleTools]] -git-tree-sha1 = "90914795fc59df44120fe3fff6742bb0d7adb1d0" -uuid = "d9ec5142-1e00-5aa0-9d6a-321866360f50" -version = "0.14.3" - [[deps.NetworkOptions]] uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" version = "1.3.0" @@ -2302,24 +900,6 @@ deps = ["Artifacts", "Libdl"] uuid = "458c3c95-2e84-50aa-8efc-19380b2a3a95" version = "3.5.1+0" -[[deps.OpenSpecFun_jll]] -deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl"] -git-tree-sha1 = "1346c9208249809840c91b26703912dff463d335" -uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" -version = "0.5.6+0" - -[[deps.Optim]] -deps = ["Compat", "EnumX", "FillArrays", "ForwardDiff", "LineSearches", "LinearAlgebra", "NLSolversBase", "NaNMath", "PositiveFactorizations", "Printf", "SparseArrays", "StatsBase"] -git-tree-sha1 = "61942645c38dd2b5b78e2082c9b51ab315315d10" -uuid = "429524aa-4258-5aef-a3af-852621145aeb" -version = "1.13.2" - - [deps.Optim.extensions] - OptimMOIExt = "MathOptInterface" - - [deps.Optim.weakdeps] - MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" - [[deps.Opus_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] git-tree-sha1 = "c392fc5dd032381919e3b22dd32d6443760ce7ea" @@ -2336,28 +916,12 @@ deps = ["Artifacts", "Libdl"] uuid = "efcefdf7-47ab-520b-bdef-62a2eaa19f15" version = "10.44.0+1" -[[deps.PDMats]] -deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "d922b4d80d1e12c658da7785e754f4796cc1d60d" -uuid = "90014a1f-27ba-587c-ab20-58faa44d9150" -version = "0.11.36" -weakdeps = ["StatsBase"] - - [deps.PDMats.extensions] - StatsBaseExt = "StatsBase" - [[deps.Pango_jll]] deps = ["Artifacts", "Cairo_jll", "Fontconfig_jll", "FreeType2_jll", "FriBidi_jll", "Glib_jll", "HarfBuzz_jll", "JLLWrappers", "Libdl"] git-tree-sha1 = "1f7f9bbd5f7a2e5a9f7d96e51c9754454ea7f60b" uuid = "36c8627f-9965-5494-a995-c6b170f724f3" version = "1.56.4+0" -[[deps.Parameters]] -deps = ["OrderedCollections", "UnPack"] -git-tree-sha1 = "34c0e9ad262e5f7fc75b10a9952ca7692cfc5fbe" -uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a" -version = "0.12.3" - [[deps.Parsers]] deps = ["Dates", "PrecompileTools", "UUIDs"] git-tree-sha1 = "7d2f8f21da5db6a806faf7b9b292296da42b2810" @@ -2423,18 +987,6 @@ git-tree-sha1 = "3faff84e6f97a7f18e0dd24373daa229fd358db5" uuid = "7f904dfe-b85e-4ff6-b463-dae2292396a8" version = "0.7.73" -[[deps.PolyaGammaHybridSamplers]] -deps = ["Distributions", "Random", "SpecialFunctions", "StatsFuns"] -git-tree-sha1 = "9f6139650ff57f9d8528cd809ebc604c7e9738b1" -uuid = "c636ee4f-4591-4d8c-9fae-2dea21daa433" -version = "1.2.6" - -[[deps.PositiveFactorizations]] -deps = ["LinearAlgebra"] -git-tree-sha1 = "17275485f373e6673f7e7f97051f703ed5b15b20" -uuid = "85a6dd25-e78a-55b7-8502-1745935b8125" -version = "0.2.4" - [[deps.PrecompileTools]] deps = ["Preferences"] git-tree-sha1 = "07a921781cab75691315adc645096ed5e370cb77" @@ -2452,12 +1004,6 @@ deps = ["Unicode"] uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7" version = "1.11.0" -[[deps.ProgressMeter]] -deps = ["Distributed", "Printf"] -git-tree-sha1 = "fbb92c6c56b34e1a2c4c36058f68f332bec840e7" -uuid = "92933f4c-e287-5a05-a399-4b506db050ca" -version = "1.11.0" - [[deps.PtrArrays]] git-tree-sha1 = "1d36ef11a9aaf1e8b74dacc6a731dd1de8fd493d" uuid = "43287f4e-b6f4-7ad1-bb20-aadabca52c3d" @@ -2487,18 +1033,6 @@ git-tree-sha1 = "8f528b0851b5b7025032818eb5abbeb8a736f853" uuid = "e99dba38-086e-5de3-a5b1-6e4c66e897c3" version = "6.8.2+2" -[[deps.QuadGK]] -deps = ["DataStructures", "LinearAlgebra"] -git-tree-sha1 = "9da16da70037ba9d701192e27befedefb91ec284" -uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" -version = "2.11.2" - - [deps.QuadGK.extensions] - QuadGKEnzymeExt = "Enzyme" - - [deps.QuadGK.weakdeps] - Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" - [[deps.REPL]] deps = ["InteractiveUtils", "JuliaSyntaxHighlighting", "Markdown", "Sockets", "StyledStrings", "Unicode"] uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" @@ -2509,22 +1043,6 @@ deps = ["SHA"] uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" version = "1.11.0" -[[deps.ReactiveMP]] -deps = ["BayesBase", "DataStructures", "DiffResults", "Distributions", "DomainIntegrals", "DomainSets", "ExponentialFamily", "FastCholesky", "FastGaussQuadrature", "FixedArguments", "ForwardDiff", "HCubature", "LazyArrays", "LinearAlgebra", "MacroTools", "MatrixCorrectionTools", "Optim", "PolyaGammaHybridSamplers", "PositiveFactorizations", "Random", "Rocket", "SpecialFunctions", "StaticArrays", "StatsBase", "StatsFuns", "TinyHugeNumbers", "Tullio", "TupleTools", "Unrolled"] -git-tree-sha1 = "dcd2f85ba9f2f7be1b1b31708db889d078b161f6" -uuid = "a194aa59-28ba-4574-a09c-4a745416d6e3" -version = "5.6.2" - - [deps.ReactiveMP.extensions] - ReactiveMPOptimisersExt = "Optimisers" - ReactiveMPProjectionExt = "ExponentialFamilyProjection" - ReactiveMPRequiresExt = "Requires" - - [deps.ReactiveMP.weakdeps] - ExponentialFamilyProjection = "17f509fa-9a96-44ba-99b2-1c5f01f0931b" - Optimisers = "3bd65402-5787-11e9-1adc-39752487f4e2" - Requires = "ae029012-a4dd-5104-9daa-d747884805df" - [[deps.RecipesBase]] deps = ["PrecompileTools"] git-tree-sha1 = "5c3d09cc4f31f5fc6af001c250bf1278733100ff" @@ -2554,53 +1072,10 @@ git-tree-sha1 = "62389eeff14780bfe55195b7204c0d8738436d64" uuid = "ae029012-a4dd-5104-9daa-d747884805df" version = "1.3.1" -[[deps.Rmath]] -deps = ["Random", "Rmath_jll"] -git-tree-sha1 = "4395a4cad612f95c1d08352f8c53811d6af3060b" -uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa" -version = "0.8.1" - -[[deps.Rmath_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl"] -git-tree-sha1 = "58cdd8fb2201a6267e1db87ff148dd6c1dbd8ad8" -uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f" -version = "0.5.1+0" - -[[deps.Rocket]] -deps = ["DataStructures", "Sockets", "Unrolled"] -git-tree-sha1 = "fe7373bf6b935c4431002fd91fa581d5eb835d09" -uuid = "df971d30-c9d6-4b37-b8ff-e965b2cb3a40" -version = "1.8.3" - -[[deps.RxInfer]] -deps = ["BayesBase", "DataStructures", "Dates", "Distributions", "DomainSets", "ExponentialFamily", "FastCholesky", "GraphPPL", "HTTP", "JSON", "LinearAlgebra", "Logging", "MacroTools", "Optim", "Preferences", "ProgressMeter", "Random", "ReactiveMP", "Reexport", "Rocket", "Static", "Statistics", "TupleTools", "UUIDs"] -git-tree-sha1 = "bcaf218d6b5329dc5e5be4299cbb5818c2c7bb19" -uuid = "86711068-29c9-4ff7-b620-ae75d7495b3d" -version = "4.6.2" - - [deps.RxInfer.extensions] - PrettyTablesExt = "PrettyTables" - ProjectionExt = "ExponentialFamilyProjection" - - [deps.RxInfer.weakdeps] - ExponentialFamilyProjection = "17f509fa-9a96-44ba-99b2-1c5f01f0931b" - PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" - [[deps.SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" version = "0.7.0" -[[deps.SciMLPublic]] -git-tree-sha1 = "ed647f161e8b3f2973f24979ec074e8d084f1bee" -uuid = "431bcebd-1456-4ced-9d72-93c2757fff0b" -version = "1.0.0" - -[[deps.ScopedValues]] -deps = ["HashArrayMappedTries", "Logging"] -git-tree-sha1 = "c3b2323466378a2ba15bea4b2f73b081e022f473" -uuid = "7e506255-f358-4e82-b7e4-beb19740aa63" -version = "1.5.0" - [[deps.Scratch]] deps = ["Dates"] git-tree-sha1 = "9b81b8393e50b7d4e6d0a9f14e192294d3b7c109" @@ -2611,17 +1086,6 @@ version = "1.3.0" uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" version = "1.11.0" -[[deps.Setfield]] -deps = ["ConstructionBase", "Future", "MacroTools", "StaticArraysCore"] -git-tree-sha1 = "c5391c6ace3bc430ca630251d02ea9687169ca68" -uuid = "efcf1570-3423-57d1-acb7-fd33fddbac46" -version = "1.1.2" - -[[deps.SharedArrays]] -deps = ["Distributed", "Mmap", "Random", "Serialization"] -uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383" -version = "1.11.0" - [[deps.Showoff]] deps = ["Dates", "Grisu"] git-tree-sha1 = "91eddf657aca81df9ae6ceb20b959ae5653ad1de" @@ -2633,12 +1097,6 @@ git-tree-sha1 = "f305871d2f381d21527c770d4788c06c097c9bc1" uuid = "777ac1f9-54b0-4bf8-805c-2214025038e7" version = "1.2.0" -[[deps.SimpleTraits]] -deps = ["InteractiveUtils", "MacroTools"] -git-tree-sha1 = "be8eeac05ec97d379347584fa9fe2f5f76795bcb" -uuid = "699a6c99-e7fa-54fc-8d76-47d257e15c1d" -version = "0.9.5" - [[deps.Sockets]] uuid = "6462fe0b-24de-5631-8697-dd941f90decc" version = "1.11.0" @@ -2654,49 +1112,12 @@ deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" version = "1.12.0" -[[deps.SpecialFunctions]] -deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] -git-tree-sha1 = "f2685b435df2613e25fc10ad8c26dddb8640f547" -uuid = "276daf66-3868-5448-9aa4-cd146d93841b" -version = "2.6.1" - - [deps.SpecialFunctions.extensions] - SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" - - [deps.SpecialFunctions.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - [[deps.StableRNGs]] deps = ["Random"] git-tree-sha1 = "95af145932c2ed859b63329952ce8d633719f091" uuid = "860ef19b-820b-49d6-a774-d7a799459cd3" version = "1.0.3" -[[deps.Static]] -deps = ["CommonWorldInvalidations", "IfElse", "PrecompileTools", "SciMLPublic"] -git-tree-sha1 = "49440414711eddc7227724ae6e570c7d5559a086" -uuid = "aedffcd0-7271-4cad-89d0-dc628f76c6d3" -version = "1.3.1" - -[[deps.StaticArrays]] -deps = ["LinearAlgebra", "PrecompileTools", "Random", "StaticArraysCore"] -git-tree-sha1 = "b8693004b385c842357406e3af647701fe783f98" -uuid = "90137ffa-7385-5640-81b9-e52037218182" -version = "1.9.15" - - [deps.StaticArrays.extensions] - StaticArraysChainRulesCoreExt = "ChainRulesCore" - StaticArraysStatisticsExt = "Statistics" - - [deps.StaticArrays.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" - -[[deps.StaticArraysCore]] -git-tree-sha1 = "192954ef1208c7019899fbf8049e717f92959682" -uuid = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" -version = "1.4.3" - [[deps.Statistics]] deps = ["LinearAlgebra"] git-tree-sha1 = "ae3bb1eb3bba077cd276bc5cfc337cc65c3075c0" @@ -2719,28 +1140,10 @@ git-tree-sha1 = "2c962245732371acd51700dbb268af311bddd719" uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" version = "0.34.6" -[[deps.StatsFuns]] -deps = ["HypergeometricFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] -git-tree-sha1 = "8e45cecc66f3b42633b8ce14d431e8e57a3e242e" -uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c" -version = "1.5.0" - - [deps.StatsFuns.extensions] - StatsFunsChainRulesCoreExt = "ChainRulesCore" - StatsFunsInverseFunctionsExt = "InverseFunctions" - - [deps.StatsFuns.weakdeps] - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" - [[deps.StyledStrings]] uuid = "f489334b-da3d-4c2e-b8f0-e476e12c162b" version = "1.11.0" -[[deps.SuiteSparse]] -deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] -uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" - [[deps.SuiteSparse_jll]] deps = ["Artifacts", "Libdl", "libblastrampoline_jll"] uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" @@ -2767,11 +1170,6 @@ deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" version = "1.11.0" -[[deps.TinyHugeNumbers]] -git-tree-sha1 = "83c6abf376718345a85c071b249ef6692a8936d4" -uuid = "783c9a47-75a3-44ac-a16b-f1ab7b3acf04" -version = "1.0.3" - [[deps.TranscodingStreams]] git-tree-sha1 = "0c45878dcfdcfa8480052b6ab162cdd138781742" uuid = "3bb67fe8-82b1-5028-8e26-92a6c54297fa" @@ -2782,29 +1180,6 @@ git-tree-sha1 = "372b90fe551c019541fafc6ff034199dc19c8436" uuid = "410a4b4d-49e4-4fbc-ab6d-cb71b17b3775" version = "0.1.12" -[[deps.Tullio]] -deps = ["DiffRules", "LinearAlgebra", "Requires"] -git-tree-sha1 = "972698b132b9df8791ae74aa547268e977b55f68" -uuid = "bc48ee85-29a4-5162-ae0b-a64e1601d4bc" -version = "0.3.8" - - [deps.Tullio.extensions] - TullioCUDAExt = "CUDA" - TullioChainRulesCoreExt = "ChainRulesCore" - TullioFillArraysExt = "FillArrays" - TullioTrackerExt = "Tracker" - - [deps.Tullio.weakdeps] - CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" - ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" - FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" - Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" - -[[deps.TupleTools]] -git-tree-sha1 = "41e43b9dc950775eac654b9f845c839cd2f1821e" -uuid = "9d95972d-f1c8-5527-a6e0-b4b365fa01f6" -version = "1.6.0" - [[deps.URIs]] git-tree-sha1 = "bef26fb046d031353ef97a82e3fdb6afe7f21b1a" uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4" @@ -2815,11 +1190,6 @@ deps = ["Random", "SHA"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" version = "1.11.0" -[[deps.UnPack]] -git-tree-sha1 = "387c1f73762231e86e0c9c5443ce3b4a0a9a0c2b" -uuid = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" -version = "1.0.2" - [[deps.Unicode]] uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" version = "1.11.0" @@ -2854,12 +1224,6 @@ git-tree-sha1 = "af305cc62419f9bd61b6644d19170a4d258c7967" uuid = "45397f5d-5981-4c77-b2b3-fc36d6e9b728" version = "1.7.0" -[[deps.Unrolled]] -deps = ["MacroTools"] -git-tree-sha1 = "6cc9d682755680e0f0be87c56392b7651efc2c7b" -uuid = "9602ed7d-8fef-5bc8-8597-8f21381861e8" -version = "0.1.5" - [[deps.Unzip]] git-tree-sha1 = "ca0969166a028236229f63514992fc073799bb78" uuid = "41fe7b60-77ed-43a1-b4f0-825fd5a5650d" @@ -3135,10 +1499,6 @@ version = "1.9.2+0" # ╔═╡ Cell order: # ╟─96547560-d294-11ef-0fa7-6b6489f7baba # ╟─af24aa27-b0a1-4c9b-aee0-0e5143d2f47e -# ╟─9654ea3e-d294-11ef-335c-657af1ceaf19 -# ╟─96552348-d294-11ef-16d8-b53563054687 -# ╟─9655959e-d294-11ef-0ca6-5f20aa579e91 -# ╟─05db1eab-5b63-4ab9-8b4a-ab2cde554295 # ╟─9655b2c2-d294-11ef-057f-9b3984064411 # ╟─403845d6-0229-4c93-bc23-3aecfd76e874 # ╟─9655d360-d294-11ef-0f06-ab58e2ad0e5f @@ -3148,10 +1508,6 @@ version = "1.9.2+0" # ╟─965606f2-d294-11ef-305b-870427879e50 # ╟─2d5cd42e-26ea-4ee4-979f-3a4e2c6271b6 # ╟─96561594-d294-11ef-1590-198382927808 -# ╟─9656cf72-d294-11ef-03aa-b715dd686c09 -# ╟─27b9e811-4542-4fe8-86a5-f8ba09074761 -# ╟─9656d850-d294-11ef-21a1-474b07ea7729 -# ╟─9658329c-d294-11ef-0d03-45e6872c4985 # ╟─f0181b53-a604-489f-a89e-db6fc58571dd # ╟─9862f675-87e6-4bba-a07d-9e51839819c7 # ╟─ea4a720f-a644-46a0-ad35-b215780e0928 @@ -3161,113 +1517,26 @@ version = "1.9.2+0" # ╟─9656fae2-d294-11ef-10d8-ff921d5956bd # ╟─b33b2aef-e672-490c-bdf4-a5f655fa4695 # ╟─f46473f8-e749-431d-8350-751115f0aaf0 -# ╟─cb3df230-6c7e-41b9-ba13-3c5f8a7fbb62 # ╟─96570d3e-d294-11ef-0178-c34dda717495 -# ╟─9657b088-d294-11ef-3017-e95c4c69b62b +# ╟─679d05b5-9fee-46e7-9627-c886c34a612d # ╟─0afe3cdc-15ed-4d9a-848a-d1977d051866 # ╟─96571c34-d294-11ef-11ef-29beeb1f96c2 -# ╟─253d4703-03d6-4961-8c3b-b70d2cbc0710 # ╟─a7b1f559-3c34-491e-83e7-ba95c8c22c80 -# ╟─70736e62-2b6c-4b3a-ab59-7e51522d620b # ╟─96575dd4-d294-11ef-31d6-b39b4c4bdea1 # ╟─5cc2016e-0383-448c-bd33-5b3a687b7436 # ╟─f65f5d0e-2583-4b88-b9f2-5fee15257c05 # ╟─7009cdc8-892c-499e-b932-b828fa300b6c -# ╟─91f81188-727c-4754-9a07-e754eef8bbe0 -# ╟─0633afea-5e92-4bad-8402-d159c534af81 -# ╟─f11564db-aafc-4df9-b494-4e5ced9bfcfe # ╟─9651f976-b834-4b81-8810-649f0290969d # ╟─96587a66-d294-11ef-2c7a-9fd7bea76582 # ╟─89e2757e-a09f-40c6-8dd7-9b4b4d232e17 # ╟─c4b5b124-e52a-41fc-b27e-a58181622e5c -# ╟─a1c957c1-69b7-4178-ab59-c0b2439bb01a -# ╟─9658c106-d294-11ef-01db-cfcff611ed81 -# ╟─96594d44-d294-11ef-22b8-95165fb08ce4 -# ╟─96597ce0-d294-11ef-3478-25c6bbef601e -# ╟─965998a8-d294-11ef-1d18-85876e3656c5 -# ╟─284a9dd5-1e26-4fd3-bb58-6e7ac0a0872f -# ╟─42dd67e6-eb0f-4368-9947-47de229f7be1 -# ╟─480165f9-33d9-4db1-bf05-8d99f0d9fb3e -# ╟─1c9c7994-672c-42a3-8ae7-8ce092ada9f0 -# ╟─7764541a-c11e-4e12-bbac-f8906cbc5dc6 -# ╟─8a2019af-9500-42c5-8408-ff93104a2d79 -# ╠═2f86a4cf-2075-45b5-bf2d-2d4d6888461a -# ╠═22f22f59-c320-4654-b472-b64cc3a001ff -# ╟─1a40ef8d-d677-4bf0-9186-18c5aa43a849 -# ╟─965a1df0-d294-11ef-323c-3da765f1104a -# ╠═fd338a30-9622-405a-96fa-caca6bd4ccfb -# ╠═1070063a-ef85-4527-ae82-1f01c1a506ff -# ╠═c03b1140-adce-467a-b953-50ad1bf3bc34 -# ╟─9431bc9a-bd83-4e4d-b64d-0571c1d01c87 -# ╠═83a70a4b-b114-4351-8fa2-dd565ebc9916 -# ╟─b3262127-69e0-4efb-875b-074d1d70437c -# ╟─fb61c774-34a3-493a-b149-c870993b6d46 -# ╟─5bcefd5f-4cd2-4cfe-8c1f-1129e5020d9a -# ╟─965a37e8-d294-11ef-340f-0930b229dd32 -# ╟─1832bffd-2729-4d3f-86f4-0e2d9ab26ba3 # ╟─4a10044c-e044-43e1-bd44-847f56019061 -# ╠═92f7bcfd-00a4-4cb7-a3eb-c1e101fdbcf6 -# ╟─965a6c20-d294-11ef-1c91-4bd237afbd20 +# ╟─1f16c8db-f2da-4a72-8a69-743ec857ac3b # ╟─fa5bdb1c-4412-48cc-950c-9ed92b4c9f76 # ╟─be670693-2036-46cb-8452-a2d0e1bf1172 -# ╟─25492eea-e649-43f9-b71f-ac6d1a80d0ee -# ╟─a5cd774f-57ad-4cb5-86c0-35987aa6e221 -# ╟─b6de3f00-d3b8-44d8-b72a-48cd5628b607 -# ╟─05375a01-4d1b-44cc-b1c4-a5eb4b6c5c5b -# ╟─206c34b3-1873-460b-911e-f2cd4f8886af -# ╟─45251c19-6eae-41e7-b0ed-8bd70a67d4e0 -# ╟─a6e155eb-7376-4e57-8e63-628934e14e78 -# ╟─9dc870d7-a5f3-447c-96ee-ad23199bc253 -# ╟─e8a35c28-6d6d-4066-8251-f091f28622a9 -# ╟─965a8a1a-d294-11ef-1d2f-65abf76665e8 -# ╟─965aa14c-d294-11ef-226f-65d587fefa64 -# ╟─56e8a1bd-ef80-4265-b926-e5e9e085b72f -# ╟─965ab77c-d294-11ef-2510-95b1a998589f -# ╟─965af708-d294-11ef-112c-f5470031dbbe -# ╟─965b11a4-d294-11ef-1d04-dbdf39ce91a3 -# ╟─965b25ac-d294-11ef-0b9a-9d5a50a76069 -# ╟─8dd6874c-12b7-47b6-b589-009849198024 -# ╟─bfbf3d09-23f5-4f54-96f6-bfe536cfc228 -# ╠═e7e4b6d0-bdf0-4a93-9a73-7971e6e33065 -# ╟─2f5415e5-70b1-47ea-9790-7ac953bca538 -# ╠═1b76ab6c-ffa2-40eb-a6c6-55d7097a5108 -# ╟─965b886e-d294-11ef-1b10-0319896874cf -# ╟─a3e11d46-5a22-4eb6-ba91-7258ba3c667e -# ╟─0efe10d8-1d0e-4a8f-8005-25ee261322b8 -# ╠═1be3121d-be18-46a1-9af9-f108a2257c22 -# ╟─e5658c95-6cd0-426f-b819-31f9f2c7eaf4 -# ╠═94ca674e-1a01-424c-8657-6510be7097c3 -# ╟─965c18f8-d294-11ef-2456-b945a46241f4 -# ╟─e0add49a-94ac-4247-8554-5a50d4abbebb -# ╟─d05277c1-fb9e-4b2b-bcbc-d8be5e63cab5 -# ╟─965c5f28-d294-11ef-324e-4df3e38b5045 -# ╠═d27f7af6-e094-44fa-8ba4-4ad2fa38f8bc -# ╟─90d62ba0-ca97-43f6-8f5a-0c1086a13f3d -# ╠═053e9dde-c088-4f15-9ca6-98b8185a8a11 -# ╠═07b09ac1-7fa7-4b62-b130-97315adb6fa7 -# ╟─defb2149-294b-47a8-99ed-1b3746b275f1 -# ╟─b3656d6c-4717-4fcd-90c6-ae4f4aa5e1be -# ╟─b15f28ce-c8c1-439b-aeca-74a58d2557e2 -# ╠═86e67c05-068d-4de4-80f3-1a20cc8a43ea -# ╠═fffa27d5-eb68-4dd3-9995-4a53fba6c1e4 -# ╟─578ec319-337d-4396-bb75-eaf99d95a38d # ╟─89da2fc0-a7c8-4a9d-82d9-622a311d010d # ╠═5a8dcadb-f0c2-4fb0-b8cd-db8cf49cc292 # ╠═965a08f4-d294-11ef-0604-1586ff37c0d4 -# ╠═981b08cc-7fb4-4880-8e8a-0b60a5dd72a2 -# ╠═2cb7d369-e7fd-4d66-8321-66a9197a26bd -# ╠═387f55ba-0fa4-4171-a405-3b1fb4e6b586 -# ╟─997235c1-08bd-4dbc-b1bc-cb10a3b83da4 -# ╠═96ef3cfb-ca18-46d6-bcac-0122c2c85fba -# ╠═aec4726a-954e-4e76-aae5-2dd6c979b12d -# ╠═79a0d02b-368f-4371-854c-cf2cea9328e5 # ╟─cb4427f9-0cb1-4393-b4de-14a4d64cc29c -# ╟─aca1f927-bc3b-48f6-af5c-12ee2ea4a49b -# ╠═99265e22-e8dc-40fe-989f-0d2a6c72faac -# ╟─f6fc4fad-70fb-432f-b77d-8e6ad42eef6c -# ╠═e20e9048-1271-41c7-97d3-635f320aa365 -# ╟─3a045b5c-9d87-46a6-a404-85c4bd77dd61 -# ╠═34ebbbe1-2a6b-422b-aeb1-cd2953acddca -# ╠═ba7a2dbd-f068-4249-bc29-77f2d0804676 # ╟─00000000-0000-0000-0000-000000000001 # ╟─00000000-0000-0000-0000-000000000002 diff --git a/mlss/Probability Theory Review.jl b/mlss/Probability Theory Review.jl index 1567806a..105ef3de 100644 --- a/mlss/Probability Theory Review.jl +++ b/mlss/Probability Theory Review.jl @@ -48,42 +48,6 @@ title("Probability Theory Review") # ╔═╡ bcb4be20-0439-4809-a166-8c50b6b9206b PlutoUI.TableOfContents() -# ╔═╡ 3e1803d0-d294-11ef-0304-df2b9b698cd1 -md""" -## Preliminaries - -##### Goal - -- Review of Probability Theory as a theory for rational/logical reasoning with uncertainties (i.e., a Bayesian interpretation) - -##### Materials - -- Mandatory - - - These lecture notes - -- Optional - - - Bishop pp. 12-24 - - - [3Blue1Brown, YouTube video on Bayes theorem (2019)](https://youtu.be/HZGCoVF3YvM?si=JaXdesPjU8B_BtrC) - - Nice animated tutorial on Bayes rule. - - - [Edwin Jaynes, Probability Theory–The Logic of Science (2003)](http://www.med.mcgill.ca/epidemiology/hanley/bios601/GaussianModel/JaynesProbabilityTheory.pdf). - - Brilliant book on the Bayesian view of probability theory. Just for fun, scan the annotated bibliography and references. - - - [Aubrey Clayton, Bernoulli's Fallacy–Statistical Illogic and the Crisis of Modern Science (2021)](https://aubreyclayton.com/bernoulli) - - A very readable account of the history of statistics and probability theory. Discusses why most popular statistics recipes are very poor scientific analysis tools. Use probability theory instead! - - - [Ariel Caticha, Entropic Inference and the Foundations of Physics (2012)](https://github.com/bmlip/course/blob/main/assets/files/Caticha-2012-Entropic-Inference-and-the-Foundations-of-Physics.pdf), pp.7-56 (ch.2: probability) - - Great introduction to probability theory, in particular w.r.t. its correct interpretation as a state-of-knowledge. - - Absolutely worth your time to read the whole chapter, even if you skip section 2.2.4 (pp.15-18) on Cox's proof. - - - [Joram Soch et al ., The Book of Statistical Proofs (2023 - )](https://statproofbook.github.io/) - - Online resource for proofs in probability theory and statistical inference. - -""" - # ╔═╡ 61713e1c-8e37-45d9-9f58-c3db69e15b66 challenge_statement("Disease Diagnosis",header_level=1) @@ -224,135 +188,6 @@ keyconcept(" ", """ ) -# ╔═╡ 3e18f18c-d294-11ef-33e4-b7f9495e0508 -md""" -## Why Probability Theory for Machine Learning? - -Machine learning concerns updating our beliefs about appropriate settings for model parameters from new information (namely a data set), and therefore PT provides the *optimal calculus for machine learning*. - -""" - -# ╔═╡ 3e1906ea-d294-11ef-236e-c966a9474170 -md""" -In general, nearly all interesting questions in machine learning (and information processing in general) can be stated in the following form (a conditional probability): - -```math -p(\texttt{whatever-we-want-to-know}\, | \,\texttt{whatever-we-do-know}) -``` - -where ``p(A|B)`` means the probability that ``A`` is true, given that ``B`` is true. - -""" - -# ╔═╡ 3e191b6c-d294-11ef-3174-d1b4b36e252b -md""" -##### Examples - - * Predictions - -```math -p(\,\texttt{future-observations}\,|\,\texttt{past-observations}\,) -``` - - * Classify a received data point ``x`` - -```math -p(\,x\texttt{-belongs-to-class-}k \,|\,x\,) -``` - - * Update a model based on a new observation - -```math -p(\,\texttt{model-parameters} \,|\,\texttt{new-observation},\,\texttt{past-observations}\,) -``` - -""" - -# ╔═╡ 3e192ef4-d294-11ef-1fc4-87175eeec5eb -md""" -## Frequentist vs. Bayesian Interpretation of Probabilities - -The interpretation of a probability as a **degree-of-belief** about the truth value of an event is also called the **Bayesian** interpretation. - -""" - -# ╔═╡ 3e19436c-d294-11ef-11c5-f9914f7a3a57 -md""" -In the **Bayesian** interpretation, the probability is associated with a **state-of-knowledge** (usually held by a person, but formally by a rational agent). - - * For instance, in a coin tossing experiment, ``p(\texttt{outcome} = \texttt{tail}) = 0.4`` should be interpreted as the belief that there is a 40% chance that ``\texttt{tail}`` comes up if the coin were tossed. - * Under the Bayesian interpretation, PT calculus (sum and product rules) **extends boolean logic to rational reasoning with uncertainty**. - -""" - -# ╔═╡ 4edf38ab-a940-4ab0-be22-fa95cf571146 -md""" -In the Bayesian interpretation, all probabilities are, in principle, conditional probabilities of the type ``p(A|I)``, since there is always some background knowledge. However, we often write ``p(A)`` rather than ``p(A|I)`` if the background knowledge ``I`` is assumed to be obviously present. E.g., we usually write ``p(A)`` rather than ``p(\,A\,|\,\text{the-sun-comes-up-tomorrow}\,)``. - -""" - -# ╔═╡ 3e194ef2-d294-11ef-3b38-1ddc3063ff35 -md""" -The Bayesian interpretation contrasts with the **frequentist** interpretation of a probability as the relative frequency that an event would occur under repeated execution of an experiment. - - * For instance, if the experiment is tossing a coin, then ``p(\texttt{outcome} = \texttt{tail}) = 0.4`` means that in the limit of a large number of coin tosses, 40% of outcomes turn up as ``\texttt{tail}``. - -""" - -# ╔═╡ 3e1964b4-d294-11ef-373d-712257fc130f -md""" -The Bayesian viewpoint is more generally applicable than the frequentist viewpoint, e.g., it is hard to apply the frequentist viewpoint to events like ``\texttt{"it will rain tomorrow"}``. - -""" - -# ╔═╡ 3e196d6a-d294-11ef-0795-41c045079251 -md""" -The Bayesian viewpoint is clearly favored in the machine learning community. (In this class, we also strongly favor the Bayesian interpretation). - -""" - -# ╔═╡ 3e198336-d294-11ef-26fd-03cd15876486 -md""" -Aubrey Clayton, in his wonderful book [Bernoulli's fallacy](https://aubreyclayton.com/bernoulli) (2021), writes about this issue: - -> “Compared with Bayesian methods, standard [frequentist] statistical techniques use only a small fraction of the available information about a research hypothesis (how well it predicts some observation), so naturally they will struggle when that limited information proves inadequate. Using standard statistical methods is like driving a car at night on a poorly lit highway: to keep from going in a ditch, we could build an elaborate system of bumpers and guardrails and equip the car with lane departure warnings and sophisticated navigation systems, and even then we could at best only drive to a few destinations. Or we could turn on the headlights.” - - -""" - -# ╔═╡ 3e198ba6-d294-11ef-3fe7-d70bf4833fa6 -md""" -In this class, we aim to turn on the headlights and illuminate the elegance and power of the Bayesian approach to information processing. - -""" - -# ╔═╡ 3e19e95a-d294-11ef-3da4-6d23922a5150 -md""" -## Variable Assignments as Propositions - - -""" - -# ╔═╡ 3e1a69f4-d294-11ef-103e-efc47025fb8f -md""" -If ``X`` is a variable, then an *assignment* ``X=x`` (where ``x`` is a value, e.g., ``X=5``) can be interpreted as an event. Hence, the expression ``p(X=5)`` should be interpreted as the *degree-of-belief of the event* that variable ``X`` takes on the value ``5``. - -""" - -# ╔═╡ 3e1a7c8e-d294-11ef-1f97-55e608d49141 -md""" -If ``X`` is a *discretely* valued variable, then ``p(X=x)`` is a probability *mass* function (PMF) with ``0\le p(X=x)\le 1`` and normalization ``\sum_x p(x) =1``. - -""" - -# ╔═╡ 3e1a8eca-d294-11ef-1ef0-c15b24d05990 -md""" -If ``X`` is *continuously* valued, then ``p(X=x)`` is a probability *density* function (PDF) with ``p(X=x)\ge 0`` and normalization ``\int_x p(x)\mathrm{d}x=1``. - - * Note that if ``X`` is continuously valued, then the value of ``p(x)`` is not necessarily ``\le 1``. E.g., a uniform distribution on the continuous domain ``[0,.5]`` has value ``p(x) = 2`` over its domain. - -""" - # ╔═╡ 3e1ab104-d294-11ef-1a98-412946949fba md""" # $(HTML("Probability Theory Calculus")) @@ -510,16 +345,6 @@ Note that marginalization can be understood as applying a "generalized" sum rule """ -# ╔═╡ 3e1bba8e-d294-11ef-1f61-295af16078ce -md""" -Of course, in the continuous domain, marginalization becomes - -```math -p(X)=\int_Y p(X,Y) \,\mathrm{d}Y -``` - -""" - # ╔═╡ 3e1bcb00-d294-11ef-2795-bd225bd00496 md""" ## $(HTML("Bayes Rule")) @@ -562,223 +387,66 @@ Bayes rule tells us how to update our knowledge about model parameters when faci # ╔═╡ 16c2eb59-16b8-4347-9aab-6e4b99016c79 keyconcept("", md"Bayes rule is the fundamental rule for learning from data!") -# ╔═╡ 3e1bffec-d294-11ef-2a49-9ff0f6331add +# ╔═╡ 3e1d33c8-d294-11ef-0a08-bdc419949925 md""" -## Bayes Rule Nomenclature - -Some nomenclature associated with Bayes rule: - -```math -\underbrace{p(\theta | D)}_{\text{posterior}} = \frac{\overbrace{p(D|\theta)}^{\text{likelihood}} \times \overbrace{p(\theta)}^{\text{prior}}}{\underbrace{p(D)}_{\text{evidence}}} -``` - -""" +## Probabilistic Inference -# ╔═╡ 3e1c0e80-d294-11ef-0d19-375e01988f16 -md""" -Note that the evidence (a.k.a. *marginal likelihood* ) can be computed from the numerator through marginalization since +**Probabilistic inference** refers to computing ```math - p(D) = \int p(D,\theta) \,\mathrm{d}\theta = \int p(D|\theta)\,p(\theta) \,\mathrm{d}\theta +p(\,\text{whatever-we-want-to-know}\, | \,\text{whatever-we-already-know}\,) ``` -""" - -# ╔═╡ 3e1c1e3e-d294-11ef-0955-bdf9d0ba3c53 -md""" -Hence, having access to likelihood and prior is in principle sufficient to compute both the evidence and the posterior. To emphasize that point, Bayes rule is sometimes written as a transformation: +For example: ```math - \underbrace{\underbrace{p(\theta|D)}_{\text{posterior}}\cdot \underbrace{p(D)}_{\text{evidence}}}_{\text{this is what we want to compute}} = \underbrace{\underbrace{p(D|\theta)}_{\text{likelihood}}\cdot \underbrace{p(\theta)}_{\text{prior}}}_{\text{this is available}} +\begin{align*} + p(\,\text{Mr.S.-killed-Mrs.S.} \;&|\; \text{he-has-her-blood-on-his-shirt}\,) \\ + p(\,\text{transmitted-codeword} \;&|\;\text{received-codeword}\,) + \end{align*} ``` -""" +This can be accomplished by repeatedly applying the sum and product rules. -# ╔═╡ 3e1c4224-d294-11ef-2707-49470aaae6eb -md""" -For a given data set ``D``, the posterior probabilities of the parameters scale relatively against each other as +In particular, consider a joint distribution ``p(X,Y,Z)``. Assume we are interested in ``p(X|Z)``: ```math -p(\theta|D) \propto p(D|\theta) p(\theta) +\begin{align*} +p(X|Z) \stackrel{p}{=} \frac{p(X,Z)}{p(Z)} \stackrel{s}{=} \frac{\sum_Y p(X,Y,Z)}{\sum_{X,Y} p(X,Y,Z)} \,, +\end{align*} ``` -Hence, all that we can learn from the observed data is contained in the likelihood function ``p(D|\theta)``. This is called the **likelihood principle**. +where the ``s`` and ``p`` above the equality sign indicate whether the sum or product rule was used. +In the rest of this course, we'll encounter many lengthy probabilistic derivations. For each manipulation, you should be able to associate an 's' (for sum rule), a 'p' (for product or Bayes rule) or an 'm' (for a simplifying model assumption like conditional independency) above any equality sign. """ -# ╔═╡ 3e1c51e2-d294-11ef-2c6d-d32a98308c6f -md""" -## The Likelihood Function vs the Sampling Distribution - -Consider a distribution ``p(D|\theta)``, where ``D`` relates to variables that are observed (i.e., a "data set") and ``\theta`` are model parameters. - -""" +# ╔═╡ 3e1b05ee-d294-11ef-33de-efed64d01c0d +keyconcept( + "", + md""" + All valid probabilistic relations can be derived from just two fundamental principles: the **sum rule** and the **product rule**. These two rules form the foundation of probability theory, from which more complex constructs such as conditional probabilities, Bayes’ theorem, and marginalization naturally follow. + + """ +) -# ╔═╡ 3e1c60ba-d294-11ef-3a01-cf9e97512857 +# ╔═╡ b176ceae-884e-4460-9f66-020c1ac447f1 md""" -In general, ``p(D|\theta)`` is just a function of the two variables ``D`` and ``\theta``. We distinguish two interpretations of this function, depending on which variable is observed (or given by other means). - +# Examples """ -# ╔═╡ 3e1c70be-d294-11ef-14ed-0d46515541c5 -md""" -The **sampling distribution** (a.k.a. the **data-generating** distribution) +# ╔═╡ ab223dea-8ba8-4d30-94f4-72c8e070aadf +θ_bond = @bind θ Scrubbable(0.0:0.02:1; format=".2f"); -```math -p(D|\theta=\theta_0) -``` +# ╔═╡ 922770f4-ddc8-4089-b378-f14088276b43 +exercise_statement("Which color does the ball have?"; prefix="Inference ") -(which is a function of ``D`` only) describes a probability distribution for data ``D``, assuming that it is generated by the given model with parameters fixed at ``\theta = \theta_0``. +# ╔═╡ 3e1de32c-d294-11ef-1f63-f190c8361404 +md""" -""" +##### Problem -# ╔═╡ 3e1c806a-d294-11ef-1fad-17e5625279f7 -md""" -In a machine learning context, often the data is observed, and ``\theta`` is the free variable. In that case, for given observations ``D=D_0``, the **likelihood function** (which is a function only of the model parameters ``\theta``) is defined as - -```math -L(\theta) \triangleq p(D=D_0|\theta) -``` - -""" - -# ╔═╡ 3e1c9184-d294-11ef-3e35-5393d97fbc44 -md""" -Note that ``L(\theta)`` is not a probability distribution for ``\theta`` since in general ``\sum_\theta L(\theta) \neq 1``. - -""" - -# ╔═╡ 3e1d33c8-d294-11ef-0a08-bdc419949925 -md""" -## Probabilistic Inference - -**Probabilistic inference** refers to computing - -```math -p(\,\text{whatever-we-want-to-know}\, | \,\text{whatever-we-already-know}\,) -``` - -For example: - -```math -\begin{align*} - p(\,\text{Mr.S.-killed-Mrs.S.} \;&|\; \text{he-has-her-blood-on-his-shirt}\,) \\ - p(\,\text{transmitted-codeword} \;&|\;\text{received-codeword}\,) - \end{align*} -``` - -This can be accomplished by repeatedly applying the sum and product rules. - -In particular, consider a joint distribution ``p(X,Y,Z)``. Assume we are interested in ``p(X|Z)``: - -```math -\begin{align*} -p(X|Z) \stackrel{p}{=} \frac{p(X,Z)}{p(Z)} \stackrel{s}{=} \frac{\sum_Y p(X,Y,Z)}{\sum_{X,Y} p(X,Y,Z)} \,, -\end{align*} -``` - -where the ``s`` and ``p`` above the equality sign indicate whether the sum or product rule was used. - -In the rest of this course, we'll encounter many lengthy probabilistic derivations. For each manipulation, you should be able to associate an 's' (for sum rule), a 'p' (for product or Bayes rule) or an 'm' (for a simplifying model assumption like conditional independency) above any equality sign. -""" - -# ╔═╡ 3e1b05ee-d294-11ef-33de-efed64d01c0d -keyconcept( - "", - md""" - All valid probabilistic relations can be derived from just two fundamental principles: the **sum rule** and the **product rule**. These two rules form the foundation of probability theory, from which more complex constructs such as conditional probabilities, Bayes’ theorem, and marginalization naturally follow. - - """ -) - -# ╔═╡ b176ceae-884e-4460-9f66-020c1ac447f1 -md""" -# Examples -""" - -# ╔═╡ e3157dc0-5a64-4479-a37a-40fe25cccc07 -code_example("Sampling Distribution and Likelihood Function for the Coin Toss") - -# ╔═╡ 7d493e09-f7cc-4e13-a506-b792edcbf390 -md""" - -and the likelihood function - -```math -L(\theta) \triangleq p(y=1|\theta) = \theta \,. -``` - -""" - -# ╔═╡ 3e1d20e0-d294-11ef-2044-e1fe6590a600 -md""" -!!! note - The (discrete) sampling distribution is a valid probability distribution. - - However, the likelihood function ``L(\theta)`` clearly isn't, since ``\int_0^1 L(\theta) \mathrm{d}\theta = 0.5 \neq 1``. -""" - -# ╔═╡ fc733d61-fd0f-4a13-9afc-4505ac0253df -f(y,θ) = θ.^y .* (1 .- θ).^(1 .- y) # p(y|θ) - -# ╔═╡ ab223dea-8ba8-4d30-94f4-72c8e070aadf -θ_bond = @bind θ Scrubbable(0.0:0.02:1; format=".2f"); - -# ╔═╡ d93f73d4-2783-4777-b0ce-cdc0444cb300 -md""" - -Consider the following simple model for the outcome ``y \in \{0,1\}`` (tail = ``0``, head = ``1``) of a biased coin toss with a real parameter $θ_bond ``= \theta \in [0,1]``: - -```math -\begin{align*} -p(y|\theta) = \theta^y (1-\theta)^{1-y}\\ -\end{align*} -``` - -Next, we use Julia to plot both the sampling distribution - -""" - -# ╔═╡ 8a7dd8b7-5faf-4091-8451-9769f842accb -let - p1 = plot( - [0,1], f([0,1], θ); - line=:stem, - marker=:circle, - xrange=(-0.5, 1.5), yrange=(0,1), - title="Sampling Distribution", - xlabel="y", ylabel=L"p(y|θ=%$θ)", label="" - ) - - _θ = 0:0.01:1 - y=1 - p2 = plot( - _θ, f(y, _θ); - ylabel=L"p(y=%$y | θ)", xlabel=L"θ", - title="Likelihood Function", label="" - ) - scatter!(p2, - [θ], [f(y, θ)]; - label=nothing, - ) - - plot(p1, p2) -end - -# ╔═╡ b7445b9b-7fbb-4560-b947-a23af0fcf101 -md""" -Click and drag this number to change ``\theta``: $θ_bond. -""" - -# ╔═╡ 922770f4-ddc8-4089-b378-f14088276b43 -exercise_statement("Which color does the ball have?"; prefix="Inference ") - -# ╔═╡ 3e1de32c-d294-11ef-1f63-f190c8361404 -md""" - -##### Problem - -- A bag contains one ball, known to be either white or black. A white ball is put in , and the bag is shaken. Next, a ball is drawn out, which proves to be white. If we now take out another ball, what is the probability it will be white? +- A bag contains one ball, known to be either white or black. A white ball is put in , and the bag is shaken. Next, a ball is drawn out, which proves to be white. If we now take out another ball, what is the probability it will be white? """ @@ -864,6 +532,477 @@ Show solution? $(@bind(show_disease_diagnosis_solution, CheckBox(default=false)) """ +# ╔═╡ 9b92fc89-2036-4525-979b-d296ab29329c +md""" +# Summary +""" + +# ╔═╡ 5da42e9a-4318-48ea-9f43-4c1e1c97bceb +keyconceptsummary() + +# ╔═╡ 03692f4d-0daf-4dfc-a7ff-6b954326e4d0 +exercises() + + +# ╔═╡ 3a1d380e-df80-4727-9772-f199214cf05d +md""" +##### The Sum Rule (**) + +Derive the general sum rule, +```math +p(A + B) = p(A) + p(B) - p(A,B) +``` +from the elementary sum rule ``p(A) + p(\bar A) = 1`` and the sum and product rules. + +""" + +# ╔═╡ 99d9099f-4908-4bb3-8d59-da9cb69af04c +hint( + md""" + Here, you may make use of the (Boolean logic) fact that ``A + B = \overline {\bar A \bar B }``. + """ +) + +# ╔═╡ 3b1b0869-b815-4697-9dba-3c4b4cb5ac47 +hide_solution( +md""" +```math +\begin{align} +p\left( A + B \right) &\underset{\mathrm{bool}}{=} p\left( \overline {\bar A \bar B } \right) \\ + &\underset{\mathrm{sum}}{=} 1 - p\left( \bar{A} \bar{B} \right) \\ + &\underset{\mathrm{prod}}{=} 1 - p\left( \bar{A} |\bar{B} \right) p\left(\bar{B} \right) \\ + &\underset{\mathrm{sum}}{=} 1 - \left( 1 - p\left(A|\bar B \right) \right) \left( 1 - p\left( B \right) \right) \\ + &= p(B) + \left( {1 - p\left( B \right)} \right)p\left( {A|\bar B } \right) \\ + &\underset{\mathrm{prod}}{=} p(B) + \left( 1 - p\left( B \right) \right) p\left( \bar{B} |A \right) \frac{ p\left( A \right) }{ p\left(\bar{B}\right)} \\ + &\underset{\mathrm{sum}}{=} p(B) + p\left(\bar{B} |A \right) p\left( A \right) \\ + &\underset{\mathrm{sum}}{=} p(B) + \left( 1 - p\left( {B|A} \right) \right) p\left( A \right) \\ + &\underset{\mathrm{sum}}{=} p\left( A \right) + p(B) - p\left( A,B \right) +\end{align} +``` +Note that, aside from the first boolean rewrite, everything follows straight application of sum and product rules. + + +""") + +# ╔═╡ 5f377237-d9a5-4778-aa4d-1c6ce109b705 +md""" +##### Apples and Oranges + +Box 1 contains 8 apples and 4 oranges. Box 2 contains 10 apples and 2 oranges. Boxes are chosen with equal probability. You pick a box and then select a fruit from that box. +- (a) (*) What is the probability of choosing an apple? +- (b) (**) If an apple is chosen, what is the probability that it came from box 1? +""" + +# ╔═╡ 5613e9b7-ff0d-435a-9de6-aaf293ebf592 +hide_solution( +md""" +The following probabilities are given in the problem statement, +```math +\begin{align} +p(b_1) &= p(b_2) = 1/2 \\ +p(a|b_1) &= 8/12, \quad p(a|b_2) = 10/12 \\ +p(o|b_1) &= 4/12, \quad p(o|b_2) = 2/12 +\end{align} +``` +(a) +```math +p(a) = \sum_i p(a,b_i) = \sum_i p(a|b_i)p(b_i)=\frac{8}{12}\cdot\frac{1}{2} + \frac{10}{12}\cdot\frac{1}{2} = \frac{3}{4} +``` +(b) +```math +p(b_1|a) = \frac{p(a,b_1)}{p(a)} = \frac{p(a|b_1)p(b_1)}{p(a)} = \frac{\frac{8}{12}\cdot\frac{1}{2}}{\frac{3}{4}} = \frac{4}{9} +``` +""" +) + +# ╔═╡ fc3151f9-e143-4e31-b7b7-3f25b4fe9dab +md""" +##### What is a Random Signal? (*) +Is a speech signal a "probabilistic" (random) or a deterministic signal? +""" + +# ╔═╡ 66ebe33c-8360-4938-9b51-625e5bed176c +hide_solution( +md""" +That depends. The term “probabilistic” refers to a state-of-knowledge (or beliefs) about something—in this case, about the values of a speech signal. The key point is that the signal itself is neither inherently probabilistic nor deterministic; these labels describe our knowledge about it. + +If you had a perfect microphone and recorded the speech signal flawlessly at its source, you would know all its values exactly—no uncertainty—so you could call it deterministic. + +However, before making the recording, how would you represent your knowledge about the signal values you are going to measure? You face uncertainty, so the appropriate description is a probability distribution over all possible signal values. +""") + +# ╔═╡ 5b681e41-ad14-4c58-8ea0-4b6d85885c51 +md""" +##### Who Speaks the Truth? (***) +The inhabitants of an island tell the truth one-third of the time. They lie with probability ``2/3``. On an occasion, after one of them made a statement, you ask another person "was that statement true?" and he says "yes". What is the probability that the statement was indeed true? + +""" + +# ╔═╡ 91dd40f0-c373-48b3-b83b-6e8df2c43e5a +hide_solution( +md""" +We use variables ``S_1 \in \{\text{t},\text{f}\}`` and ``S_2 \in \{\text{y},\text{n}\}`` for statements 1 and 2 and shorthand "y", "n", "t" and "f" for "yes", "no", "true" and "false", respectively. The problem statement provides us with the following probabilities, +```math +\begin{align} +p(S_1=\text{t}) &= 1/3 \\ +p(S_1=\text{f}) &= 1 - p(S_1=\text{t}) = 2/3\\ +p(S_2=\text{y} | S_1=\text{t}) &= 1/3 \\ +p(S_2=\text{y} | S_1=\text{f}) &= 2/3 +\end{align} +``` +We are asked to compute ``p(S_1=\text{t} | S_2=\text{y})``. Use Bayes rule, +```math +\begin{align} +p(S_1=\text{t} | S_2=\text{y}) &= \frac{p(S_1=\text{t},S_2=\text{y})}{p(S_2=\text{y})} \\ +&= \frac{\overbrace{p(S_2=\text{y}|S_1=\text{t})p(S_1=\text{t})}^{\text{both speak the truth}}}{\underbrace{p(S_2=\text{y}|S_1=\text{t})p(S_1=\text{t})}_{\text{both speak the truth}}+\underbrace{p(S_2=\text{y}|S_1=\text{f})p(S_1=\text{f})}_{\text{both lie}}}\\ +&= \frac{\frac{1}{3}\cdot\frac{1}{3}}{\frac{1}{3}\cdot\frac{1}{3}+\frac{2}{3}\cdot\frac{2}{3}} = \frac{1}{5} +\end{align} +``` +""") + +# ╔═╡ a8d4a517-84a7-426e-a49e-482c5fd047ae +md""" +##### The Likelihood Function is a Function of What? (*) + +When considering the distribution ``p(D|\theta)``, is it more correct to speak about the likelihood of the model parameters ``\theta`` than about the likelihood of the observed data set ``D``. And why? + +""" + +# ╔═╡ d3b003c6-70ca-419f-a343-e35b266323f3 +hide_solution( +md""" +Yes, it’s more correct to speak about the likelihood of the model parameters, not of the observed data set. Once ``D`` has been observed, it is no longer a random variable; it’s just a fixed outcome. What varies is ``\theta``, so ``L(\theta) = p(D|\theta)`` is a function of the parameters, not of the data. + +Saying “likelihood of the data” is misleading because it confuses likelihood with the sampling distribution ``p(D|\theta)`` seen as a function of ``D`` (where ``\theta`` is fixed). The latter is a probability distribution over possible data sets before observing them. +""") + +# ╔═╡ be66b697-f920-4361-9ff2-b12cc50ae8c9 +md""" +# Optional Slides +""" + +# ╔═╡ 3e1803d0-d294-11ef-0304-df2b9b698cd1 +md""" +## Preliminaries + +##### Goal + +- Review of Probability Theory as a theory for rational/logical reasoning with uncertainties (i.e., a Bayesian interpretation) + +##### Materials + +- Mandatory + + - These lecture notes + +- Optional + + - Bishop pp. 12-24 + + - [3Blue1Brown, YouTube video on Bayes theorem (2019)](https://youtu.be/HZGCoVF3YvM?si=JaXdesPjU8B_BtrC) + - Nice animated tutorial on Bayes rule. + + - [Edwin Jaynes, Probability Theory–The Logic of Science (2003)](http://www.med.mcgill.ca/epidemiology/hanley/bios601/GaussianModel/JaynesProbabilityTheory.pdf). + - Brilliant book on the Bayesian view of probability theory. Just for fun, scan the annotated bibliography and references. + + - [Aubrey Clayton, Bernoulli's Fallacy–Statistical Illogic and the Crisis of Modern Science (2021)](https://aubreyclayton.com/bernoulli) + - A very readable account of the history of statistics and probability theory. Discusses why most popular statistics recipes are very poor scientific analysis tools. Use probability theory instead! + + - [Ariel Caticha, Entropic Inference and the Foundations of Physics (2012)](https://github.com/bmlip/course/blob/main/assets/files/Caticha-2012-Entropic-Inference-and-the-Foundations-of-Physics.pdf), pp.7-56 (ch.2: probability) + - Great introduction to probability theory, in particular w.r.t. its correct interpretation as a state-of-knowledge. + - Absolutely worth your time to read the whole chapter, even if you skip section 2.2.4 (pp.15-18) on Cox's proof. + + - [Joram Soch et al ., The Book of Statistical Proofs (2023 - )](https://statproofbook.github.io/) + - Online resource for proofs in probability theory and statistical inference. + +""" + +# ╔═╡ 3e18f18c-d294-11ef-33e4-b7f9495e0508 +md""" +## Why Probability Theory for Machine Learning? + +Machine learning concerns updating our beliefs about appropriate settings for model parameters from new information (namely a data set), and therefore PT provides the *optimal calculus for machine learning*. + +""" + +# ╔═╡ 3e1906ea-d294-11ef-236e-c966a9474170 +md""" +In general, nearly all interesting questions in machine learning (and information processing in general) can be stated in the following form (a conditional probability): + +```math +p(\texttt{whatever-we-want-to-know}\, | \,\texttt{whatever-we-do-know}) +``` + +where ``p(A|B)`` means the probability that ``A`` is true, given that ``B`` is true. + +""" + +# ╔═╡ 3e191b6c-d294-11ef-3174-d1b4b36e252b +md""" +##### Examples + + * Predictions + +```math +p(\,\texttt{future-observations}\,|\,\texttt{past-observations}\,) +``` + + * Classify a received data point ``x`` + +```math +p(\,x\texttt{-belongs-to-class-}k \,|\,x\,) +``` + + * Update a model based on a new observation + +```math +p(\,\texttt{model-parameters} \,|\,\texttt{new-observation},\,\texttt{past-observations}\,) +``` + +""" + +# ╔═╡ 3e192ef4-d294-11ef-1fc4-87175eeec5eb +md""" +## Frequentist vs. Bayesian Interpretation of Probabilities + +The interpretation of a probability as a **degree-of-belief** about the truth value of an event is also called the **Bayesian** interpretation. + +""" + +# ╔═╡ 3e19436c-d294-11ef-11c5-f9914f7a3a57 +md""" +In the **Bayesian** interpretation, the probability is associated with a **state-of-knowledge** (usually held by a person, but formally by a rational agent). + + * For instance, in a coin tossing experiment, ``p(\texttt{outcome} = \texttt{tail}) = 0.4`` should be interpreted as the belief that there is a 40% chance that ``\texttt{tail}`` comes up if the coin were tossed. + * Under the Bayesian interpretation, PT calculus (sum and product rules) **extends boolean logic to rational reasoning with uncertainty**. + +""" + +# ╔═╡ 4edf38ab-a940-4ab0-be22-fa95cf571146 +md""" +In the Bayesian interpretation, all probabilities are, in principle, conditional probabilities of the type ``p(A|I)``, since there is always some background knowledge. However, we often write ``p(A)`` rather than ``p(A|I)`` if the background knowledge ``I`` is assumed to be obviously present. E.g., we usually write ``p(A)`` rather than ``p(\,A\,|\,\text{the-sun-comes-up-tomorrow}\,)``. + +""" + +# ╔═╡ 3e194ef2-d294-11ef-3b38-1ddc3063ff35 +md""" +The Bayesian interpretation contrasts with the **frequentist** interpretation of a probability as the relative frequency that an event would occur under repeated execution of an experiment. + + * For instance, if the experiment is tossing a coin, then ``p(\texttt{outcome} = \texttt{tail}) = 0.4`` means that in the limit of a large number of coin tosses, 40% of outcomes turn up as ``\texttt{tail}``. + +""" + +# ╔═╡ 3e1964b4-d294-11ef-373d-712257fc130f +md""" +The Bayesian viewpoint is more generally applicable than the frequentist viewpoint, e.g., it is hard to apply the frequentist viewpoint to events like ``\texttt{"it will rain tomorrow"}``. + +""" + +# ╔═╡ 3e196d6a-d294-11ef-0795-41c045079251 +md""" +The Bayesian viewpoint is clearly favored in the machine learning community. (In this class, we also strongly favor the Bayesian interpretation). + +""" + +# ╔═╡ 3e198336-d294-11ef-26fd-03cd15876486 +md""" +Aubrey Clayton, in his wonderful book [Bernoulli's fallacy](https://aubreyclayton.com/bernoulli) (2021), writes about this issue: + +> “Compared with Bayesian methods, standard [frequentist] statistical techniques use only a small fraction of the available information about a research hypothesis (how well it predicts some observation), so naturally they will struggle when that limited information proves inadequate. Using standard statistical methods is like driving a car at night on a poorly lit highway: to keep from going in a ditch, we could build an elaborate system of bumpers and guardrails and equip the car with lane departure warnings and sophisticated navigation systems, and even then we could at best only drive to a few destinations. Or we could turn on the headlights.” + + +""" + +# ╔═╡ 3e198ba6-d294-11ef-3fe7-d70bf4833fa6 +md""" +In this class, we aim to turn on the headlights and illuminate the elegance and power of the Bayesian approach to information processing. + +""" + +# ╔═╡ 3e19e95a-d294-11ef-3da4-6d23922a5150 +md""" +## Variable Assignments as Propositions + + +""" + +# ╔═╡ 3e1a69f4-d294-11ef-103e-efc47025fb8f +md""" +If ``X`` is a variable, then an *assignment* ``X=x`` (where ``x`` is a value, e.g., ``X=5``) can be interpreted as an event. Hence, the expression ``p(X=5)`` should be interpreted as the *degree-of-belief of the event* that variable ``X`` takes on the value ``5``. + +""" + +# ╔═╡ 3e1a7c8e-d294-11ef-1f97-55e608d49141 +md""" +If ``X`` is a *discretely* valued variable, then ``p(X=x)`` is a probability *mass* function (PMF) with ``0\le p(X=x)\le 1`` and normalization ``\sum_x p(x) =1``. + +""" + +# ╔═╡ 3e1a8eca-d294-11ef-1ef0-c15b24d05990 +md""" +If ``X`` is *continuously* valued, then ``p(X=x)`` is a probability *density* function (PDF) with ``p(X=x)\ge 0`` and normalization ``\int_x p(x)\mathrm{d}x=1``. + + * Note that if ``X`` is continuously valued, then the value of ``p(x)`` is not necessarily ``\le 1``. E.g., a uniform distribution on the continuous domain ``[0,.5]`` has value ``p(x) = 2`` over its domain. + +""" + +# ╔═╡ 3e1bffec-d294-11ef-2a49-9ff0f6331add +md""" +## Bayes Rule Nomenclature + +Some nomenclature associated with Bayes rule: + +```math +\underbrace{p(\theta | D)}_{\text{posterior}} = \frac{\overbrace{p(D|\theta)}^{\text{likelihood}} \times \overbrace{p(\theta)}^{\text{prior}}}{\underbrace{p(D)}_{\text{evidence}}} +``` + +""" + +# ╔═╡ 3e1c0e80-d294-11ef-0d19-375e01988f16 +md""" +Note that the evidence (a.k.a. *marginal likelihood* ) can be computed from the numerator through marginalization since + +```math + p(D) = \int p(D,\theta) \,\mathrm{d}\theta = \int p(D|\theta)\,p(\theta) \,\mathrm{d}\theta +``` + +""" + +# ╔═╡ 3e1c1e3e-d294-11ef-0955-bdf9d0ba3c53 +md""" +Hence, having access to likelihood and prior is in principle sufficient to compute both the evidence and the posterior. To emphasize that point, Bayes rule is sometimes written as a transformation: + +```math + \underbrace{\underbrace{p(\theta|D)}_{\text{posterior}}\cdot \underbrace{p(D)}_{\text{evidence}}}_{\text{this is what we want to compute}} = \underbrace{\underbrace{p(D|\theta)}_{\text{likelihood}}\cdot \underbrace{p(\theta)}_{\text{prior}}}_{\text{this is available}} +``` + +""" + +# ╔═╡ 3e1c4224-d294-11ef-2707-49470aaae6eb +md""" +For a given data set ``D``, the posterior probabilities of the parameters scale relatively against each other as + +```math +p(\theta|D) \propto p(D|\theta) p(\theta) +``` + +Hence, all that we can learn from the observed data is contained in the likelihood function ``p(D|\theta)``. This is called the **likelihood principle**. + +""" + +# ╔═╡ 3e1c51e2-d294-11ef-2c6d-d32a98308c6f +md""" +## The Likelihood Function vs the Sampling Distribution + +Consider a distribution ``p(D|\theta)``, where ``D`` relates to variables that are observed (i.e., a "data set") and ``\theta`` are model parameters. + +""" + +# ╔═╡ 3e1c60ba-d294-11ef-3a01-cf9e97512857 +md""" +In general, ``p(D|\theta)`` is just a function of the two variables ``D`` and ``\theta``. We distinguish two interpretations of this function, depending on which variable is observed (or given by other means). + +""" + +# ╔═╡ 3e1c70be-d294-11ef-14ed-0d46515541c5 +md""" +The **sampling distribution** (a.k.a. the **data-generating** distribution) + +```math +p(D|\theta=\theta_0) +``` + +(which is a function of ``D`` only) describes a probability distribution for data ``D``, assuming that it is generated by the given model with parameters fixed at ``\theta = \theta_0``. + +""" + +# ╔═╡ 3e1c806a-d294-11ef-1fad-17e5625279f7 +md""" +In a machine learning context, often the data is observed, and ``\theta`` is the free variable. In that case, for given observations ``D=D_0``, the **likelihood function** (which is a function only of the model parameters ``\theta``) is defined as + +```math +L(\theta) \triangleq p(D=D_0|\theta) +``` + +""" + +# ╔═╡ 3e1c9184-d294-11ef-3e35-5393d97fbc44 +md""" +Note that ``L(\theta)`` is not a probability distribution for ``\theta`` since in general ``\sum_\theta L(\theta) \neq 1``. + +""" + +# ╔═╡ e3157dc0-5a64-4479-a37a-40fe25cccc07 +code_example("Sampling Distribution and Likelihood Function for the Coin Toss") + +# ╔═╡ d93f73d4-2783-4777-b0ce-cdc0444cb300 +md""" + +Consider the following simple model for the outcome ``y \in \{0,1\}`` (tail = ``0``, head = ``1``) of a biased coin toss with a real parameter $θ_bond ``= \theta \in [0,1]``: + +```math +\begin{align*} +p(y|\theta) = \theta^y (1-\theta)^{1-y}\\ +\end{align*} +``` + +Next, we use Julia to plot both the sampling distribution + +""" + +# ╔═╡ 7d493e09-f7cc-4e13-a506-b792edcbf390 +md""" + +and the likelihood function + +```math +L(\theta) \triangleq p(y=1|\theta) = \theta \,. +``` + +""" + +# ╔═╡ b7445b9b-7fbb-4560-b947-a23af0fcf101 +md""" +Click and drag this number to change ``\theta``: $θ_bond. +""" + +# ╔═╡ 3e1d20e0-d294-11ef-2044-e1fe6590a600 +md""" +!!! note + The (discrete) sampling distribution is a valid probability distribution. + + However, the likelihood function ``L(\theta)`` clearly isn't, since ``\int_0^1 L(\theta) \mathrm{d}\theta = 0.5 \neq 1``. +""" + +# ╔═╡ fc733d61-fd0f-4a13-9afc-4505ac0253df +f(y,θ) = θ.^y .* (1 .- θ).^(1 .- y) # p(y|θ) + +# ╔═╡ 8a7dd8b7-5faf-4091-8451-9769f842accb +let + p1 = plot( + [0,1], f([0,1], θ); + line=:stem, + marker=:circle, + xrange=(-0.5, 1.5), yrange=(0,1), + title="Sampling Distribution", + xlabel="y", ylabel=L"p(y|θ=%$θ)", label="" + ) + + _θ = 0:0.01:1 + y=1 + p2 = plot( + _θ, f(y, _θ); + ylabel=L"p(y=%$y | θ)", xlabel=L"θ", + title="Likelihood Function", label="" + ) + scatter!(p2, + [θ], [f(y, θ)]; + label=nothing, + ) + + plot(p1, p2) +end + # ╔═╡ ef264651-854e-4374-8ea8-5476c85150c4 md"# Moments and Transformations" @@ -1244,150 +1383,6 @@ p_y(y) &= p_x(g(y)) \cdot g^\prime(y) \\ In the statistics literature, ``y = \frac{x-\mu}{\sigma}`` is called the **standardized** variable since it transforms a general normal variable into a standard normal one.) """) -# ╔═╡ 9b92fc89-2036-4525-979b-d296ab29329c -md""" -# Summary -""" - -# ╔═╡ 5da42e9a-4318-48ea-9f43-4c1e1c97bceb -keyconceptsummary() - -# ╔═╡ 03692f4d-0daf-4dfc-a7ff-6b954326e4d0 -exercises() - - -# ╔═╡ 3a1d380e-df80-4727-9772-f199214cf05d -md""" -##### The Sum Rule (**) - -Derive the general sum rule, -```math -p(A + B) = p(A) + p(B) - p(A,B) -``` -from the elementary sum rule ``p(A) + p(\bar A) = 1`` and the sum and product rules. - -""" - -# ╔═╡ 99d9099f-4908-4bb3-8d59-da9cb69af04c -hint( - md""" - Here, you may make use of the (Boolean logic) fact that ``A + B = \overline {\bar A \bar B }``. - """ -) - -# ╔═╡ 3b1b0869-b815-4697-9dba-3c4b4cb5ac47 -hide_solution( -md""" -```math -\begin{align} -p\left( A + B \right) &\underset{\mathrm{bool}}{=} p\left( \overline {\bar A \bar B } \right) \\ - &\underset{\mathrm{sum}}{=} 1 - p\left( \bar{A} \bar{B} \right) \\ - &\underset{\mathrm{prod}}{=} 1 - p\left( \bar{A} |\bar{B} \right) p\left(\bar{B} \right) \\ - &\underset{\mathrm{sum}}{=} 1 - \left( 1 - p\left(A|\bar B \right) \right) \left( 1 - p\left( B \right) \right) \\ - &= p(B) + \left( {1 - p\left( B \right)} \right)p\left( {A|\bar B } \right) \\ - &\underset{\mathrm{prod}}{=} p(B) + \left( 1 - p\left( B \right) \right) p\left( \bar{B} |A \right) \frac{ p\left( A \right) }{ p\left(\bar{B}\right)} \\ - &\underset{\mathrm{sum}}{=} p(B) + p\left(\bar{B} |A \right) p\left( A \right) \\ - &\underset{\mathrm{sum}}{=} p(B) + \left( 1 - p\left( {B|A} \right) \right) p\left( A \right) \\ - &\underset{\mathrm{sum}}{=} p\left( A \right) + p(B) - p\left( A,B \right) -\end{align} -``` -Note that, aside from the first boolean rewrite, everything follows straight application of sum and product rules. - - -""") - -# ╔═╡ 5f377237-d9a5-4778-aa4d-1c6ce109b705 -md""" -##### Apples and Oranges - -Box 1 contains 8 apples and 4 oranges. Box 2 contains 10 apples and 2 oranges. Boxes are chosen with equal probability. You pick a box and then select a fruit from that box. -- (a) (*) What is the probability of choosing an apple? -- (b) (**) If an apple is chosen, what is the probability that it came from box 1? -""" - -# ╔═╡ 5613e9b7-ff0d-435a-9de6-aaf293ebf592 -hide_solution( -md""" -The following probabilities are given in the problem statement, -```math -\begin{align} -p(b_1) &= p(b_2) = 1/2 \\ -p(a|b_1) &= 8/12, \quad p(a|b_2) = 10/12 \\ -p(o|b_1) &= 4/12, \quad p(o|b_2) = 2/12 -\end{align} -``` -(a) -```math -p(a) = \sum_i p(a,b_i) = \sum_i p(a|b_i)p(b_i)=\frac{8}{12}\cdot\frac{1}{2} + \frac{10}{12}\cdot\frac{1}{2} = \frac{3}{4} -``` -(b) -```math -p(b_1|a) = \frac{p(a,b_1)}{p(a)} = \frac{p(a|b_1)p(b_1)}{p(a)} = \frac{\frac{8}{12}\cdot\frac{1}{2}}{\frac{3}{4}} = \frac{4}{9} -``` -""" -) - -# ╔═╡ fc3151f9-e143-4e31-b7b7-3f25b4fe9dab -md""" -##### What is a Random Signal? (*) -Is a speech signal a "probabilistic" (random) or a deterministic signal? -""" - -# ╔═╡ 66ebe33c-8360-4938-9b51-625e5bed176c -hide_solution( -md""" -That depends. The term “probabilistic” refers to a state-of-knowledge (or beliefs) about something—in this case, about the values of a speech signal. The key point is that the signal itself is neither inherently probabilistic nor deterministic; these labels describe our knowledge about it. - -If you had a perfect microphone and recorded the speech signal flawlessly at its source, you would know all its values exactly—no uncertainty—so you could call it deterministic. - -However, before making the recording, how would you represent your knowledge about the signal values you are going to measure? You face uncertainty, so the appropriate description is a probability distribution over all possible signal values. -""") - -# ╔═╡ 5b681e41-ad14-4c58-8ea0-4b6d85885c51 -md""" -##### Who Speaks the Truth? (***) -The inhabitants of an island tell the truth one-third of the time. They lie with probability ``2/3``. On an occasion, after one of them made a statement, you ask another person "was that statement true?" and he says "yes". What is the probability that the statement was indeed true? - -""" - -# ╔═╡ 91dd40f0-c373-48b3-b83b-6e8df2c43e5a -hide_solution( -md""" -We use variables ``S_1 \in \{\text{t},\text{f}\}`` and ``S_2 \in \{\text{y},\text{n}\}`` for statements 1 and 2 and shorthand "y", "n", "t" and "f" for "yes", "no", "true" and "false", respectively. The problem statement provides us with the following probabilities, -```math -\begin{align} -p(S_1=\text{t}) &= 1/3 \\ -p(S_1=\text{f}) &= 1 - p(S_1=\text{t}) = 2/3\\ -p(S_2=\text{y} | S_1=\text{t}) &= 1/3 \\ -p(S_2=\text{y} | S_1=\text{f}) &= 2/3 -\end{align} -``` -We are asked to compute ``p(S_1=\text{t} | S_2=\text{y})``. Use Bayes rule, -```math -\begin{align} -p(S_1=\text{t} | S_2=\text{y}) &= \frac{p(S_1=\text{t},S_2=\text{y})}{p(S_2=\text{y})} \\ -&= \frac{\overbrace{p(S_2=\text{y}|S_1=\text{t})p(S_1=\text{t})}^{\text{both speak the truth}}}{\underbrace{p(S_2=\text{y}|S_1=\text{t})p(S_1=\text{t})}_{\text{both speak the truth}}+\underbrace{p(S_2=\text{y}|S_1=\text{f})p(S_1=\text{f})}_{\text{both lie}}}\\ -&= \frac{\frac{1}{3}\cdot\frac{1}{3}}{\frac{1}{3}\cdot\frac{1}{3}+\frac{2}{3}\cdot\frac{2}{3}} = \frac{1}{5} -\end{align} -``` -""") - -# ╔═╡ a8d4a517-84a7-426e-a49e-482c5fd047ae -md""" -##### The Likelihood Function is a Function of What? (*) - -When considering the distribution ``p(D|\theta)``, is it more correct to speak about the likelihood of the model parameters ``\theta`` than about the likelihood of the observed data set ``D``. And why? - -""" - -# ╔═╡ d3b003c6-70ca-419f-a343-e35b266323f3 -hide_solution( -md""" -Yes, it’s more correct to speak about the likelihood of the model parameters, not of the observed data set. Once ``D`` has been observed, it is no longer a random variable; it’s just a fixed outcome. What varies is ``\theta``, so ``L(\theta) = p(D|\theta)`` is a function of the parameters, not of the data. - -Saying “likelihood of the data” is misleading because it confuses likelihood with the sampling distribution ``p(D|\theta)`` seen as a function of ``D`` (where ``\theta`` is fixed). The latter is a probability distribution over possible data sets before observing them. -""") - # ╔═╡ dd31ec7c-708d-4fd7-958d-f9887798a5bc md""" # Code @@ -1430,16 +1425,6 @@ result = (sensitivity * prevalence) / (sensitivity * prevalence + (1 - specifici # ╔═╡ 4a81342c-17c7-4eb9-933b-edb98df7b9c4 n(x; digits=2) = @sprintf("%.*f", digits, x) -# ╔═╡ 079157c9-5d97-4dbc-8c47-afa8b661db06 -let - θ_str = n(θ) - @mdx """ - ```math - p(y|\\theta=$(n(θ))) = \\begin{cases} $(n(1-θ)) & \\text{if }y=0 \\\\ $(n(θ)) & \\text{if } y=1 \\end{cases} - ``` - """ -end - # ╔═╡ 2156f96e-eebe-4190-8ce9-c76825c6da71 if show_disease_diagnosis_solution @mdx(""" @@ -1466,6 +1451,16 @@ Many people have trouble distinguishing ``p(A|B)`` from ``p(B|A)`` in their head """) end +# ╔═╡ 079157c9-5d97-4dbc-8c47-afa8b661db06 +let + θ_str = n(θ) + @mdx """ + ```math + p(y|\\theta=$(n(θ))) = \\begin{cases} $(n(1-θ)) & \\text{if }y=0 \\\\ $(n(θ)) & \\text{if } y=1 \\end{cases} + ``` + """ +end + # ╔═╡ 00000000-0000-0000-0000-000000000001 PLUTO_PROJECT_TOML_CONTENTS = """ [deps] @@ -2761,7 +2756,6 @@ version = "1.9.2+0" # ╔═╡ Cell order: # ╟─3e17df5e-d294-11ef-38c7-f573724871d8 # ╟─bcb4be20-0439-4809-a166-8c50b6b9206b -# ╟─3e1803d0-d294-11ef-0304-df2b9b698cd1 # ╟─61713e1c-8e37-45d9-9f58-c3db69e15b66 # ╟─3e185ab0-d294-11ef-3f7d-9bd465518274 # ╟─840ab4dc-0d2e-4bf8-acc7-5f1ee2b0dcaf @@ -2772,6 +2766,51 @@ version = "1.9.2+0" # ╟─3e18d2ea-d294-11ef-35e9-2332dd31dbf0 # ╟─4abbb3de-3b21-4c31-b015-e16c466a20aa # ╟─3e18e4bc-d294-11ef-38bc-cb97cb4e0963 +# ╟─3e1ab104-d294-11ef-1a98-412946949fba +# ╟─fea8ae4c-8ef9-4b74-ad13-1314afef97de +# ╟─3e1b4b1c-d294-11ef-0423-9152887cc403 +# ╟─3e1b5c9c-d294-11ef-137f-d75b3731eae4 +# ╟─3e1b7d14-d294-11ef-0d10-1148a928dd57 +# ╟─5377c5a4-77c4-4fa7-9f84-0c511e3bf708 +# ╟─3e1b8bf4-d294-11ef-04cc-6364e46fdd64 +# ╟─3e1b9ba8-d294-11ef-18f2-db8eed3d87d0 +# ╟─3e1babca-d294-11ef-37c1-cd821a6488b2 +# ╟─3e1bcb00-d294-11ef-2795-bd225bd00496 +# ╟─3e1bdd02-d294-11ef-19e8-2f44eccf58af +# ╟─3e1bf116-d294-11ef-148b-f7a1ca3f3bad +# ╟─16c2eb59-16b8-4347-9aab-6e4b99016c79 +# ╟─3e1d33c8-d294-11ef-0a08-bdc419949925 +# ╟─3e1b05ee-d294-11ef-33de-efed64d01c0d +# ╟─b176ceae-884e-4460-9f66-020c1ac447f1 +# ╟─ab223dea-8ba8-4d30-94f4-72c8e070aadf +# ╟─922770f4-ddc8-4089-b378-f14088276b43 +# ╟─3e1de32c-d294-11ef-1f63-f190c8361404 +# ╟─4c639e65-e06b-4c5e-b6e7-aabed6b6c0b4 +# ╟─ff9142ba-3a85-48cf-8b78-07e0b554e280 +# ╟─3e027ede-4ac1-4521-a026-dc00bfca4adf +# ╟─3e1e2b96-d294-11ef-3a68-fdc78232142e +# ╟─727dc817-0284-4c0f-9a92-21dcbea50807 +# ╟─fae6f2ce-ac8f-4ea6-b2cf-38b30a7e20d4 +# ╟─ea803646-b18e-4f13-ab6a-fc5080d44e92 +# ╟─3e1d6d00-d294-11ef-1081-e11b8397eb91 +# ╟─9b85a92b-c56b-48a3-97c1-6b1882e33a22 +# ╟─2156f96e-eebe-4190-8ce9-c76825c6da71 +# ╟─9b92fc89-2036-4525-979b-d296ab29329c +# ╟─5da42e9a-4318-48ea-9f43-4c1e1c97bceb +# ╟─03692f4d-0daf-4dfc-a7ff-6b954326e4d0 +# ╟─3a1d380e-df80-4727-9772-f199214cf05d +# ╟─99d9099f-4908-4bb3-8d59-da9cb69af04c +# ╟─3b1b0869-b815-4697-9dba-3c4b4cb5ac47 +# ╟─5f377237-d9a5-4778-aa4d-1c6ce109b705 +# ╟─5613e9b7-ff0d-435a-9de6-aaf293ebf592 +# ╟─fc3151f9-e143-4e31-b7b7-3f25b4fe9dab +# ╟─66ebe33c-8360-4938-9b51-625e5bed176c +# ╟─5b681e41-ad14-4c58-8ea0-4b6d85885c51 +# ╟─91dd40f0-c373-48b3-b83b-6e8df2c43e5a +# ╟─a8d4a517-84a7-426e-a49e-482c5fd047ae +# ╟─d3b003c6-70ca-419f-a343-e35b266323f3 +# ╟─be66b697-f920-4361-9ff2-b12cc50ae8c9 +# ╟─3e1803d0-d294-11ef-0304-df2b9b698cd1 # ╟─3e18f18c-d294-11ef-33e4-b7f9495e0508 # ╟─3e1906ea-d294-11ef-236e-c966a9474170 # ╟─3e191b6c-d294-11ef-3174-d1b4b36e252b @@ -2787,20 +2826,6 @@ version = "1.9.2+0" # ╟─3e1a69f4-d294-11ef-103e-efc47025fb8f # ╟─3e1a7c8e-d294-11ef-1f97-55e608d49141 # ╟─3e1a8eca-d294-11ef-1ef0-c15b24d05990 -# ╟─3e1ab104-d294-11ef-1a98-412946949fba -# ╟─fea8ae4c-8ef9-4b74-ad13-1314afef97de -# ╟─3e1b4b1c-d294-11ef-0423-9152887cc403 -# ╟─3e1b5c9c-d294-11ef-137f-d75b3731eae4 -# ╟─3e1b7d14-d294-11ef-0d10-1148a928dd57 -# ╟─5377c5a4-77c4-4fa7-9f84-0c511e3bf708 -# ╟─3e1b8bf4-d294-11ef-04cc-6364e46fdd64 -# ╟─3e1b9ba8-d294-11ef-18f2-db8eed3d87d0 -# ╟─3e1babca-d294-11ef-37c1-cd821a6488b2 -# ╟─3e1bba8e-d294-11ef-1f61-295af16078ce -# ╟─3e1bcb00-d294-11ef-2795-bd225bd00496 -# ╟─3e1bdd02-d294-11ef-19e8-2f44eccf58af -# ╟─3e1bf116-d294-11ef-148b-f7a1ca3f3bad -# ╟─16c2eb59-16b8-4347-9aab-6e4b99016c79 # ╟─3e1bffec-d294-11ef-2a49-9ff0f6331add # ╟─3e1c0e80-d294-11ef-0d19-375e01988f16 # ╟─3e1c1e3e-d294-11ef-0955-bdf9d0ba3c53 @@ -2810,9 +2835,6 @@ version = "1.9.2+0" # ╟─3e1c70be-d294-11ef-14ed-0d46515541c5 # ╟─3e1c806a-d294-11ef-1fad-17e5625279f7 # ╟─3e1c9184-d294-11ef-3e35-5393d97fbc44 -# ╟─3e1d33c8-d294-11ef-0a08-bdc419949925 -# ╟─3e1b05ee-d294-11ef-33de-efed64d01c0d -# ╟─b176ceae-884e-4460-9f66-020c1ac447f1 # ╟─e3157dc0-5a64-4479-a37a-40fe25cccc07 # ╟─d93f73d4-2783-4777-b0ce-cdc0444cb300 # ╟─079157c9-5d97-4dbc-8c47-afa8b661db06 @@ -2821,19 +2843,6 @@ version = "1.9.2+0" # ╟─b7445b9b-7fbb-4560-b947-a23af0fcf101 # ╟─3e1d20e0-d294-11ef-2044-e1fe6590a600 # ╠═fc733d61-fd0f-4a13-9afc-4505ac0253df -# ╟─ab223dea-8ba8-4d30-94f4-72c8e070aadf -# ╟─922770f4-ddc8-4089-b378-f14088276b43 -# ╟─3e1de32c-d294-11ef-1f63-f190c8361404 -# ╟─4c639e65-e06b-4c5e-b6e7-aabed6b6c0b4 -# ╟─ff9142ba-3a85-48cf-8b78-07e0b554e280 -# ╟─3e027ede-4ac1-4521-a026-dc00bfca4adf -# ╟─3e1e2b96-d294-11ef-3a68-fdc78232142e -# ╟─727dc817-0284-4c0f-9a92-21dcbea50807 -# ╟─fae6f2ce-ac8f-4ea6-b2cf-38b30a7e20d4 -# ╟─ea803646-b18e-4f13-ab6a-fc5080d44e92 -# ╟─3e1d6d00-d294-11ef-1081-e11b8397eb91 -# ╟─9b85a92b-c56b-48a3-97c1-6b1882e33a22 -# ╟─2156f96e-eebe-4190-8ce9-c76825c6da71 # ╟─ef264651-854e-4374-8ea8-5476c85150c4 # ╟─3e1e4dda-d294-11ef-33b7-4bbe3300ca22 # ╟─3e1e5a5a-d294-11ef-2fdf-efee4eb1a0f2 @@ -2869,21 +2878,7 @@ version = "1.9.2+0" # ╟─db73766d-643c-41d7-a1eb-f376c657f860 # ╟─3e1fb370-d294-11ef-1fb6-63a41a024691 # ╟─317707a3-9ef1-4c67-b451-6adcfcff50f0 -# ╟─9b92fc89-2036-4525-979b-d296ab29329c -# ╟─5da42e9a-4318-48ea-9f43-4c1e1c97bceb -# ╟─03692f4d-0daf-4dfc-a7ff-6b954326e4d0 -# ╟─3a1d380e-df80-4727-9772-f199214cf05d -# ╟─99d9099f-4908-4bb3-8d59-da9cb69af04c -# ╟─3b1b0869-b815-4697-9dba-3c4b4cb5ac47 -# ╟─5f377237-d9a5-4778-aa4d-1c6ce109b705 -# ╟─5613e9b7-ff0d-435a-9de6-aaf293ebf592 -# ╟─fc3151f9-e143-4e31-b7b7-3f25b4fe9dab -# ╟─66ebe33c-8360-4938-9b51-625e5bed176c -# ╟─5b681e41-ad14-4c58-8ea0-4b6d85885c51 -# ╟─91dd40f0-c373-48b3-b83b-6e8df2c43e5a -# ╟─a8d4a517-84a7-426e-a49e-482c5fd047ae -# ╟─d3b003c6-70ca-419f-a343-e35b266323f3 -# ╟─dd31ec7c-708d-4fd7-958d-f9887798a5bc +# ╠═dd31ec7c-708d-4fd7-958d-f9887798a5bc # ╠═eeb9a1f5-b857-4843-920b-2e4a9656f66b # ╠═5394e37c-ae00-4042-8ada-3bbf32fbca9e # ╠═b305a905-06c2-4a15-8042-72ef6375720f diff --git a/mlss/archive/Regression.jl b/mlss/archive/Regression.jl index ec8e3ffa..908a0bed 100644 --- a/mlss/archive/Regression.jl +++ b/mlss/archive/Regression.jl @@ -39,36 +39,6 @@ title("Regression") # ╔═╡ 66998cd5-78d6-4b22-a9c9-886436cba4dd PlutoUI.TableOfContents() -# ╔═╡ 234b8c8e-d294-11ef-296a-3b38564babc4 -md""" -## Preliminaries - -##### Goal - -* Introduction to Bayesian (Linear) Regression - -##### Materials - -* Mandatory - - * These lecture notes - -* Optional - - * [Bishop PRML book](https://www.microsoft.com/en-us/research/wp-content/uploads/2006/01/Bishop-Pattern-Recognition-and-Machine-Learning-2006.pdf), pp. 152-158 - - * Matrix Calculus - * In this and forthcoming lectures, we will make use of some elementary matrix calculus. Please see the [Formula Cheatsheet](https://github.com/bmlip/course/blob/main/assets/files/5SSD0_formula_sheet.pdf) for formulas that will be made available to you at the written exam. - - * [RxInfer Bayesian Linear Regression example](https://examples.rxinfer.com/categories/basic_examples/bayesian_linear_regression/) - * A tutorial on Bayesian linear regression with RxInfer. - - * Jaynes (1990), [Straight Line Fitting - A Bayesian Solution](https://github.com/bmlip/course/blob/main/assets/files/Jaynes-1990-straight-line-fitting-a-Bayesian-solution.pdf) - * A fully Bayesian solution on straight line fitting with uncertainties in both ``x`` and ``y`` coordinates. - - -""" - # ╔═╡ 234ba8c2-d294-11ef-36f6-b1f61f65557a @@ -401,12 +371,69 @@ md""" See the [Mini about Basis Functions](https://bmlip.github.io/course/minis/Basis%20Functions.html) to learn more! """ -# ╔═╡ 142f4700-ccf4-4019-b3a9-57035c458276 -σ_basis² = 0.01; +# ╔═╡ 22e76656-b9f4-463e-9bf8-bd383e92948b +const Layout = PlutoUI.ExperimentalLayout + +# ╔═╡ 9fd4a9b4-3296-4fe3-931f-17744bc4df81 +Layout.vbox([ + secret_function_bond, + N_bond, + σ_noise_bond, +]) + +# ╔═╡ 338ee8e9-b786-48e1-b084-a0e8a6d12118 +baseplot(args...; kwargs...) = plot(args...; size=(650,400), xlim=(-0.0, 1.0), ylim=(-1.2,1.2), kwargs...) + +# ╔═╡ 0e9435fc-3206-4249-b5ff-42cc35c98d47 +function plot_data!(D) + plot!(; legend=:bottomleft) + plot!(secret_function; + label="True function", + color=3, + lw=3, + linestyle=:dash, + ) + scatter!( + D; + label="Observations", + color=1, + # markerstrokewidth=0, + ) +end + +# ╔═╡ 88a2bd82-6663-48cc-a535-5b3e47d814a9 +const deterministic_randomness = MersenneTwister + +# ╔═╡ 68141653-e444-4e29-bbec-4cd7359cb84c +σ_data_noise² = σ_data_noise^2 + +# ╔═╡ 72fcb6a3-36ee-4840-bdc3-ddb743e5c149 +D = let + xs = rand(deterministic_randomness(19), Uniform(0,1), N) + + ys_exact = secret_function.(xs) + + rng = deterministic_randomness(37) + ys = [ + rand(rng, Normal(y, sqrt(σ_data_noise²))) + for y in ys_exact + ] + + collect(zip(xs, ys)) +end + +# ╔═╡ 5d48e25f-9a98-43ce-8f23-d9ab28f69996 +let + baseplot() + plot_data!(D) +end # ╔═╡ 70ca3a3f-ee1c-4f3d-9d77-bf55e8e808c1 μ_basis = range(0.0, 1.0; length=10); +# ╔═╡ 142f4700-ccf4-4019-b3a9-57035c458276 +σ_basis² = 0.01; + # ╔═╡ 3a3b7ff2-68aa-411c-b7fb-c6cd00d0dd7b ϕ(μ, x) = exp(-(x - μ)^2 / σ_basis²); @@ -417,11 +444,38 @@ function f(w, x) end end; -# ╔═╡ 4d2be102-8849-4b8d-9962-8b45099ab8f2 -md""" -#### Bayesian inference -We have a closed-form solution for the posterior: -""" +# ╔═╡ 8a2730b8-3262-48cc-81f0-777cf85b9836 +# This is called the "design matrix" +Φ = [ + ϕ(μ, datum[1]) + for datum in D, μ in μ_basis +]; + +# ╔═╡ 98d729ab-79f8-4a0f-9db5-387f488fc19d +weights_posterior = MvNormalCanon( + # Posterior potential vector + Φ' * last.(D) / σ_data_noise², + # Posterior precision matrix (inverse covariance) + Φ' * Φ / σ_data_noise² + I / σ_prior² +); + +# ╔═╡ f9a5c91e-12be-4e8b-930d-74e46e39ea58 +let + baseplot() + if true + for i in 1:40 + w = rand(weights_posterior) + plot!( + x -> f(w, x); + opacity=.3, + color=2, + label=i==1 ? "Posterior samples" : nothing, + ) + end + end + + plot_data!(D) +end # ╔═╡ ec0ccf94-e12e-422d-b4d2-dcb933453146 md""" @@ -638,19 +692,6 @@ md""" _Reading this lecture online? Click **"View code"** in the top right to read the implementation of this visualisation._ """ -# ╔═╡ 234f5d32-d294-11ef-279f-f331396e47ad -md""" - -## Uncertainty About Inputs? - -In this lesson, we focused on modelling the map from given inputs ``x`` to uncertain outputs ``y``, or more formally, on the distribution ``p(y|x)``. - -What if you want to fit the best curve through a data set ``\{(x_1,y_1),\dotsc,(x_N,y_N)\}`` where both variables ``x_n`` and ``y_n`` are subject to errors? In other words, we must now also fit a model ``p(x)`` for the inputs, leading to a generative model ``p(y,x) = p(y|x) p(x)``. - -While this is a very common problem that occurs throughout the sciences, a proper solution to this problem is still hardly covered in statistics textbooks. Edwin T. Jaynes (author of the brilliant book [Probability Theory: The Logic of Science](https://bayes.wustl.edu/etj/prob/book.pdf)) discusses a fully Bayesian solution in his 1990 paper on [Straight Line Fitting - A Bayesian Solution](https://github.com/bmlip/course/blob/main/assets/files/Jaynes-1990-straight-line-fitting-a-Bayesian-solution.pdf). (Optional reading). - -""" - # ╔═╡ 9577225a-9ce3-4cf2-ac63-499ae8e905bd md""" # Summary @@ -807,96 +848,6 @@ md""" # Code """ -# ╔═╡ 22e76656-b9f4-463e-9bf8-bd383e92948b -const Layout = PlutoUI.ExperimentalLayout - -# ╔═╡ 9fd4a9b4-3296-4fe3-931f-17744bc4df81 -Layout.vbox([ - secret_function_bond, - N_bond, - σ_noise_bond, -]) - -# ╔═╡ 338ee8e9-b786-48e1-b084-a0e8a6d12118 -baseplot(args...; kwargs...) = plot(args...; size=(650,400), xlim=(-0.0, 1.0), ylim=(-1.2,1.2), kwargs...) - -# ╔═╡ 0e9435fc-3206-4249-b5ff-42cc35c98d47 -function plot_data!(D) - plot!(; legend=:bottomleft) - plot!(secret_function; - label="True function", - color=3, - lw=3, - linestyle=:dash, - ) - scatter!( - D; - label="Observations", - color=1, - # markerstrokewidth=0, - ) -end - -# ╔═╡ 88a2bd82-6663-48cc-a535-5b3e47d814a9 -const deterministic_randomness = MersenneTwister - -# ╔═╡ 68141653-e444-4e29-bbec-4cd7359cb84c -σ_data_noise² = σ_data_noise^2 - -# ╔═╡ 72fcb6a3-36ee-4840-bdc3-ddb743e5c149 -D = let - xs = rand(deterministic_randomness(19), Uniform(0,1), N) - - ys_exact = secret_function.(xs) - - rng = deterministic_randomness(37) - ys = [ - rand(rng, Normal(y, sqrt(σ_data_noise²))) - for y in ys_exact - ] - - collect(zip(xs, ys)) -end - -# ╔═╡ 5d48e25f-9a98-43ce-8f23-d9ab28f69996 -let - baseplot() - plot_data!(D) -end - -# ╔═╡ 8a2730b8-3262-48cc-81f0-777cf85b9836 -# This is called the "design matrix" -Φ = [ - ϕ(μ, datum[1]) - for datum in D, μ in μ_basis -]; - -# ╔═╡ 98d729ab-79f8-4a0f-9db5-387f488fc19d -weights_posterior = MvNormalCanon( - # Posterior potential vector - Φ' * last.(D) / σ_data_noise², - # Posterior precision matrix (inverse covariance) - Φ' * Φ / σ_data_noise² + I / σ_prior² -); - -# ╔═╡ f9a5c91e-12be-4e8b-930d-74e46e39ea58 -let - baseplot() - if true - for i in 1:40 - w = rand(weights_posterior) - plot!( - x -> f(w, x); - opacity=.3, - color=2, - label=i==1 ? "Posterior samples" : nothing, - ) - end - end - - plot_data!(D) -end - # ╔═╡ 00000000-0000-0000-0000-000000000001 PLUTO_PROJECT_TOML_CONTENTS = """ [deps] @@ -2179,7 +2130,6 @@ version = "1.9.2+0" # ╔═╡ Cell order: # ╟─234b77a8-d294-11ef-15d5-ff54ed5bec1e # ╟─66998cd5-78d6-4b22-a9c9-886436cba4dd -# ╟─234b8c8e-d294-11ef-296a-3b38564babc4 # ╟─234ba8c2-d294-11ef-36f6-b1f61f65557a # ╟─e248df9b-0c48-4803-8d1e-b466ab07692e # ╟─f1bf64f6-09f9-45a3-acd1-7975ab9e79fc @@ -2204,19 +2154,23 @@ version = "1.9.2+0" # ╟─234d6dd8-d294-11ef-3abf-8d6cb00b1907 # ╟─ab3baeb4-51d7-4f50-9e06-ed00bb783ebd # ╟─fb113692-f00c-4b48-85cc-d7bba88c7099 -# ╟─f600c228-e048-42aa-b79a-60592b367dec -# ╟─c0c57aa6-155a-49a9-9ed2-d568de1b5be2 -# ╟─9fd4a9b4-3296-4fe3-931f-17744bc4df81 -# ╟─b48b93c3-1ff2-4be0-8fad-181035f3e50e -# ╟─3fe01c67-6f95-4f6d-8c7f-5a389272ff65 -# ╟─f9a5c91e-12be-4e8b-930d-74e46e39ea58 -# ╟─018b6c7b-36bc-4867-a058-3802b43fd1eb -# ╟─90cb881a-7b5d-44e3-a7d1-bb93bef4a82b -# ╠═142f4700-ccf4-4019-b3a9-57035c458276 +# ╠═f600c228-e048-42aa-b79a-60592b367dec +# ╠═c0c57aa6-155a-49a9-9ed2-d568de1b5be2 +# ╠═9fd4a9b4-3296-4fe3-931f-17744bc4df81 +# ╠═b48b93c3-1ff2-4be0-8fad-181035f3e50e +# ╠═3fe01c67-6f95-4f6d-8c7f-5a389272ff65 +# ╠═f9a5c91e-12be-4e8b-930d-74e46e39ea58 +# ╠═018b6c7b-36bc-4867-a058-3802b43fd1eb +# ╠═90cb881a-7b5d-44e3-a7d1-bb93bef4a82b +# ╠═22e76656-b9f4-463e-9bf8-bd383e92948b +# ╠═338ee8e9-b786-48e1-b084-a0e8a6d12118 +# ╠═0e9435fc-3206-4249-b5ff-42cc35c98d47 +# ╠═88a2bd82-6663-48cc-a535-5b3e47d814a9 +# ╠═68141653-e444-4e29-bbec-4cd7359cb84c # ╠═70ca3a3f-ee1c-4f3d-9d77-bf55e8e808c1 +# ╠═142f4700-ccf4-4019-b3a9-57035c458276 # ╠═3a3b7ff2-68aa-411c-b7fb-c6cd00d0dd7b # ╠═290bc994-d0f9-4af3-bd63-78de1640c85c -# ╟─4d2be102-8849-4b8d-9962-8b45099ab8f2 # ╠═98d729ab-79f8-4a0f-9db5-387f488fc19d # ╠═8a2730b8-3262-48cc-81f0-777cf85b9836 # ╟─ec0ccf94-e12e-422d-b4d2-dcb933453146 @@ -2239,7 +2193,6 @@ version = "1.9.2+0" # ╟─b6443a13-9301-4559-a5c3-396bae2a27b9 # ╟─234ef126-d294-11ef-17a9-3da87a7e7d0a # ╟─e9804f92-29b0-4463-bf37-872183061ee2 -# ╟─234f5d32-d294-11ef-279f-f331396e47ad # ╟─9577225a-9ce3-4cf2-ac63-499ae8e905bd # ╟─70ddbde1-5de4-4bc0-ace0-9414ac616888 # ╟─8e2b2c1d-81f3-4283-ae2e-d8b3e9c201b3 @@ -2252,10 +2205,5 @@ version = "1.9.2+0" # ╠═f8c69b91-4415-454e-a50d-c4a37ada89d1 # ╠═33ca4c67-d96f-457f-bc19-171f4b4b03c6 # ╠═3ff2bd04-1490-4be6-8b26-b82d1902bb07 -# ╟─22e76656-b9f4-463e-9bf8-bd383e92948b -# ╟─338ee8e9-b786-48e1-b084-a0e8a6d12118 -# ╟─0e9435fc-3206-4249-b5ff-42cc35c98d47 -# ╟─88a2bd82-6663-48cc-a535-5b3e47d814a9 -# ╟─68141653-e444-4e29-bbec-4cd7359cb84c # ╟─00000000-0000-0000-0000-000000000001 # ╟─00000000-0000-0000-0000-000000000002