Submind YouTube summaries
Thumbnail for Probabilistic ML - 23 - Variational Inference

Probabilistic ML - 23 - Variational Inference

Watch on YouTube

Video summary

Variational Inference represents a significant evolution in probabilistic machine learning, building upon classic algorithms like K-means and Expectation-Maximization while aiming to induce structural efficiency for faster convergence. Rooted in early 20th-century physics and the work of figures such as Richard Feynman and Edward Teller, this method addresses scenarios where computing exact posteriors or Evidence Lower Bounds is computationally intractable. To overcome these challenges, Variational Inference imposes a factorization assumption on the approximating distribution, effectively treating variables as independent through the Mean Field approximation. This approach transforms complex optimization problems into an iterative process where each variable's approximate distribution is updated based on the expected energy contributions from all other variables, creating a self-consistent loop that simplifies the mathematical landscape. The historical workflow of developing these models involved manually deriving algorithms by defining joint distributions over data, latent cluster assignments, and parameters using conjugate priors like Dirichlet and Gaussian-Inverse-Gamma distributions. By computing expected values under the mean-field factorization assumption, the mathematics naturally induced a specific structure where latent variables became independent and parameter posteriors emerged as conjugate distributions. Although this manual derivation process was tedious, it yielded highly efficient and bug-free algorithms capable of automatically discovering the optimal number of clusters in a dataset. This technique served as a precursor to modern attention mechanisms, demonstrating how inducing factorization could lead to powerful, automatic model configurations without relying on brute-force computation. Over time, the machine learning community shifted away from these rigorous derivations in favor of simpler gradient descent methods commonly used in deep learning, inadvertently losing the valuable "induced factorization" structure that made earlier algorithms so effective. The speaker argues that this transition resulted in a loss of functionality, as the concept was only later rediscovered within deep learning architectures under the guise of attention mechanisms. By abandoning the structured approaches of Variational Inference for generic optimization, the field missed opportunities to leverage these inherent mathematical properties for more efficient and interpretable models. Ultimately, the demonstration of the derived Variational Inference algorithm shows how it iteratively updates responsibilities and parameters to converge on an optimal model configuration, bridging the gap between classical probabilistic methods and modern neural network architectures. This journey highlights the importance of understanding the underlying mathematical structures that drive learning algorithms, rather than treating them as black boxes. The legacy of these early derivations persists today, reminding practitioners that the elegant solutions found in physics and statistics can still inform and improve contemporary machine learning practices, ensuring that future developments continue to benefit from induced structure and efficient convergence.
Read the full video transcript
[Music] On Tuesday we spoke about K means and then about EM two classic algorithms. K means is 100 years old. EM 70 years old. We're now making our way to 40 years, 50 years old, 20 years old and then uh next Tuesday we'll be at five years old. The idea behind these fancy mathematical ideas that I have had to compress in these difficult equations is on the high level that inducing structure in a probabilistic model can be used to create algorithmic efficiency to construct programming structures in your code that allow much faster convergence of much more expressive models. So the first instance of this from 100 years ago hubush dynos is k means which I had up there. This is what I spoke about on Tuesday. It's an algorithm to fit a mixture model specifically a gausian mixture model where the gaussian clusters have no coariance or standard unit coariance but just a mean that gets fitted. And we found that this algorithm which has an iterative structure alternating between assigning um points to clusters and then fitting clusters to the points assigned to them and then repeating that that amounts to maximizing a log likelihood of such a gausian mixture model. actually maximizing the log likelihood of the parameters of the model which are the assignment of uh points to clusters and the parameters of the clusters. In the case of K means by the way this lock probability for the data given the clusters is just a gausian with a mean and a standard unit coariance. But this could be something else as well for other data types. And then we made a step and said, "Oh, actually there's a refinement of this algorithm that doesn't maximize the complete data log likelihood." So the likelihood involving both theta and zed, but it maximizes the expected complete data log likelihood. And that's called the EM algorithm. So we maximize this expression but under an expectation over the latent quantity. And then it turned out that the right thing to take the expectation over is the posterior distribution over zed because then at the time when we set this distribution this posterior to the to the postivia at theta hat and then optimize with respect to this theta we are actually closing a g a b a the gap between a lower bound on the thing we're trying to maximize the log evidence for x given the parameters having marginalized over zed and the so-called elbow the evidence lower bound which is actually this expression here on the left and then because we close this gap and then raise it we are uh guaranteed to increase the marginal lo likelihood of the model in every step and then therefore get a good model in the end. So that is an algorithm. So this you know um I've introduced it as as having been invented by Hamato Hartley and um and Dempster and others in the 50s. What I want to talk about today is the the direct continuation of this idea in in statistics at least throughout the second half of this of the 20th century. And um what I want to highlight is also its sort of story of how it influenced the development of the Nassen field that we now call machine learning. So let's to to do that let me briefly remind you of how the EM algorithm worked. So um our goal is we we someone has given us a model. We want this model involves a probability for some observations given some parameters. We'd like to know what a good what what good parameters are for this model. We want to do maximum likelihood inference. So to do that we invent some new variables called latent variables zed that we think somehow um simplifies the computation um and there is usually some ingenuity involved in coming up with those then that means we want to maximize this expression. This is called the complete data log likelihood because there's now this new object that we've invented. Um and the way we're going to do this is that we first compute a postivia over these latent variables. Hopefully we've chosen them such that we can actually do that. Um then given a particular current estimate for the parameters theta, we set a current estimate for the distribution over zed as the posterivia given those parameters. That means that when we look at this object called the evidence lower bound, we're effectively making this a tight lower bound on this by setting the K divergence between our approximation and this particular posterior to zero. And then we raise this expression as a function of theta but only this theta. So the theta in Q is kept constant. Um, and that's often possible analytically. And then we're guaranteed to raise the marginal likelihood. The first thing I want to point out is that I said maximum likelihood inference and this is a probabilistic learning course. So you probably should think ah but didn't shouldn't we do postivia rather than likelihood because likelihoods can be dangerous. That's completely true. It might be better to do this in a postivia sense rather than in a likelihood sense. But actually doing that is really straightforward. So if you want to do maximum apostroori inference rather than maximum likelihood inference, the only thing you need to change on on this slide even is just this red bit. So you just add a prior for the parameters theta. Then you're effectively maximizing a posterior rather than a likelihood. And um this is doing this is very easy because you can see that this theta does not show up in here. So this integral is actually just I mean if you think about what what what we add here then we just that there's just a term plus integral over q of z p theta dz p theta doesn't depend on zed. So the integral commutes out we just have an p of theta a log p of theta* integral over q of zed which is one. So just log p of theta so we just add a regularizer basically to the optimization. Then you keep this running until it converges. Now you have a maximum posterior estimate. And I motivated all of this with this f um like uh well sad story about um brilliant mathematicians in well central eastern Europe in the 1930s. And um actually that the path of these people now leads us as science tends to be because science tends to actually happen in very small communities to across time and space to a different place where the mathematical story of the idea actually continues. So this is a picture from um um a colloquium. So today we would say like a weekly talk um at the Los Alamos um laboratory during the Manhattan project in 1946. So the second world war had already ended. The nuclear bomb on Hiroshima Nagasaki had already been dropped. And you see here in the um front row I actually made some notes so I'm not going to tell you anything wrong. Um the uh in the first row Norris Bradbury he was the director of of uh Los Alamos. The this person you probably don't know this is a noble laureate Enrico Farmy he's uh probably the father of the nuclear bomb um as such. Does anyone reci recognize anyone else in the picture? Maybe the guy with the black jacket. Yes. Yes, this is J. Robert Oenheimer, now Hollywood famous. In this picture, he's 42 years old. Um, clearly still a young man. And, um, they are actually watching a talk on what they called the super bomb. Now we call them H bombs. Um, hydrogen fusion bombs. Probably being that talk is probably given by Edward Teller who was also an immigr. Uh he studied in um he was Hungarian, studied in Germany, immigrated, he was Jewish, he had to immigrate in 1934 having just finished his degree with his wife Mitsy Teller. Um and then he got to work quickly at Los Alamos and um sort of had this crazy idea that it should be possible to build a fusion bomb. He worked on this together with a colleague who we remember and met already in Ruf or Rof or Lmberg last lecture. Stanislafulam who was a part of the Scottish cafe. They uh together did computations using early Monte Carlo methods for which Mitsy teller also invented some part of the algorithm essentially Metropolis Hastings but they did the computations with an analog computer called the ENIAC. was like a little machine that you had to move around on a drawing to simulate the path of particles being scattered inside of a geometry of a nuclear bomb. They tried to find the geometry such that fusion would ignite and they would get very high um concentration of um of particles. So there is Monte Carlo methods being used well for impact let's say because they created I mean overall this was the outcome of all of this but in this picture so this was the very first hydrogen bomb shot this is IV Mike um in this picture there's another person a young guy a group leader today we would say an independent group leader at one of our institutions right next to Oppenheimer does anyone guess who this is Richard Fineman, later Noble laurate and he made his name with work on what we now call variational inference in physics though he already gave at this point he's 28 years old he already has his PhD he's basically a posttock um he just gave a talk recently like a few years ago to people like Albert Einstein and John Fon Noman who were quite intrigued by his work he was supposedly an absolutely brilliant theorist and he's extending work that was introduced as the under the name of the calculus of variations my uh various well old mathematicians from a long time ago so some names associated with this are also Isaac Newton very early version then leonard oiler got interested but lrange actually provided the first formal definition of this idea of the so-called calculus of variations that oiler then directly took over the idea of the calculus of variations very vaguely spoken without showing you any math but we'll go into the math quickly enough is that just like we can do derivatives on functions and then follow derivatives and gradients in multivariate calculus to find modes of functions right so to find to maximize the value of a function it should also be possible to construct and the equivalent of a derivative in spaces of functions rather than in real valued vector spaces in spaces of functions to maximize what is the corresponding concept to maximize a functional. So a function of a function and for that we need something like a derivative in function space that we can then follow to maximize such a functional an expression and LR was like totally excited about this idea and fineman got excited about it as well because it was a crucial ingredient for the theory he was working on quantum field theory where the expressions you want to maximize for are functionals of incoming probability distributions over particles. There is by the way just before we move on just to not leave it out also a corresponding equivalent to um well sort of the idea that you also know from lagrange multipliers. So finding directions along which a functional is invariant right? So just as just as you can follow a gradient to find a mode, you can also think about spaces identified by zero gradient which are equipotential lines. So spaces along which something is invariant and does not change. This is called um the the theory of of uh invariances in variations and it's associated with um the work of emintor who uh like was probably the like most famous mathematician of all times which also made ma major contributions to physics identifying conservation laws rather than optimal choices. But today we're going to focus on the optimization problem on the the idea of how to how to maximize a functional. Our elbow is such a functional. It takes in a distribution Q, right? That's a function of Z. And then we want to find the function Q which makes this expression maximal. So this expression is clear. Oh, I forgot the differential. So it's a derivative in integral over dz, right? So it's something that takes in the input over zed and marginalizes it out, integrates over it. And now we have a number or actually we have a function of x that we can try to maximize for a particular x. So what does that mean? So this is the this is the physics motivation. I'll talk more about that in a second. Where does this like what does this come in? What does this have to do with machine learning? What does this have to do with with like with inference and approximation? Well, machine learning is the continuation of physics into everything else. It's the application of the of mechanisms from science into every other aspect of human activity. That's maybe why it's so exciting and why it's so influential at the moment. So to try let me try and generalize the ideas both from the previous lecture and from this physics background in a few slides. So let's and for a moment let me try to be more concrete. So let's say we're trying to do inference in a model that assigns probability to some data points X. So X is the things we get to see and our model will involve some parameters theta and it involves some latent quantities which we need to describe what's going on. We could think of them as nuisance variables or maybe we think of them as the things we're particularly interested in the variables of the model zed. So in the algorithms we've encountered so far in K means and EM we had to make very specific choices for zed. So the the the name of the game was let's find the right latent variables such that this posterior over zed given x becomes tractable and we can just write it out for a particular value of theta. When we can do that then we our our elbow h well then we need two things actually. We need to be able to get the the posterior attractable. So we need to be able to evaluate this function. Not just to evaluate it, it has to have a closed form. And then in addition, we need to be able to compute this quantity where Q is now set to P. So let me go back up. Why did we do it? Well, I mean, yeah. So we need to plug for Q. We need to plug in this postivia and then be able to do this integral. Then we can do EN. But maybe we can't do that integral or we can't evaluate this posterior. That seems like a strong constraint and make and makes em a sort of a very specific tool that you can only apply when you can do a when you actually can compute the posterior and the associated integral for the elbow. Maybe sometimes we can get away with just saying somehow will allow Q to be something more like more restricted, something that like that some space of functions in which we can do tractable computations and then we just maximize this elbow this evidence lower bound within the space of things in which we can do this within the space of probability distributions Q in which we can compute and when we do that we'll then find some Q that minimizes KL divergence right because raising the evidence lower bound means we get closer to the evidence and as we do that we narrow the gap the gap is given by the KL divergence between our Q and the posterior over Z given X and then this Q that we find within the space of tractable distributions will be the best possible tractable approximation to the posterior that we might want to have and that's actually the entire idea. So this is an optimization process within some space of probability distributions Q and that idea is called variational inference. Well because because variational is the word for optimization or sort of the generalization of derivatives to function spaces. And of course that raises the question of what do I mean by a space of tractable distributions? I'll tell you in a second because it's actually surprisingly fluid this definition of of tractable. But before I do that, let me first make the connection to physics because I said this all comes that's all associated with people like fineman and so on. So what that you understand why people might be interested in that and actually also to make my case that machine learning is just a continuation of uh physics to everything else. So here's how this works in physics. We um uh um so these problems show up in what's called statistical mechanics or thermodynamics. So this these are this is the study of systems that involve many many many particles a large set of particles or large set of things that uh all have their degrees of freedom. Let's say their degrees of freedom are x and zed whatever this is like the position the momentum of all the particles in a gas. spin in some quantum gas and so on. These are called the micro state of a system. So it's the set of parameters that or variables that describe everything. But we assume that there's many of them like 10 to the 23, right? So a very very large number of them and therefore we can't actually see them. But we can write down how energy arises in this system. In principle it's just too many numbers to track. So we write down an energy function a total energy E that's just a function of X and Z and then we that energy as you now know already is associated with some probability distribution. So particles just take a configuration according to e to the minus energy. That's actually sort of a circular statement. The energy is basically this right? So the energy is the logarithm of the ne negative logarithm of the probability. So the typical challenge now in this theory like the statistical theory of mechanics or thermodynamics is that you can't actually go in and measure the state of every single particle or every single element of your system because it's large. There's like 10 to the 23 molecules and they all behave in a very complicated way. What you can measure though is what's called macro states. So you can measure the temperature, you can measure the pressure, um maybe you can measure some chemical potential, something like this. These are called macro states. And the way you usually measure those is by taking the system and letting it do something. You're letting it expand to measure its temperature or its pressure. You measure the force it exerts onto walls to measure pressure, something like this. And um so this means you're using some energy that's available in the system to measure the macro state in a very abstract fashion. And that's why the the this uh like the the the quantities that we tend to measure in terms of macro states are associated with something called the free energy. So free energy is really just a word that physicists and chemists actually use to describe properties of a of a microscopic system of a statistical system that one can use to do work. And uh that allows us to now separate the x's and the zs basically right. So so sort of as a in a sort of backward definition we could think of the um the variables x in our model as the things that we the aspects of the system that we get to measure in as a function of that that we get get to measure a function of when we measure macro states like temperature um or uh pressure or some work that the system can do. And the zetss are the other variables that somehow only indirectly affect the macroates. And so we can um define something called the free energy as the log negative logarithm of a marginal over our probability distribution under set. So it's the negative logarithm of the distribution of x. And um so x might actually be a very small set of of variables and zed might be a huge um remaining set of variables and there are actually now um well okay so I'll talk about free energy in a second again. So what what you now want might want to do is you might want to figure out what the remaining degrees of freedom are what the zs are that you may still be interested in for some reason but you can't measure them directly. So you're going to find some probability distribution that approximates the correct postivity on Z, the conditional you would get under all of the observations by fitting some parameters of a model. And that involves writing down something that um in the in uh variational theory of physics is called variational free energy. And that happens to be the exa exactly the negative of our elbow because you're minimizing it. So then what the physicist does is they try to find an expression that minimizes variational free energy which is the same as maximizing the elbow. So they're effectively doing approximate inference to find an approximation to the posterior on a on uh the distribution over some interesting latent quantities of a system that can't be directly measured through workable energy through free energy. And then a large part of the game in uh thermodynamics is that people have come up with various different ideas for free energy. And there's all these different words that if you've taken a a thermodynamics class or a statistical physics class, you've encountered like the helmholds free energy or the gibs free energy. Um, and they associated with these people. So Rudolfph Claus was a German uh physicist who coined the term entropy first. It came from him. He was looking for a word that sounded a bit like energy because he was thinking about the ability of a system, the differential of a heat in a system that can be used to return heat back and do work with it. He initially wanted to call it the favandlunks energy, but that was a weird word that nobody understood. So he was looking for a Greek word because he thought ancient languages survive for everyone and everyone uses the same ancient languages at least in the western world in you know 18 something. Um so he came up with entropy which is the Greek word for um transformation basically. So favandlong's energy um and as a contemporary said he that he succeeded in finding a word that meant the same for everyone namely nothing. And let's be honest to this day nobody really knows what entropy actually is beyond an expectation of a logarithm. Then he had was the first maybe to come up with this idea of free energy. He defined a certain type of free energy that involves entropy of a system and it's all about what kind of work the system can do and how you measure it. Then Josiah Gibbs had um also a very like influential English physicist in in thermodynamics had a a new way of expanding the idea of um free energy. So he added an extra term to it. That's pressure times volume. He calls this um u by the way is the potential energy plus pressure times volume. He calls this the enalpy. Somehow funny that they all came up with something similar to energy, right? Entropy, energy, enthalpy. And then I'll do a jump just to connect the dots to a contemporary colleague Dave Dave Bllye from Columbia University now in in New York. And to my knowledge, he was the first to came up with the word elbow for evidence lower bound. yet another word with an E. Um, and he's a statistician or machine learner we would say these days. Um, and uh um really like phrased this very clearly as the variational energy being the negative of a lower bound on the evidence in a general statistical model that we can apply way beyond physics and maximize. He also came up with some early very structured algorithms for language modeling long before large language models. So how does this optimization in the space of functions actually work? So it's going to work into the following. I now have several slides that try to explain in like as directly as I managed to because it's actually quite challenging how this idea of optimizing for a a functional for a function might actually work in practice. And you can maybe imagine that this is going to be hard. Why? Because well for functions over a vector space we know how to compute gradients and then the gradients are vectors in the vector space. So these vectors are tractable objects right? They're just a bunch of numbers that you can follow. You can write out an algorithm that just follows a gradient. But if you think about a functional that takes in a function. So a functional is an integral and we need to optimize this with respect to a function. Functions are infinite dimensional objects. So we can't expect that we're going to compute some infinitely long vector that we can somehow follow on a computer. Right? So the process is not going to be that we will compute a derivative in function space and then do a few steps to optimize the values of the function. That doesn't work. And instead we're going to find implicit statements that characterize a mode in the space of functions. And then those implicit statements will surprisingly sometimes actually tell us what the optimal function is. So here this is how is this how this is going to work. We start again with our model. We invent some latent variable Z XR are actually the given data and we have some parameters theta. We want to find a good approximation to the posterior over the latent variable Z given X. If you have a model then the posterior is well defined. It exists but it might be intractable. we might not be able to write it down. So what we're now going to do is we invent some family of distributions Q of Z. Um and then well we uh that's not a that's a stupid sentence. Why did I add this there? No, that's not actually what I want to write. This might be again one of the things that I fixed. Let me just do this. Well, yeah. Or maybe it's just confusing and it's what I wanted to write initially. But let's actually keep it simpler and just write it like this. So we invent some family of distributions Q of Z. And then we want to optimize the elbow which is this thing as a function of as a functional of Q of zed. I mean so Q of Z might well depend on our parameters. Sure. But the goal will be to find a Q a function of zed that maximizes this expression. And this expression is an expectation over this function under Q of zed. So it's an integral that gets rid of zed basically. So here it is again when we do that if we manage to do that somehow we don't yet know how we're going to minimize KL divergence to the true posterior. So what we find will be a good approximation in the sense of this divergence. So if we believe that KL divergence is a good measure of distance in the space of probability distributions then this process will find the best approximation to the posterior within the space of Q's that we consider. Now if you if you don't expect any if you don't impose any any restrictions on Q if you just say within the space of all probability distributions find the optimal Q we know what the extremum is going to be the optimal value it's going to be this probability distribution the posterivia but we know that we can't compute that or we assume that this is not tractable because if it were we would just compute it so we'll probably have to somehow restrict rict the family of distributions that we want to consider. And you might think that this restriction will be something like, oh, we're just going to assume that they're all gausian or something like that. And that's actually an option. But it turns out that sometimes you can be much more uh subtle. You can impose very weak restrictions on the space of functions that aren't even of a parametric form and you still get an actual tractable distribution out. And I'll show you how this works after the break. But first there was a question. Was it done? >> Ah yes. So we are not so this is this is the thing like let me let me just make sure you're not not getting confused. We're not yet going to assume that we will that Q of zed will be a parametric family of distributions. It doesn't have to be an exponential family even. It could in principle just be anything for the moment. So the restrictions we're going to put in will be something very subtle. It will not be oh it's a gausian with a known variance and we just have to optimize for the mean. That's actually also something we can do. This is this is uh sometimes called fixed form variational inference and I'm thinking that I might introduce it at the beginning of the next lecture also in another historical anecdote. But um yeah so well like for now let's actually keep it at that. So if there are no more questions then let's take our customary five minutes. There was a good question during the break which was also maybe triggered by my mistake of like the spelling mistake here or writing mistake here on the on the slide. Namely what about this theta? on Tuesday we spoke about optimizing the elbow for theta and that's the at the heart of the EM algorithm in in EM this maximization of the elbow in terms of the distribution Q of zed sort of happens almost implicitly right we just set Q of Z to the to correct posterior because we can track the full posterior and then we optimize in theta so in the setting that we're now discussing we can still optimize for theta that's still fine but this is not actually what I'm talking about right Now, so the what what we're now considering is this sort of maybe like to make this aspect of finding a good distribution over zed an explicit part of the thought process in in EM finding the best approximation Q is straightforward. It's just P of Z given theta sorry P of Z given X and theta the true posterior. But we're now saying well what if we can't track the full posterior? What is another thing we could find? So we're looking for a function of zed. That function of zed might still also be a function of theta. And then we can also optimize for theta. Fine. But it's just not important for now, right? So we'll keep those parameters like out for a moment and care about this thing as a function of zed. And we're looking for the best possible function of zed. And we know that the one best possible function of zed is the postivia of zed given x. That's the thing we'd like to compute because we want to do Beijian inference, right? But let's just say this is not something we can write down. It's just too complex. I'll come back to the example of Gaussian mixture models later on in a few minutes and then we'll see that there might be generalizations of that model where this Q of Z given X is intractable. And now I said this ominous vague thing is that we're somehow going to put some restrictions on the space of functions probability distributions Q of zed. And again it might seem like we would have to put some parametric descriptions. We'd like to say we are going to write this function in terms of a finite number of features. But actually just like in the gausian process world we got away sometimes without actually imposing a functional form on the space of functions. but instead using some kernel function. We're now going to get away without imposing a functional form but by imposing a factorization. So a factorization is a statement like the one in the first row or above that says so we're looking for this probability distribution over these latent variable set. So of these sets let's say there is like a 100 right or I don't know n of them and what we're going to have a weird clicker problem. Um, and what we're going to do is we're going to impose that Q of Z factorizes into a bunch of terms that depend on well, let's say some subsets of Z. So the Z I might still be multivaried subsets or they might just be individual scalar variables. It doesn't really matter actually. It's just some new variables, some subp part of this entire set of variables set I. And you know that factorization amounts to independence assumptions, right? So maybe by imposing independence under the our approximation to the posterior so our Q is an approximation to P of Z given X maybe the computation becomes tractable and we are not going to impose anything else we won't say what the functional form of Q is just that it has to have a product over the Z I so then what happens is the following so let's look at our variational bound our elbow that we're trying to maximize I'm going to put it in um in its definition. So that's the full definition. But now I've plugged in the fact that we assume that Q of Z factoriizes. So in our variational bound, remember that the variational bound is an expectation over Q over log of joint over Q, right? Log of P of X and Z divided by Q of zed. So if you take the logarithm, we can take this into a difference. And now we impose the fact that Q factorizes. So we have a product here and a sum here. And now we're going to do a thing namely that we're going to pick out one particular set of variables set J. So there's this product over all the Z I. We're going to pick one of those. Let's call that set J. I've actually highlighted this in this slide in red, but I realize now that you can barely see the red, right? So this set J is going to be for a moment the thing we focus on. Then of course what's going to happen on the computer is we're going to iterate over the J's, right? Going to do this for one J and then for all the other ones. So if we do that then we pick out this one term from the product. Um and then in here let's of move the rest of the terms of the product into the into a set of brackets. And then this joint integral over all the sets will be an outer integral over Z J. and then all the inner integrals over z the all the z i's minus um over here we'll get one particular term because here's a sum over the q i z i right so in one of those terms set j is in there so we'll get an integral over q of zj j log q of zj j that's the entropy of our distribution q of zj j or partial distribution or marginal over q of zj j plus something I'll call a constant. So why why is that a constant? Well, because it involves the sum over the logarithm of q of i. Um when you integrate out all the set i's that are in j and then there's an integral over q of zj j which is just a one. And so it's just a number that doesn't depend on zj j. Okay? Because we want to optimize for q of zj j. That number doesn't matter. So now we're going to look at this expression and we realize that it's sort of in a functional form. We have marginalized out the set I we've taken integrals over Z I all the I's that aren't J. So what we're left with here if you think of this as a function then it's a function of X and of Z J because all the other Z I not J have been absorbed. They've all been integrated out. So we're left with a function that we can think of as the logarithm of an unnormalized probability distribution over x and just z j. So it's a function of only set j but not the other variables minus uh sorry plus the entropy of uh q of zj j plus constants. So this is so if uh for a particular function zj this is just this is just a constant. So what we have left here is an x uh sorry uh is of course not right it's just an entropy but what we have here is um this expression which we can uh think of as an unnormalized probability distribution over zj which is equal to the expectation of the log joint under all the other approximations. So now if you think for a moment in terms of an algorithm what we're going to have in the end is some some for loop that in each loop goes through all of our Z eyes I and when it considers a particular set J this object here will be the expectation under all the other approximations under all the other Q of Q of Z I not J. So if we find an algorithm that sort of inductively always constructs approximations over all the Z I's then whenever we consider one of the one of them Z J all the other ones will be assigned a distribution and if we're lucky we might be able to compute this expected value. Now of course we don't know yet whether we're going to be able to do that. We just hope right and now let's just let's just go in and see if this might work. How would we then find a Q a Q of Z j which maximizes this expression? So here we now sort of instead of trying to comput a gradient which we can't in function space we remember that maximizing the elbow amounts to minimizing KL divergence. And here we now have a reduced elbow. We have an elbow over Q of Z J. And that means when we maximize this elbow, we're going to be minimizing the K divergence between Q of Z J and this unnormalized well now hopefully also normalized probability distribution implied posterior distribution over Z J which arises from integrating out all the other Q of Z I and that object right we know how to minimize the K divergence we just set it to zero by setting this equal to this. So if we have this log P of X and Z I this thing sorry Z J this thing then we can look at this expression and it might be of a functional form that we can identify as a function of Z of J for which we might know the normalization constant and then we're done we have a probability distribution over Z set J that minimizes the KL divergence maximizes the elbow and therefore is the best approximation within this loop that seems crazy at first let me show you an example and and you'll maybe get the structure. So we'll find an algorithm that iterates repeatedly in a for loop and inside one iteration of the for loop. It'll iterate over our variables set i 1 to n in each case setting the approximate distribution over z j to this implied posterior distribution which we get by computing this expected value of the logarithm of p of x under and zed under all the other set js. and then looking at the expression and finding hopefully that it happens to be the logarithm of a distribution that we can write down. Why might we be able to write down such a distribution? Well, because we're going to put in some assumptions about what P of X and Z actually looks like, right? We're still going to define our generative model for the data for X and Z in some nice structured form. For example, we might choose an exponential family. So if p this joint distribution of a latent variables and data happens to be an exponential family parameterized by some natural parameters then as you know we can write it like this that's the standard form of an exponential family with sufficient statistics natural parameters log partition function and then the thing that we need to compute this expected value under our approximations for the logarithm of the joint distribution log of p of x and zed will be this expression, right? So now we look where the zs are, right? So the zetss don't show up in here. That's just a log normalization constant. Capital zed and lowerase zed are different things. So we're left with just an expectation over the sufficient statistics times natural parameters minus log of zed. And now we know that we want to set this expression to zero with respect to the parameters. Um and that means that we're going to have to set the gradient of the log partition function to the expan um to the ex to this thing right so in terms of theta we're going to have to set it to this and we know from our lecture on exponential families that this expected value of the sufficient statistics sorry that this gradient of the lo partition function is equal to the expected value of the sufficient statistics under our uh p of x and um zed. And that might actually be something we can look up on Wikipedia if it's a known exponential family. And then we have an expression here that we can just plug in there and then we're done. So here's a summary and then I'll show you an example. So this so-called free form variational inference is a framework to construct approximating probability distributions that approximate a posterior P of Z given X that might not be analytically tractable. And we do that by iteratively finding um minimizers of the Kel D or overall finding a minimizer of the Kel divergence between the approximation and the true postivia which we do in an iterative sense by going through subsets of the variables and imposing some factorization. This by the way is called the mean field. I mention I actually had this on the slide but I sort of rushed over it. The reason why it's called the mean field approximation again comes from physics. That's also a little bit why I gave you this physics introduction. Namely, you can think of this factorization as saying there's this system that has contains lots and lots and lots of particles. There's too many of them for me to track all of them. And they all interact with each other all the time across time and space. And now what I'm going to do is I'm just going to pretend I'm sitting on one of these particles that J. And I want to know what happens to this thing that sits within the system and I want to want to see its own distribution. What's the best local distribution for Q of Z J that I can think of? Well, for that I need to consider the expected value of this overall probability log of P of X and Z. The overall systems energy. Remember that log of P of XZ is a negative energy under the effect of all the other particles. if everyone else in the system also pretends that they are alone. And then the way that they interact with each other is solely through this update step which you can think of as well as we just derived an expected value of the overall systems energy under everyone's individual contribution and that is a function of myself said J and that sort of mediation of everyone else on average to me that's the mean field right the average contribution ution from everyone else. That's why this is called mean field theory. So in practice this means we'll define a probability distribution. Then we impose a factorization. We do a lot of computations to hopefully find an expected value of the log of P of X and Z under all of these approximations. Look at it. Find that it has an algebraic form that we might be able to track. And if we can do this across all of the variables in our model, then we're done. And then we can just write an algorithm that iterates. So I'll show you now for the rest of the lecture how this works for our Gaussian mixture model in an extension of the Gaussian mixture model. But before I do that, I want to tell you a little bit also why I'm doing this in this lecture. Because every time I bring this up, I realize that what I'm going to do now next is going to be going to seem tedious and complicated and like mentally painful. And that's maybe actually true. But there is a very interesting very high level lesson hidden behind that. And what I want you to do now in the next 20 odd minutes is to maybe first keep the following story in mind that I'm telling you really like because I think it contains an interesting insight and it sounds like old men talking about war. So when when I did my PhD, I already mentioned this. What I'm going to show you now was really the daily work of machine learning engineers and statisticians. Namely, it was a very manual process. We would write down these graphical models that I'm going to actually I'm going to show it to you now. Graphical models that look like this. This would happen on a blackboard. You would stand together with your colleagues, write down models and lot think a lot about algebraic expressions. And there was a whole lot of theory about these graphical models. This lecture also used to be a lot about these graphical models. There were all sorts of concepts like induced factorizations and um um induced covariances explaining away structures, V structures, colliders and so on. We had words for all of these with the plates and so on because the work largely involved writing down math on a on a piece of paper or a blackboard. then doing the very tedious thing that we are now going to do together and then finding an algorithm arising from all of this that you then had to implement on a computer that was very hard to get right. There was usually a lot of bugs so they were they were quite hard to fix. Um and then you had these updates that we're going to see emerging now from the computations that would give an algorithm that actually runs that if you once you hide it all right worked really beautifully and then you could write a cool paper or build some product um and actually also build products. So this is how some early tools for like for machine learning were actually built like how uh online recommener systems were built for for large platforms like Xbox Live and so on. Um and um that felt like a really tedious process and it used to be a joke that that was like how how PhD students would spend their time and someone actually said well why don't we just always use mark of J Monte Carlo and a running g that by the time you had implemented this algorithm your mark of J Monte Carlo algorithm had converged and you could have just taken the samples from that so the value of course was in doing the derivations once and then having a very efficient algorithm and because it was so tedious very precisely I think actually because it was so tedious and you're going to see in a moment how tedious it was people stopped doing it by the time deep learning arrived and they got excited about gradient descent and said how about we just get rid of all of this nasty derivations and just do gradient descent from now on I'm going to highlight this on Tuesday next week and interestingly In doing that, the whole community actually moved to a different type of variational bounds which don't require these nasty derivations which are just gradient descent. And in doing that, we actually lost the plot for a little bit. It lost some functionality. And I'm going to show you on Tuesday an example, a high level example of how I think this happened. And actually none of us really realized I think for a while up until someone discovered that you need the structure again in some deep learning architectures and it can be very powerful and that idea was called attention. So I'll try and highlight that story and I'm telling this to you now giving this spoiler because I want you to like just keep your attention for 10 minutes as I'm doing these nasty derivations. I realize that they are tedious, but there is a lesson hidden in them if you squint through the math for why people might have stopped doing them. So on Tuesday, I showed you this um Gausian mixture model. Here's our data. Remember the scatter plot with the dots. There's like we think that there are clusters in there. We think that they are described by some Gausian cluster which have a mean and a variance. And um each cluster has a weight that is um parameterized by some probability distribution, a discrete distribution that uh is represented by a vector of numbers called pi k that contains numbers from between zero and one that sum to one. I showed you on Tuesday how we can optimize those parameters mu sigma pi by inventing these latent variables set n and k which then induce a factorization once we know what the set nks are then it's easy to optimize for mu and sigma and pi and then when we have mu and sigma and pi we can again compute a posterior over z and k that was called the em algorithm so we iterate between computing a postivia over zed optimizing the parameters given that postivia by optimizing the ex the the elbow the expected um complete data log likelihood and then setting again the posterior over uh z to the the resulting posterior if you fix mu sigma and pi we now realize that there was a special case of variational inference where we did completely free form variational inference basically we found the optimal approximation to the distribution over Z, namely the true posterior if we fix the parameters. Now what we're going to do is to do generalize even further to a point where the true posterior is not tractable anymore. So why might you want to do that? Well, one problem that this model has is that we have to say how many clusters there are. We have to say there are K clusters and what we don't know what K is, right? 2 3 4 5 6 7 10 I don't know. You have to look at a data set and that's difficult to do. I mean it's easy to do in two dimensions but difficult in 100 dimensions. So what you might want to do as a general concept um is a fully basian gausian mixture model. This was something that people like Neil Lawrence and Chris Bishop and David Mai were very interested in around the sort of early 2000s. Um and they wrote several papers about them which I forgot to site here. Maybe I'll do that next Tuesday. So we're going to do a general Beijian version of this model. What does it mean to be basian? Well, we're just going to replace these parameters pi, mu and sigma with actual variables. We'll say we'll be uncertain about what pi mu and sigma are and compute posteriors over them. And how would we do that? Well, first let's think about the algebraic aspect. So instead of a p of x and z given mu sigma and pi we're now going to have a joint distribution over p of x over x z pi mu and sigma by um the uh product rule we can write that as the likelihood we already know p of x and z given pi mu and sigma. So that's an object we already had in our gaussian mixture model times a prior over those three parameters P of pi, mu, and sigma. But we're already going to impose some generative assumptions. We're going to say this factorizes into a prior over pi, a prior over mu given sigma, and a prior over sigma. Why? Why? Why might might you want to do that? Because you've done a homework exercise about conjugate prior inference in Gausians and you know that there is a conjugate prior for Gausians with unknown mean and variance called the gamma inverse vishar and you did some of you did the homework on this with um Bernardo and we also know what the conjugate prior for a for a um categorical probability distribution is does someone know it's actually on the slide so category ical distributions are an exponential family. They have a conjugate prior and the conjugate prior is the it's the multivariate generalization of the thing I use to infer who is wearing glasses. It's on the slide so you can just point it out. It's the durishlay distribution. So the conjugate prior for a categorical distribution is called the durishlay and it's here. So I've just written it down. So we say and this is sort of what would use what used to happen in front of a blackboard, right? Someone said, "Oh, I've drawn this graph." So this is going to be a variable. Oh, we've got mu and sigma. They're going to a gausian. Oh, I know what the conjugate prior for a gausian with unknown mean and variance is. It's a gaus inverse var looks like this. We're going to put independent priors over every single component with a vishar distribution over the inverse of the coariance matrix the position matrix and a gausian prior over the mean and then someone knew from memory that you have to do this product form where you have a vishar prior over the precision and then a gausian prior given the precision that's how you get a conjugate prior and then pi is a probability distribution so everyone knows what the conjugate prior for a a categorical probability distribution is it's a dishly distribution which has this particular form that's a dishly distribution and those things have parameters right they are themselves exponential families because the conjugate prior to an exponential family are themselves exponential families they have par parameters they're called alpha for the dish a vector of numbers and w which is a symmetric positive definite matrix for the vishar and mu which is a count variable you know you might remember that conjugate prior have these count variables and um um m and beta for this Gaus inverse var and you just happen to know that if you lived in that time and now we realize okay that posterior over pi mu and sigma is not going to be tractable. Actually, we used to be able to read this off a graph. And if if you would have taken this class like four years ago, you would now be able to read it off the graph as well because you can see that there's a collider structure here coming in all these arrows pointing in. This is a variable. It's sorry, it's an observable. We've seen this. So therefore, all of these things become dependent on each other in the graph and we cannot expect to find a joint posterior over them. You don't know this anymore because I didn't do three lectures on graphical models and instead talked about deep learning and atlas approximations. But you know I can tell you that that's something you what you would be able to read off if you are an expert in graphical models. So now what people used to do is to say oh we're going to build a variational bound and we're going to do a free form variational bound. So the only thing we're going to do is we're going to impose some factorization over this joint distribution over these four things. And the factorization we probably want to have is that the variables sorry the parameters of what used to be our Gaussian mixture model will be one set of variables and these um cluster labels these set n case ks those will be the other variable why well because that's how em works right it was somehow very useful to have this distribution over set. Our approximation that we're going to impose is that whatever our whatever our approximation to the postivio over z pi mu and sigma will be it will have to factorized like this. There'll have to be one distribution for who belongs to which cluster and then one distribution for all the variables. That's it. That's the only only constraint we're going to put in. And that seems like a very weak thing to do, right? Just one little dot between the distributions. And now we do variational inference. And this variational inference with this really like cranking a big squeaking handle on a big cog. So you sit down with a with a piece of paper. Sometimes was like an A3 piece of paper in like landscape mode and then you start with the derivation. You say what was the rule to derive a variational bound? Well, we will need under this mean field assumption, this factorization, we will need to compute the the well to find this object, we will need to compute the expected value of the joint distribution the loin under all the other variables in our distribution. So if we want to find the the approximation to zed then we need to take an expected value of the log joint with respect to pi mu and sigma. And then afterwards hopefully we might be able to get an to get a distribution on pi mu and sigma by taking an expectation over the same thing under this distribution over zed that we've now just found. And then we can iterate. So actually it might be a good idea to start with zed. Why? Well, because the zetss are they are binary variables, right? So, we're going to clearly because they're just a set of binary variables, we will learn probabilities for them. I mean, the most general distribution over a set finite set of binary variables is a probability distribution, a discrete distribution that just contains numbers. Okay, so let's do those first. So our goal is going to be let's compute an expected value under whatever the other distributions are. We don't know yet what they are over this log joint. And now we look at the expression for the log joint. So we go back up again. We look at it. Okay. This is the thing we need to plug in, right? That's the object. And these have actually forms, right? So P of pi is the dish. P of mu and sigma is this thing down here which has this factorization structure inside. And this thing p of x and z given pi mu and sigma. Well, it's up there. this this this thing we're we're just going to plug that in with a logarithm in front. So now we find that this when we plug this in we'll get one term over p of zed given pi. Clearly this only contains pi. There is no mu and sigma. So our expectation over pi mu and sigma actually already simplifies. It's just going to be an expectation over pi. And then there's a term over p of x given z mu and uh sigma. So that contains no pi only mu and sigma. So that will simplify the computation as well. We'll only need to take an expected value over mu and sigma. And then if you go one slide back, there's a lot of terms here in the back that don't actually depend on zed at all. So we just can forget about them. Okay? So they won't contribute to our approximation to q of zed. So those are the two terms that we need to compute. And now we plug them in. So what was p of zed given pi? Well, it was actually happened to be a product over the individual Z N and then a product over the individual component K times pi K raised to the Z N power. So if you take the logarithm we get a sums over N and K Z N K expectation of log of pi K. So if you take the expectation under this joint distribution over pi, mu and sigma, well mu and sigma don't matter and pi is the only thing left and it goes through all of these sums and pass the z and k and we're just left with an expectation over log of pi k under whatever our approximating distribution on pi will be. We don't know yet what that is but we know now that we will need to be able to compute the expectation of log of pi k under this distribution. Okay, let's go to the second term. We just write this in. We leave it like this. Now we go to the second term. Here we need to take an expected value over the log of P of X given zed. So on the previous slide we had log of um we we can actually if we know zed right then we know what um log of P of X given zed is it's just a product over the individual XN over gausian and specifically the Gausians from cluster K indicated to us by Z NK. So we go through a sum over both the n and the k if you want to do if if we did the logarithm right all the products turn into sums over zn and k which is dv these 01 variables that tell us which cluster we are in that basically select an entry in the sum over k and then once we know which term in the sum we select we just have to write down the logarithm of a gausian so that's this right so gausians are exponentials of minus quadratic form with mean and coariance and then there's the log normalization constant that contains sigma. We need to take the expected value over this thing under whatever our approximation to mu and sigma is which we don't know yet but that's the thing we're going to need to compute. So now we know something we'll need to compute. And now we do this cool variational inference trick and we say what do we actually have here as a function of zed. So this expression here this will be up to normalization our the logarithm of our distribution over zed and that's a function of zed and we see that what this is is well it's a sum over n and k ah that means that there'll be a product in our probability distribution that factorizes over all the n over all the individual samples and then there'll be a z and k and some number that we have to compute. We don't know yet what that number is. It'll have to be computed from our approximations over q of distributions over pi and mu and sigma. But it's just going to be a number. So let's give a name to that number. Let's call it the logarithm of row n k. This is a real number. So it's it can be the logarithm of something. And then what we have here if we stare at this and like squint with our eyes and think of it as a distribution over zed we can take the exponential and we see that what we're going to get is a product over every individual datim all the individual components all the clusters times the number row n k raised to the z nk power that's not just a discrete distribution it's a factorizing discrete distribution. ution. Uhhuh. And now the only thing we need to make sure is that it's actually a distribution. So that means the uh row and k have to sum to one. And we can do that through the normalization. So by using this constant over here which we haven't used for anything yet to normalize and set um the actual parameters of our distribution to RNK the responsibilities from EM um which are now going to be equal to row NK divided by well the row sum of the row NK so that's a um dot sum Xis equal to one basically operation which will then give us our approximating distribution to Q of to the posterior over zed which is a fully factorizing probability distribution condition categorical distribution over the set n case parameterized by a responsibility rnk. Uhhuh. And so what's just happened is we now have we started by assuming that there will be some approximating distribution over the zs and one over the pi mu and sigma. And what the variational inference has automatically told us without us imposing it is that this bound actually factorizes over all the n. Every single x will have its own z n which contains k entries and they'll be independent of each other. We didn't impose this. It just happened because our model has this structure. This is called induced factorization. And now we know that this is a discrete distribution. Discrete distributions are very easy to work with. the expected value of a discrete distribution um sorry the expected value of the random variable under a discrete distribution is just the probability for um um like that's in the vector right so the the the the expected frequency with which you're going to see cluster K in for sample number n is literally just R and K. It's just a property of discrete distributions. Here it is. Yeah. So discrete distributions and now you could go to like the Wikipedia page for discrete distributions and redouble the properties of it like it's for example that it's uh it's independent between the Z KJ and it's like you know like everything is here okay no so but so far we only have an approximating distribution for the zed we need to to be able to build it we need an approximating distribution to pi mu and sigma and then we need to be able to evaluate this thing this expected value of this nasty expression here. So now let's turn our attention to pi, mu and sigma and see if we can find approximations for these. Now knowing that our approximating distribution to zed will be a fully factorizing categorical distribution. Okay, here I'm going to spare you some nasty computation and I'm just going to say uh I'm going to motivate some convenient notation. I'm going to introduce this NK which are the uh column-wise sums over the responsibilities. It's basically the expected number of data points in each cluster. We're going to compute a weighted um cluster mean weighted by the responsibilities and we're also going to wait uh compute a weighted uh empirical estimate for the co-variance in each cluster. We call that SK. That's just a convenient thing to have because otherwise the slide becomes so big that I can't use normal font size anymore. And now we're going to try and identify our approximating distribution over pi mu and sigma. Remember that we made no further assumptions. We just said that there will have to be a distribution over pi mu and sigma. And the way to find that is to find its logarithm actually is to compute the expected value of the log joint under q of z. That's our operator operational rule for variational inference. So let's go back again a few slides. Write down the log joint that here it is. I've just copied it over from the previous slide. We will need to take the expected value of this thing under our categorical distribution over Q over Z. And now we look at this expression and again we see okay this stuff does not depend on zed. There's no zed in here. So that's going to be easy. We just leave it in there. It's still I mean it's a function of pi and mu and k. So we'll need to keep it under consideration. We can't call it a constant, but it we don't need to take an expected value over Q of zed because zed doesn't show up here. Zed shows up here and there. Okay, so let's move the expected value inside. And then rearrange some terms because we also see that there is a pi here, but not here and not there, but there's a pi here. And there's a mu and sigma here, but not here and there. And it's a mu and sigma here and there. So let's rearrange the sum. Let's move this bit over here and this bit over there. Now um we are going to have a term that uh of a bunch of terms that depend on pi. One of them only depends on pi or nothing else. One of them contains zed. So we need to take an expected value over zed under our approximation which we now have. And um then there will be all the terms that depend on mu. One of them does not depend on mu. One of them depends h sorry one of them does not depend on zed. One of them depends on zed. So let's take an expected value over zed in here. But zed actually just shows up linearly in here in a big sum. So to take the expected value over zed n k we literally just replace this with our r n k our responsibility that we just found. Okay. So I've just copied over this expression from above. It's the exact same line again. And what we now see is if you stare at this, oh hang on, there is something over pi here. And then there is a sum over the k and a sum over the n. So this sum over k shows up both here and there. So we can take it outside. And in this q of zed we um have a factorization over the n's and the k's as well. So let's plug all those sums in. We do that here and we will get I've just rearranged terms. Um ah let's actually do one one after the other. Let's first do the one on on pi. So just this term. Let's just think about it. This will be a distribution over pi. Let's write this down. So it's the logarithm of a distribution over pi. So it's this term which just contains pi. That's our derivishly prior log derishly prior plus an expected value of under the discrete distribution over log of p of zed given pi. Remember what is what is p of z given pi? Well, it's just so the expected sorry the logarithm of p of z n k given pi is just pi k. So let's plug that in here. And now we have from the division prior this term from this term this we can rearrange the sums. You can put there's a sum over k everywhere. So we put it outside. There's one sum over n left. The log of pi shows up everywhere here and there. So we can take it outside of brackets. And our log probability distribution over pi will be a sum over k times a bunch of numbers. So that means as a distribution over pi we have to take the exponential. Our distribution will be a product over k pi k raised to some power up to normalization. And if you've paid attention in the exponential family class lecture, you know that that's a dish distribution. That's the exact form of a dishly distribution. So without imposing without requiring that our algorithm uses a dishly distribution, it turns out that our approximating distribution is a dishly just from the factorization assumption. The algebraic form emerges naturally. We'll have a dishly approximation which contains parameters. So numbers alpha plus nk where the nk are these things. It's just a name for them. And then um so you can look up on a and on on Wikipedia everything about dishlay distributions. And on Wikipedia you can click on this link. You'll find lots and lots of properties of dish distributions. In particular you'll find something about the expected value under a dish distribution of the log of pi. Ooh. And this is where the loop closes. So now we go back and we said well our problem was initially that to compute our expected value over uh so to compute our probability distribution over zed our approximation we would need an expected value over the logarithm of pi of k under whatever our approximation on pi was. So now we know our approximation on pi will be a dishly and we know that this expected value of the logarithm of pi under a dishly distribution happens to be something we can implement. It's that this nasty thing it's called the digama function and it's available in scypi. So now we're going to be able to write a piece of code that does the first step to compute the approximation to q of zed. And now we need to do actually the same thing again for mu and k. So um I can I'm not going to do this in detail here. Right. So so this is the stuff from before. I've copied it over again. Let's now focus on mu and sigma. Let's uh find the uh approximating distributions for mu and sigma. Well here we discover that there's a lot of sums in here again. So everything will become a product over independent distributions over mu k and sigma k. Um so there's again an induced factorization and if you do the math it's the basically the same as for the durish layer you find that there is a very tedious long expression that tells us what the functional form of the distribution over mu and sigma will be it'll factoriize over the clusters k and the individual distributions happen to be a product of gaus inverse gamas so they are again conjugate posteriors And um then you can go on Wikipedia and look up see if you can find. So we the last slide um whether we can find the terms we need to compute to close the loop to close the inductive loop and find the terms that we need to construct our approximation on um on the postivio on zed Q of zed. So for that from the previous slide we know that we need to be able to compute expected values of log of pi we know how to do that it's a division and then we also need expected values of the log inverse determinant and a quadratic form and for that back in the day people had big books for example uh there were a set of cool notes by someone called um Sam Roise that you could read up on and they would or you would go to your colleague down the corridor Tom Minka somewhere like that and they would tell you a yeah I know how to compute the expected values of of quadratic forms other gausian there are these nasty expressions here they are you can read them up and you just implement them and they involve functions that we know they're nasty expressions but they're are functions that we know so that means we can now do an an update loop that initially initializes a dish distribution for pi a gauss inverse gamma distribution for mu and sigma and the categorical distribution for zs at random and then it iteratively updates it first construct constructs a new approximation for the the categorical distribution over Z and K by computing these numbers which are proportional to this object which we have read off from the previous slides and these are all variables that we have because they are parameters of the approximating distributions for mu and sigma and pi for the dishlay and the gaus inverse gamma and then we do the update that was on the previous uh slide to compute conjugate posteriors given the Zs for pi, mu and sigma. And we keep doing that. If you want to see how this works, I actually uploaded the code to beat up Jupyter Hub, but it's so long that it doesn't make sense to look at it now. Instead, I can show you as a concluding picture how this algorithm works. So, here is our data set again in black dots. This this this um old faithful geyser data set. And now what I've done is I've initialized this um uh what's called a dirishly gausian mixture model or beijian gausian mixture model with nine clusters. I said like just try just try nine clusters see what happens. And this is the initial cut um um configuration. That's the expected load of each cluster. It's it's a draw from a sparse division distribution. And that those are the initial clusters in alpha value. to see their weight. So the the faint ones are very low weight and the thick ones are high weight and there are nine of them. They're randomly initialized. And now I just run the algorithm for a few steps. Every time it updates the Rs and then the parameters of the conjugate posterior to mu, pi, and sigma. And what this algorithm is now able to do is it discovers that they're just two two clusters. And back then people were super excited about this. You learn how many things there are in a data set. Isn't that a cool thing to be able to do? You show the algorithm the world and it decides for itself how many things are in your world. In this case, two and also what they look like. Namely, they have two means and two coariances. So this framework called variational inference is maybe it's a more civilized weapon from a more civilized age. It used to be the way that we built machine learning. But this is actually what we used to call machine learning building these very structured algorithms that required a lot of tedious work. In this case you first write down your generative model P. Then you decide you're going to do variational inference. You impose a factorization on your approximating distribution over the unknown latent variables. In the case of the Gaussian mixture model, there are z, pi, mu and sigma. Then you go through the math on a piece of paper like or on a slide and derive the variational updates. You find some forms that you know how to compute. It involves lots of known integrals, lots of oh, I know how to compute this under an exponential family. I know how to compute that expected value. And then you implement this as an iterative loop that keeps doing these updates. In fact, you can even evaluate the variational bound the elbow as you do and you can see the elbow rise. That used to be a buck fixing test. People did just compute the elbow and see whether it actually goes up. If it doesn't go up, you made a mistake. And um that used to be very tedious. It was very difficult to automate. It required a lot of knowledge and finger uh you know feeling um and elbow grease like the ability to sit down for several hours to do the derivations and it's what used to be a lot of work for uh research engineers and on Tuesday I will briefly tell you the thing that you probably expect that people stopped doing this because auto was cool and then we'll rediscover that in doing that we actually lost a bit of structure namely this induced factorization idea that really helps with algorithmic efficiency and it was rediscovered arguably that's going to be my point by some people in the deep learning community in the form of attention. So with that um please give feedback. I'll see you on Tuesday for the last content lecture of the of term and then afterwards there'll just be a summary. Thank you very much.