Prof Chris Holmes | Bayesian fitting and evaluation of complex models arising in...

[00:00:07]yeah thank you to the organizers for inviting me so I kind of have two hats I have one at the University of Oxford I have a chair and biostatistics that's a joint position between the medical school and the computational statistics and machine learning group in the Department of Statistics and I'm also scientific director for health at the Alan Turing Institute which is just down the road here and I would say and it's a genuine call out if you want to come and have a look and have a look around you know it's a really interested in space it's the National Institute for data science and machine learning we have a very we're in the British Library it's a really nice office space we have about 300 researchers on space what kind of 300 research all working on the kind of methodology and theoretical side of machine learning and data science and my interest here and then reflected in this talk I'm really into this in interesting this interface of kind of statistical machine learning and applications for health hence why I'm here and so I in my talk I said I was going to talk about Bayesian fitting and model evaluation and I'm actually not going to have the time to talk about the evaluation side although I think that's that's interesting in its own right but so I mean I mainly focus on fitting Bayesian models in the context of a robustness argument that says well you don't believe you actually have the true model yeah so and of course we never have had the true model so I'll talk about some recent results that we've been working on some of it's a bit quite theoretical but hopefully give you some insight into the methodological implications of that so I'm very interested in working in Bayesian models in high dimensional spaces involving different data modalities which I think we're kind of challenging from a foundational perspective for Bayesian inference because Bayesian inference assumes you somehow have captured the true sampling distribution I'm going to talk about uncertainty quantification for parameter values and this kind of wouldn't same kind of new approach but this approach where we use in BAE nonparametric models to Train parametric models so if you think about like you take a generalized linear model logistic regression so rather than fit in that using a conventional Bayesian approach yeah you're using wind bugs or something what I'm going to do is I'm going to take first of all kind of train a nonparametric model to then train the parametric model and I'm going to try and convince you that that's not a ridiculous thing to do and I'll have some examples a small example at the end from a taken out of a genome-wide Association study actually focus region of a genome-wide Association say assume a little bit of kind of Bayes 101 so Bayesian statistics is kind of founded in decision theory an optimal decision-making and uncertainty principle following the work of Savage and Bayesian statistics essentially uses probability as a vehicle I think it's a really important thing that sometimes there's a confusion about you know what's random or or what's not to a Bayesian that doesn't matter and it doesn't have to be a notion of a kind of true stochasticity to something for use to use probability so we're using probability just as a vehicle yeah we're just using it as a vehicle to express degrees of belief in something that we don't know about and then we're going to use probability calculus Bayes rule to update that and the nice thing about that is shown by savage is that if you take this approach of using probability as a vehicle to express uncertainty an updated using Bayes rule then you obtain all sorts of coherency that you might not do if you just tried to do it in a different manner and it also based on statistics allows us to apply probabilistic modeling to a wider class of problems than would be done using non-bayesian approaches you know for a you know a teaching example is the number of coins in my pocket now you know doesn't you can't that that's not really a random variable in any sense but from a Bayesian perspective it doesn't matter yeah all that matters is that you have a state of uncertainty in the number of coins in my pocket and you will use probability to express that uncertainty use probability as a vehicle and then if I gave you some information you would then use probability calculus and Bayes rule to update your beliefs so I think that's the kind of an important distinction but there's not a kind of notion of randomness there we're just using probability as a tool for information propagation okay and at the heart of course of Bayesian inference the way it gets its name although and you know you don't actually have to do this but at the heart of Bayesian inference is Bayes rule or Bayes updating that says that the posterior distribution beliefs your beliefs about some unknown theta that you care about having gathered some information are proportional to the prior believes times the likelihood and so that's that kind of prescriptive updating rule its prescriptive it says that if you have prior beliefs and you have some data and a likelihood function you must update in this format otherwise you could be led to states of in coherence which would mean that if used this side of the room has the same prior beliefs as this side of the room R on some unknown quantity of interest and I started feeding new data one at a time you know got n observations I feed you one at a time in your updating and this side I give you all of the data in one block and you update in one block now for coherency you would want for you to agree it would you know you've got the same information you had the same prior and if you wish to remain coherent in any decision-making that you then add rep that you would then address you have to use this that's the prescriptive nature the axiomatic nature of Bayesian updating as in in Savage okay so in Bayesian statistics that likelihood function actually has another special it's called the likelihood but it's a it's a bit more specific it's actually a sampling distribution which means it has to be normalized yeah so it's actually so the likelihood function doesn't necessarily have to be normalized but the sampling distribution in the Bayesian model does and it represents the joint probability of generating the data that was observed and it's written you know as a function of the thing you care about the theta so for example you know if observations are discrete then the likelihood assigned to a parameter theta is just the probability of generating the observed data from the model given that value of theta ok so formally that like I said what's the probability if I was to generate data over and over again from this likelihood function with this value of theta that it would exactly match the observation and just as an aside that it's this constructive argument that underlies things like approximate Bayesian computation which is a kind of a new approach to working in in challenging model spaces however and this is kind of part of my interest in a bit more towards the kind of foundations is of course well not of course first of all this the statement that really resonated with me since I first studied Bayesian statistics is all of Bayesian statistics is model-based yeah and that's not true of non basic everything that you wish to do in Bayesian statistics involves building models yeah so if you were testing yeah so non basics is if you were testing or doing variable selection or do an inference you might adopt different approaches whereas in a Bayesian framework you don't you build joint probability models yeah and and of course the other the flip side of that though is we know all models are wrong and so if all models are wrong and all the basing statistics is model-based where does that leave us in trying to interpret this object so kind of formally this object is the posterior distribution you've run a Bayesian analysis you've run it in wind bugs it's kicked you out a characterization of the posterior distribution what actually you know what what is that object that's coming out for my analysis how should you interpret probability statements and I've had people included referees actually on some early papers say we'll just you know why are you doing this just keep calm and kind of carry on other things so so but if you do you know if the models wrong and I take decisions using that posterior you know well first of all what is it represent and secondly can I trust this object to plug into a decision analysis yeah now of course you know careful model check-in model validation kind of alleviates a lot of those concerns so you wouldn't say well we all know we don't ever really capture the true likelihood function so we would want to do careful posterior analysis checking residuals looking for the perhaps structures transformations of variables to look at that but I was cut you know it's still kind of any interest in to to address this issue of the modernist specification inherent modeled Mis specification so I guess to emphasize this that you know mathematically we write it at this is that you know if you can entertain that nature here is f naught so nature is give it is giving you samples yeah so nature's providing new samples out of F naught which is the true sampling distribution and but you're using a sample additional likelihood function f theta her model yeah and this just says that for no value of theta does that match nature yeah so formally Bayesian statistics assumes that MIT nature's true data generating mechanism is contained under the support of the prior that's why the prior integrates the one you know you put of course zero measure outside of the support of your prior ie with probability one nature is in is in the model class okay so but increasingly that's kind of being challenged especially in in the world of data science and modern applications you know how can I start to think about a true generative model probabilistic description especially when we're trying to combine different data modalities so you know people working on UK biobank might be you know taking half a million genomes or genotype data on half a million individuals trying to combine that with a hundred thousand of those who got imaging you know it's increasingly challenging to think that you could actually either do that in a you know we you know in an accurate manner but you know as I said models are just simply models of course and they're used to kind of distilled salient features that's why we do model in is to so we can focus in on the salient features of the problem to allow us to kind of make good decisions but you know formally Bayesian statements have predicated on the model space character okay so where does that leave us so what happens when you fit a false model to data so consider one notion of a best fitted model so here I've got a likelihood function f theta and I'm just going to maximize it so you could translate this this equation here into I've called GL m let's let's take a generalized linear model take a logistic likelihood function and let's say I've just called G 11 on some data and it's given me the estimate of the parameter the maximum likelihood estimate so the odd max just means find the value of theta that maximizes that score yeah ie the MLE okay and and in that sense you can just think of this object here the log likelihood function is just a score it's just a kind of scoring mechanism or a scoring rule and then it's perfectly valid to consider the properties of you know running a GLM model on data outside of the fact the data might not really come from a price on distribution or really come from banila with logistic probability yeah and and once we do that we can then start to think about well what happens you know you know what happens to the output or what's happening to the output model you know as you get more more data so let's have a quick kind of divergence into asymptotics well it does actually there's no real asymptotic yeah but it says what would happen to think about it I've got my GLM model you know I send it in my hundred data points suppose you had two hundred three hundred four hundred a billion a trillion you know the output is going to move somewhere where is it moving - yeah so what's the output of a say GLM trained with an infinite sample size and it turns out to be to look like this yeah so what it says is that the output from your GLM will be the the thing that maximizes his expected log likelihood yeah so f naught remember I said that's nature yeah so that's not unknown and never known yeah for data arriving after out of F naught of X and moreover your Bayesian models gonna head there - doesn't matter at this point if I give you an infinite sample size the prior unless it's you know if it has a reasonable support you know over the data is that's going to get washed out yeah so the prize gonna get washed out and your Bayesian posterior is just gonna head to a spike yeah and it's going to head to a spike at this point so they're both the Bayesian and non Bayesian models are heading towards theta naught and perhaps an easier way of thinking about that is is in the following case that I've got X I these the observations that nature's given me and it's going to give me an infinite number of observations and I'm just going to call GLM you know with an infinite number of authorizations and and this is actually kind of an interesting point spec point in the model space so you can see if you look at it for a little bit of time that that value that the place you're heading to and she minimizes the cool back lie below divergence between your between your model F theta and nature's F naught so you're not heading you know arbitrarily and you head into a well-defined point in model space yeah and that's that's irrespective of whether the modern is true or not yeah if the model is true the cool back liability I've urgent's to be zero yeah because your head to the true sampling distribution but if it's not your head to the point in the model space that's closest to nature's f naught that's outside the model space and that's been well known since White's paper and the robustness of maximum likelihood estimation okay and so so we've kind of solved one thinking about it this way we've actually solved one of the the initial problems which is now if you're going to put a prior so if you suppose your goal good basin's you're sitting down and eliciting subjective priors whereas with in conventional basing update you'd be kind of forced to think about the true value of theta now here we're just going to conceptualize thinking about where this value is yes so now your prior is about you know where would this GLM model end up given an infinite data and that's what I'm going to mean is captured in pthey tonight yes it's going to be the prior on where is this value of this kind of in some sense optimal value theta naught okay so it's where's the posterior going to converge to as I give it more more data okay okay of course we don't have an infinite data yeah so we have finite day to end data points and so what I would like to know yeah is given the output of my GLM or a finite data point can I update my beliefs it feels like there should be a belief update about theta naught yeah so given what I've currently got as an output for the maximum likelihood estimate can I say anything now about where theta naught may or may not be remember the prior was specified about that position we're in the model space we'll take the closest in KL so you know is there a posterior update now you could just simply do that using Bayesian updating you could say well I've got my prior I'm just going to put it through the likelihood function but as I said before that that makes a kind of an assumption that the data is really coming from the likelihood for some value of theta and and what we've been looking at is how can you provide an update a Bayesian update in the absence of that so thinking that actually in the absence of the sampling distribution being true and the kind of the high-level view is that you know if I think about this theta naught is the sense of the output under an infinite sample size the posterior under certainty in theta flows directly from uncertainty in nature's F naught yeah if you knew the sampling if you know if I knew what ethanol was nature's F naught I could just I could just generate an infinite sample fit my GL m and I'd be done yes so the uncertainty in the model comes from an uncertainty you know in a wider sense from from f naught yeah from nature's ethanol and so if you think about it that way if you're committed to using this model and the uncertainty in this sentence in the optimal value flows from uncertainty in the sampling distribution why don't we just attack that directly so so if F naught is unknown well I could be Bayesian about it yeah because I say well if a new F naught I'd know the value of theta I don't know F naught so let's just tackle that directly and that's gonna that enters into this area of called Bayesian nonparametric s-- so what we're going to do is we're going to try and put a prior to correctly on earth yeah and this hints at the thing I want to write it we gonna talk I said I'm gonna train a nonparametric model yeah to train a parametric model this is what I mean I'm gonna I'm gonna use a nonparametric update for F and then use that to train my my parametric model and see what happens so we place a prior on the space of distribution functions and then learn about theater that way and is the essence of basic nonparametric learning okay I think I've said this so if I could generate suppose wind bugs could give me an output on a distribution function then I could then train my model against this because if you give me a sample now of the distribution function there's no uncertainty in theta so if you imagine if you were using something like Gibbs sampling where you say draw F draw a sample for the-- to draw after us there's no the posterior for theta given F is a Dirac that's there's no uncertainty so that's just going to be a maximization step rather than the sampling because I said if I draw fi from the posterior on the space of distribution functions that are know precisely how to sum well there's no sampling variation in theater I okay and repeating this operation would then give me a bag of Monte Carlo samples of the theatres that characterize the marginal posterior distribution that comes that flows from the uncertainty in F and it it turns out that there's a really interesting and useful class of nonparametric models in the bayesian literature called dirichlet processes which are very amenable to this problem and for technical reasons and I could see no at least one or two people in the audience who are going to know why we're going to take something called a degenerate directly process and the reason why we're going to take a degenerate directly process is it's going to but all of its support on the observed data so think about the empirical think about an empirical likelihood function I puts the support on the data we're going to get something a kind of a little bit similar to that so we're going to take a degenerate datura say process which is going to put all of its support on the observed data and the reason for that is going to give us a really simple update a surprisingly simple update for this problem that we've just characterized okay because if I draw this F from something called the Bayesian bootstrap which goes back to Rubin yeah again it's this is very similar to kind of bootstrapping if we're non-bayesian yeah if I draw a death from this duration a process what happens is that this object that minimizes this KL under an infinite sample size of this F just turns out to be a weighted log likelihood yeah weighted maximum log likelihood from so the MLE would have weights one yeah this sample has a weights and moreover those weights have a very particular form there there's something called der is a 1 1 1 so the Duras lay is 1 1 1 ones is a uniform distribution on the simplex what does that mean it means that there they have to sum to 1 the weights sum to 1 and then they're uniformly distributed so any point in the space of positive weights that sum to 1 has equal measure or any region so it's a kind of a flat prior over the simplex over the uniform simplex but if you look at this this is now looking really nice because what this says here in order to do that kind of bit of a fancy math thing as I was talking about for all you need to do is draw these stochastic weights and maximize a weighted log likelihood and moreover this is a Monte Carlo drawer it's not a Markov chain Monte Carlo drawer there's no serial dependence I can draw these F's independently in parallel yeah so so what does this say it says if I want to fit do this kind of strange thing of fitting my Bayesian generalized linear model using nonparametric s-- I'd say to you well how many samples do you want you'd say I don't know 10,000 I'd say right just draw me 10,000 stochastic weight sets in parallel and push those into your GLM function maximize this weighted log likelihood not the yeah and then kick those out and you're done yeah and what I need to persuade you is is that the things that this kicks out are useful yeah it's not just kind of yeah it has to have this this particular form yeah so this is just a simple simply a weighted maximum likelihood estimate with random weights and if we repeatedly redraw the weights we end up with a bag of samples from in this case okay this isn't actually this procedure again for people in the audience someone already knows is was originally specified in a paper or RSS be read paper by Newton and Raftery in 1984 so everything I've said to you so far has already been done yeah so you know you might think why not you know that does that mean it's no good you know it was a red paper why doesn't people pick it up and I think there were two reasons first of all this thing on the left hand side is actually not necessarily a good approximation to your original Bayesian model and I can characterize that yeah but that original Bayesian model assumed the model was true so when people looked at this paper they said well you know I've got my true Bayesian model if I if I look at this object here does it match and they said no it doesn't you know doesn't really get the same kind of coverage so we don't like it yeah it wasn't picked up on and the other thing you have to remember that it was just as MCMC was kicking off so this big push into Markov chain Monte Carlo yeah and sample sizes were relatively small too the thing that we're dealing with now so there wasn't so much of an issue in terms of runtime and convergence and burning for the MCMC yeah but it turns out and this is I guess part of the work that we've done is that if if you're anyway outside of the modern space if the model isn't exactly true yeah if the model is precisely true do Bayesian inference do standard basing updated if there's any sort of model missed specification this object on the left hand side is much better than the original bit much better than so I'm getting over is has certain properties which the original Bayesian approach wouldn't have so it turns out that those bag of samples from maximizing this weighted log likely had our particular interest in objects and I'm not going to bore you with with the asymptotics here but essentially some of you might know I've heard about this kind of sandwich covariance matrix estimation is that because we're using a bootstrap the posterior coverage of the bait of this bayesian model using this maximization obtains the right frequentist kind of risk whereas with Bayesian models you've got no guarantees and perhaps the simplest example to this is to think about a Poisson regression so suppose you're fit in or even just fitting the Poisson distribution suppose you're fitting a price on distribution as a Bayesian and the model was actually Kate was over disbursed yeah so they're actually the data has got less information than you think it has and what's going to happen is your posterior interval is going to be too tight yeah your modern term is specified and you think you've got more information than you actually have and vice-versa if it's from an under Sun to disburse price on you've got more information in the data than you actually think so there you're credible intervals are going to be too wide yeah now you don't I mean again careful model checking in calibration yeah but this bootstrapping approach kind of automatically corrects that so you get the kind of correct in some sense a correct learning rate outside of the the modeled specification and yeah I guess the other interesting thing is that predictively it dominates in a way if you define predictive risk in terms of like a callback lie blur an information measure of the risk prediction when you make it then by doing this bootstrapping you get you've you get better lower risk in in your predictions okay but so far is that's the other reason it wasn't like people said there's nothing bayesian about it yeah so it's actually it's really interested even if you want to know something about the you know the people know the RSS be read papers and notoriously they invite you to give a paper and then people attack it and then that's the kind of tradition but it's worth reading and of course you know people looked at this thing and they said well there's no prior yeah so how can this be a Bayesian model so that and that's in a sense is where whites are really that in its form is a bad kind of it's not necessarily a good approximation to a Bayesian model and so we would like to incorporate if we're Beijing's we'd like to incorporate information and the simple approach that we're going to use is use synthetic data and we use synthetic data in in the following math manner if you have a prior on theta because your base E and your likelihood is a is a sampling distribution it's a full sampling distribution I can just play the following game I'm gonna draw a theta from your prior yeah I'm going to draw a synthetic data set of size T from the sampling distribution from the likelihood function and and then get a I'm going to then concatenate it I'm going to add my synthetic data are going to now have an additional training set which I'm just going to combine with my observed yeah and so then I've got a data set generated under the prior and my actual observation and again I'm going to draw stick a stick weights and I'm going to do this maximization trick and we're going to have a wait in here a see in front of the synthetic data which is akin to a prior sample size so if you set C to be n then it's like saying you know I've got just as much as information in my prior as I have in my real data you know if you wanted to do kind of unit kind of a unit information Pryor would put this as one yeah and then you just randomize the weights the weights are drawn to rich lay and then you maximize independently in parallel so this is what the kind of algorithm would look like I've said all of this so I don't need to go over it again do these B's drawers but again there's B drawers can be done in parallel draw synthetic data add it on to my real data draw a set of stochastic weights maximize and then return the samples and and here's it showing what this does on a on a simple linear model so and thanks to my PhD student for doing this no idea thank you I would say an Edwin go and just do this you know a movie of this it should be ok so here I'm just fitting a simple linear model and and what you so what what am i showing here's the data up in there so the points of data points you know the MLE is this is the line that doesn't move yeah I'm fitting this model on the top right and then here are the stochastic weights yeah that are kind of moving up in it if you didn't bootstrap in yeah the kind of conventional bootstrap these weights would be like steps of 1 over N yeah and about a third of the points would have weights 0 yeah but here we've got it it's like a continuous bootstrap so here we've got the random weights and you can see if I of course we know linear regression the leverages and or in the points at the end so the best thing is to wait and if you see a big weight there you should see it starting to track down on to this point and equally and and this is the posterior output from a pretty small run but the thing to show is this point yeah yeah because it's an independent maximization it just jumps around yeah and the other thing this the simulation isn't really capturing is that this looks like MCMC but it's not we just do it in parallel so if you give me you know P cause I can deliver P samples in one iteration so it's it's trivially paralyse Abul under that and you get the usual things you get posterior credible intervals joint distributions on on your parameters okay so i've said that it's true paralyze wall it's also because you're doing these independent maximizations from random starting point say it's good at capturing multimode multimodal posteriors and so again this is just a kind of toy example with a something called a mixture model Gaussian mixture model and that's known to be a kind of a bit of a beast for for MCMC because and it's often used because we've got known symmetries we've got label switching in the parameter space I can I can pick up these parameters I can switch them over you've got exactly the same model yeah so the likelihoods unidentifiable yeah and that means that if you're MC MC is running properly you should see all of these posterior modes but because we're doing these random restarts you just capture them so here's the output from this kind of weighted likelihood bootstrap what you're seeing is you're looking at the uncertainty in centering of these scarcity of mixtures and it just picks up all of the modes and then this is an output from the no u-turn sampler so these are kind of Gibbs samplers that are run and and incense as expected because you're doing a serial update conditional where you currently are it finds it much more kind of challenging to kind of move around the space once how much time they might at a time five minutes okay so here's the genetic Association example so we wanted to do sparse logistic regression so there's the likelihood function and the prior that we're going to use and which goes into the loss it looks at it's a student it's a I actually generalized student but it's something akin to a student distribution so here's the penalty on the coefficients in the logistic regression model and you've got a scaling parameter C that kind of says how much do you want to squeeze those parameters towards zero and how much do you want to give them kind of free reign and the interesting thing in in genetics is you get really high correlation structure that comes around because of recombination so recombination breaks down genomes very slowly over generations about one recombination per chromosome per generation so you are mosaics of your grand paternal and grand maternal genomes but big blocks so your chromosome one will be split in half knew your mother's side and from your father's so what that means is over there times you still see in like European well in most all populations these recombination blocks so then you can think about this is the kind of variance covariance matrix so it's really blocking and that normally means that you have trouble but but because again we can do these random restarts this is looks like a lasso type plot but it's these are kind of credible intervals and what they show is that if you start relaxing the penalty you get more more kind of variables into the model and what I've shown here is the the lines of the medians yeah of the posterior and when you first look at this you go hold on a second something's a bit weird here because you're getting really kind of lot of jiggle in the medium but if you if you zoom in on that you can see exactly why in the next plot so this is now our kind of this is now looking down on that parameter as a function of the amount of shrinkage and it's because it's bimodal yeah so I've kind of said well it's true but shouldn't be in the model or it should and again you kind of you capture that so this is just from the posterior output so this is my final slide traditional Bayesian approaches formally assume that the sampling distribution is contained within the truth this nonparametric approach doesn't it's important to note that the nonparametric approach wouldn't give you the same as the parametric kind of conventional Bayesian they're targeting the same parameters they're both learning about exactly the same parameter so whether you like it or not the Bayesian we're both heading in the same direction as you gather more data but they're conditioning on different states of knowledge so this right hand side is conditioning on the day the true data generating mechanism being within the model space and just to note that it's scalable because for iid likelihoods you can just paralyze it and thank you and if you wanted to learn more we've got a paper in and kind of statue journals but then also one in like the nips and just recently in ICML thank you [Applause] you