DATASCI 447 Lecture 4: What Probabilistic Methods Enable (e.g. solving the most important integral)

Administrative Stuff

A Little Review

Suppose we have a model where we are trying to learn a set of unknowns, \(\boldsymbol \theta\), given a model class, likelihood, and training data.

• What is the difference between the MLE and MAP estimator of the unknowns?

• When are these two estimators equivalent to one another?

• How does the MLE estimator fit with respect to the data? What about the MAP estimate with a non-degenerate prior?

A Little Review

We showed that Ridge regression as you have seen it previously is actually the MAP estimate assuming a prior structure on coefficients.

• What was that prior structure?

• How does that prior structure relate to the hyperparameter in Ridge regression?

• What did we mean by the idea that our prior wasn’t really on the coefficients but on functional form?

A Little Review

We can express Bayesian Linear Regression with Normal Priors as a learning problem:

\[ P(\boldsymbol \beta | \mathbf X , \mathbf y) = \frac{\mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I)}{P(\mathbf y | \mathbf X)} \]

• Which part of Bayes’ theorem did we not need to think about when we were just thinking about MAP estimates?

A Little Review

Since \(\boldsymbol \beta\) is a continuous parameter vector, the posterior is over a continuous quantity. Therefore, the more common written version of the denominator is:

\[ P(\boldsymbol \beta | \mathbf X , \mathbf y) = \frac{\mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I)}{\int \limits_{\boldsymbol \beta} \mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I) d \boldsymbol \beta} \]

• Looks complicated (and it kinda is)

• What is the marginal probability that I would see our observed vector \(\mathbf y\) given \(\mathbf X\) weighted by my prior over \(\boldsymbol \beta\)?

Where We Are

From the past two lectures:

• Likelihood = probabilistic loss

• Priors = probabilistic regularization

• MAP = posterior mode/penalized point estimate

We can sum all of this up using (part of) Bayes’ theorem:

\[ P(\theta | \mathbf X) \propto P(\mathbf X | \theta) P(\theta) \]

Where We Are

\[ \hat{\bbeta} = \min_\bbeta[-\log P(\mathbf y | \bbeta, \mathbf X) - \log P(\bbeta)] \]

If:

\[ P(\bbeta) = \prod \limits_{j = 1}^P \mathcal N(\beta | 0 , \tau^2) \]

Then:

\[ \hat{\bbeta} = \min_\bbeta[-\log P(\mathbf y | \bbeta, \mathbf X) - \lambda \sum \limits_{j = 1}^P \beta_j^2] \]

• Our L2 penalty is a normal prior on coefficients

• Regularization strength is the prior variance

Why Worry?

The probabilistic framework gives us something that hand wavy heuristics cannot:

• Assumptions are explicit, consequences are easy to track

• Explains why regularization works the way it does in a way other than “having a small coefficient budget is good and bias/variance tradeoff”

• Creates explicit methods for understanding the geometry of the hypothesis space. See the prior distribution over functions by drawing from the prior and the learned distributions over functions (functional uncertainty) by taking draws from a posterior.

Why Worry?

This is neat and all

But I really haven’t shown you anything that Bayes’ can do that you couldn’t already do

• Why spend time on Bayes’ in a deep/generative learning class if it just re-describes existing practice?

Why Worry?

The promise:

Bayesian Methods are going to give us a way to unify the procedures we use to optimize to the penalized likelihood given a model and the procedures we use to choose the best model across models.

• It creates a one stop shop

• And it’s going to solve the generalization problem with an analytical solution

Conditional Inference

When we do statistical learning, we have some data and we would like to learn a function of the data that will achieve some goal in a generalizable way

• For prediction: Have low prediction error for unseen feature/outcome pairs

• For generation: Produce high clarity and realism images of frogs that were not in the original training data

Conditional Inference

Conditional Inference

But, MAP/MLE with a L2 penalty gives us:

\[ \hat{\btheta} = \min_\btheta \log P(\mathbf y | \btheta) + \lambda \|\btheta\|^2 \]

Conditional Inference

In statistical learning, what methods do we use to choose values for hyperparameters and select between models?

Conditional Inference

The standard ML recipe:

1. Choose a set of models that you want to test (separate hyperparameters inclusive)

2. For each model, find \(\hat{\btheta}\) conditional on the observed training data

3. Generate an estimate of the true model risk

4. Compare across everything and select the one with the lowest true model risk

Conditional Inference

CV methods provide good estimates of the true risk for a model class!

• They are really good at creating a score for each potential model and ranking them correctly w.r.t. to the true risk

Conditional Inference

Conditional Inference

Model class:

\[ \hat{\bbeta}_{\lambda} = \| y - \mathbf X \bbeta \|^2 + \lambda \|\boldsymbol \beta\|^2 \]

Across models:

Choose \(\lambda\) to minimize generalization score

Conditional Inference

Define the LOOCV score for a model class as:

\[ \text{LOOCV}(\lambda) = \frac{1}{N} \sum \limits_{i = 1}^N - \log P(y_i | \mathbf x_i, \hat{\bbeta}_{-i}) - \log P (\bbeta_{-i}) \]

where

\[ \hat{\bbeta}_{-i} = \min_\bbeta P(\mathbf y_{-i} | \mathbf X_{-i}, \bbeta) P(\bbeta) \]

• The average penalized negative log likelihood for each observation finding the MAP estimate using the other \(N - 1\) observations

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

Conditional Inference

CV works really well in a lot of settings!

• Each fold removes one subset of data points (or no folds if single holdout)

• Slightly different training sets. slightly different parameter estimates

• Good models don’t have any folds where the prediction is unusually bad because we left it out! Minimize this behavior.

• Average predictive loss does this because averages are sensitive to extreme values

Conditional Inference

Cross validation methods create a valid score that can be used to assess the generalizability of models

• Statistical theory shows that this quantity converges to the true risk with large enough N

• Variance increasing in number of folds; Bias decreasing in number of folds

Decades of work show that this does what we want it to do!

When CV Fails

If we can use validation sets/CV to choose models from a pool of potential models, we should!

• But what if we can’t?

When CV Fails

First situation: the hyperparameter tuning space explodes

Have any of y’all worked with XGBoost before?

When CV Fails

When CV Fails

When CV Fails

When CV Fails

When CV Fails

If we have legit like 1 million parameters to tune and we want to choose the best ones w.r.t generalization error, CV is going to be untenable:

• Run the model with each combo in the grid

• Evaluate with K-fold or single heldout set

• Repeat 1 bajillion times

When CV Fails

Second situation: evaluating the predictive outcome \(N_{valid}\) times is prohibitively expensive

• Training Llama 405B parameter model required 16,000 H100 GPUs for months and cost 100 million + in compute at current Google cloud rates for a single final run.

• A heldout set would decrease model performance drastically while refitting over multiple hyperparameter sets would be impossible - by the time it finished (estimated by Gemini to be around 1.25 years for 5-fold CV over just one hyperparameter setting), the hardware would already be obsolete!

When CV Fails

The third:

In some cases, we don’t even know what we would iterate over to create a CV score!

• And what a CV score would even measure…

When CV Fails

Frog generative model:

What is the score? Goal: good pictures that don’t exist.

• What is the iteration axis? What is the evaluation axis?

When CV Fails

When CV Fails

Prompt Engineering:

• “Summarize this text” vs. “You are a helpful assistant, summarize this text”

Tune the prompt to give the best answer in the general case.

• Discrete hyperparameter space of the entire English language

The Probability Pivot

We need something else for modern problems!

The Probability Pivot

Let’s set up a slightly different CV scheme

Sequential Validation

Assume that our data points are iid/exchangeable (e.g. order doesn’t matter) and that we have placed a prior structure on our unknowns

Step 1: Predict \(y_1\) using the prior only

\[ P(y_1 | \mathbf x_1) = \int \limits_{\bbeta} P(y_1 | \mathbf x_1, \bbeta) P(\bbeta) d\bbeta \]

• We integrate over the prior to account for our prior structural contraints on the model

• If it helps, you can think of this process as taking the average likelihood that we observe over a discrete set of possible values of \(\boldsymbol \beta\) (our prior could be discrete and uniform!)

The Probability Pivot

Step 2: Predict \(y_2\) using \(y_1\) and the prior

• We’re pretending that the data arrived in a streaming fashion - use what we saw previously to predict what we see next!

\[ P(y_2 | y_1, \mathbf x_1, \mathbf x_2) = \int \limits_{\bbeta} P(y_2 | y_1, \mathbf x_1, \mathbf x_2, \bbeta) P(\bbeta) d\bbeta \]

The Probability Pivot

Steps 3 - N: Continue sequentially until all of the data are observed

\[ P(y_i | y_{<i}, \mathbf x_i, \mathbf x_{<i}) = \int \limits_{\bbeta} P(y_i | y_{<i}, \mathbf x_i, \mathbf x_{<i}, \bbeta) P(\bbeta) d\bbeta \]

• This method yields a sequence of predictive likelihoods.

• And these predictive likelihoods are going to look pretty familiar…

The Probability Pivot

Before we define exactly where we’re going to go with this, I’m going to show this process for Ridge regression.

• The integral over the coefficients is the hard part here

• Where we’re going to go is going to eliminate the need for this integral (sorta)

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

The Probability Pivot

Sequential cross validation is pretty much the same as LOOCV!

• The score (the average over all steps) measures how quickly the model learns - the faster we get to low loss, the lower the average

• How effectively it uses past information - is it over sensitive to new data (bad fit)? Under sensitive (bad fit)?

• How well it predicts genuinely unseen data.

It responds a little differently to the effect of sample size than LOOCV

• We aren’t always respecting training size \(N\)

• Different definition

The Probability Pivot

A model that generalizes well:

• Performs reasonably well early

• Improves and data accumulates

• Does not require seeing all data before behaving sensibly

• Does not oscillate between good predictions and bad predictions - it hits a loss and sorta stabilizes

The Magic of Bayes

Okay… That still seems like CV with extra steps.

The Magic of Bayes

Our sequential score:

\[ \sum_{i = 1}^N \log P(y_i | y_{-i}, \mathbf x_{-i}, \mathbf x_i) \]

The Magic of Bayes

Our sequential score:

\[ \sum_{i = 1}^N \log P(y_i | y_{-i}, \mathbf x_{-i}, \mathbf x_i) \]

is equivalent to:

\[ P(\mathbf y | \mathbf X) \]

Recognize this?

The Magic of Bayes

Bayes’ Theorem:

\[ P(\boldsymbol \beta | \mathbf X , \mathbf y) = \frac{\mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I)}{P(\mathbf y | \mathbf X)} \]

In its more useful form:

\[ P(\boldsymbol \beta | \mathbf X , \mathbf y) = \frac{\mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I)}{\int \limits_{\boldsymbol \beta} \mathcal N_N(\mathbf y | \mathbf X \boldsymbol \beta, \hat{\sigma}^2 \mathbf I) \mathcal N_P(\boldsymbol \beta | \mathbf 0 , \tau^2 \mathbf I) d \boldsymbol \beta} \]

The Magic of Bayes

The broader formulation

We see data \(\mathcal D\) and have some unknowns. Some will be optimized conditional on hyperparameters (e.g. MAP estimates), others will be hyperparameters.

Let \(\btheta\) be the set of parameters that will be estimated via MAP and let \(\mathbf M\) be the collection of model parameters (think ridge parameters, K in KNN, depth of a tree, tuning rates, etc.)

\[ P(\btheta | \mathcal D, \mathbf M) = \frac{P(\mathcal D | \btheta, \mathbf M)P(\btheta)}{P(\mathcal D|\mathbf M)} \]

The denominator of this Bayes’ theorem formula is equivalent to how well this model generalizes

The Magic of Bayes

In our current formulation, we’re thinking of likelihoods, so we want to make this denominator as big as possible

• Maximize the probabilility that we would see the data we did conditioned on our model and integrating over the things that we learned directly from our data

This can seem like a kinda tricky concept, but remember that \(P(\mathcal D)\) is a probability

• Over the distribution of all data that could have been generated by our model given our prior, what is the probability that we would’ve observed ours?

The Magic of Bayes

There’s a really subtle logic going on here related to the prior.

By integrating over the prior, we are saying:

My MAP estimate was this, but I have some uncertainty about the estimates as a function of my prior and lack of data

• If I perturbed the plug-in MAP parameter values just a little bit (w.r.t. the posterior probability distribution around them), how much would my predictions change?

The Magic of Bayes

• Models with tight posteriors (confident about parameters) → predictions don’t change much → more concentrated probability mass → higher marginal likelihood

• Models with diffuse posteriors (uncertain about parameters) → predictions spread out → probability mass diluted across many possible datasets → lower marginal likelihood

This integration marginalizes out our parameter uncertainty to evaluate: “How probable is this exact dataset across all plausible parameter values my prior allows?”

The Magic of Bayes

For Ridge, the denominator is of the form:

\[ P(\mathbf y | \mathbf X, \lambda) \]

To make this tractable, we’ll need to find this using the marginalized integral form:

\[ P(\mathbf y | \mathbf X, \lambda) = \int \limits_\bbeta P(\mathbf y | \mathbf X, \lambda, \bbeta) P(\bbeta) d\bbeta \]

The Magic of Bayes

We’ve actually already seen this! In our sequential CV section, we said that this had a particular form. When we roll everything up (e.g. add the log conditionals):

\[ \log P(\mathbf y | \mathbf X, \lambda) = -\frac{N}{2}\log(2\pi) - \frac{1}{2}\log|\mathbf K| - \frac{1}{2}\mathbf y^T \mathbf K^{-1} \mathbf y \]

where

\[ \mathbf K = \hat{\sigma}^2_{MLE} \mathbf I + \frac{1}{\lambda}\mathbf X \mathbf X^T \]

and \(| \mathbf K |\) is the determinant of \(\mathbf K\).

• This is a case where we can evaluate this integral exactly

The Magic of Bayes

The Magic of Bayes

If we can compute this quantity directly, then the learning problem turns into a two step optimization problem:

\[ \hat{\lambda} = \max_\lambda P(\mathbf y | \mathbf X, \lambda) \]

equivalently

\[ \hat{\lambda} = \min_\lambda - \log P(\mathbf y | \mathbf X, \lambda) \]

Then:

\[ \hat{\bbeta} = \min_\bbeta P(\mathbf y | \mathbf X, \hat{\lambda}, \bbeta) P(\bbeta) \]

The Magic of Bayes

Or treat it as a single optimization problem that can happen in one single go:

\[ \{\mathbf M, \btheta \} = \max_{\mathbf M} \int \limits_{\btheta} P(\mathcal D | \btheta , \mathbf M) P(\btheta) d \btheta \]

• It’s the golden goose of probabilistic machine learning!

• Fit the model parameters and the hyperparameters in one single swoop

• No CV, no inflated compute time

• Just raw generalization

The Magic of Bayes

The Magic of Bayes

Why is this the first you’re hearing of this?

• Any thoughts?

It Feels Like I Wasted Six Years of My Life Sometimes

This

Why This Integral is Hard: Dimensionality

Consider what we’re actually computing:

\[ P(\mathcal D | \mathbf M) = \int \limits_{\btheta} P(\mathcal D | \btheta, \mathbf M) P(\btheta) d\btheta \]

This integral lives in parameter space, not data space!

• Ridge regression with 9 polynomial features: 10-dimensional integral (including intercept)
• ResNet-50: ~25 million dimensional integral
• GPT-3: 175 billion dimensional integral
• Llama 3.1: 405 billion dimensional integral

Why This Integral is Hard: Dimensionality

The curse of dimensionality: volume of space grows exponentially with dimensions

• To accurately approximate an integral in \(d\) dimensions with grid-based methods requires \(N^d\) evaluations
• Even 100 points per dimension → \(100^{1000}\) evaluations for a modest neural network

When Can We Solve This Integral?

Tractable cases (where we can compute it analytically):

1. Conjugate priors with simple likelihoods
   • Normal-Normal (Ridge regression, Gaussian Processes)
   • Beta-Binomial, Gamma-Poisson
   • We saw this with Ridge: closed form solution using determinants
2. Very low dimensional problems
   • Single parameter problems
   • Some 2-3 dimensional problems with numerical integration

When Can We Solve This Integral?

Intractable cases (where we need approximations):

1. Deep neural networks - billions of parameters, non-linear interactions
2. Non-conjugate priors - when prior and likelihood don’t “play nice”
3. Complex likelihood functions - classification, structured outputs, GANs
4. The entire modern ML landscape - which is why the rest of this course exists!

Why We Need New Bayesian Machinery: VAEs

Variational Autoencoders (VAEs) are a perfect example of where classical Bayesian methods fail

The generative model: \[ p(\mathbf x, \mathbf z) = p(\mathbf x | \mathbf z) p(\mathbf z) \]

• \(\mathbf z\) is a latent variable (e.g., the “true” features of an image)
• \(\mathbf x\) is the observed data (e.g., the pixels)
• \(p(\mathbf z)\) is a simple prior (usually \(\mathcal N(\mathbf 0, \mathbf I)\))
• \(p(\mathbf x | \mathbf z)\) is a neural network decoder

Why We Need New Bayesian Machinery: VAEs

To train this model, we need the marginal likelihood: \[ p(\mathbf x) = \int p(\mathbf x | \mathbf z) p(\mathbf z) \, d\mathbf z \]

This integral is intractable! The decoder is a complex neural network, and \(\mathbf z\) is high-dimensional.

Variational Inference: The ELBO Solution

Key idea: If we can’t compute \(p(\mathbf z | \mathbf x)\), let’s approximate it with \(q_\phi(\mathbf z | \mathbf x)\)

The Evidence Lower Bound (ELBO) decomposition: \[ \log p(\mathbf x) = \mathcal L(\theta, \phi; \mathbf x) + D_{KL}(q_\phi(\mathbf z | \mathbf x) \, \| \, p(\mathbf z | \mathbf x)) \]

where the ELBO is: \[ \mathcal L(\theta, \phi; \mathbf x) = \mathbb{E}_{q_\phi(\mathbf z | \mathbf x)}[\log p(\mathbf x | \mathbf z)] - D_{KL}(q_\phi(\mathbf z | \mathbf x) \, \| \, p(\mathbf z)) \]

Variational Inference: The ELBO Solution

Why this works:

• \(D_{KL} \geq 0\) always, so \(\mathcal L(\theta, \phi; \mathbf x) \leq \log p(\mathbf x)\) (it’s a lower bound)
• Maximizing the ELBO is equivalent to maximizing \(\log p(\mathbf x)\) and minimizing \(D_{KL}\)
• No intractable integrals! Both terms can be estimated with Monte Carlo sampling

Why We Need New Bayesian Machinery: Diffusion Models

Diffusion models take a different approach to the intractability problem

The generative process works in reverse:

1. Forward process (fixed, no learning): Gradually add noise to data \[q(\mathbf x_t | \mathbf x_{t-1}) = \mathcal N(\mathbf x_t | \sqrt{1-\beta_t} \mathbf x_{t-1}, \beta_t \mathbf I)\]
   • Start with real image \(\mathbf x_0\)
   • After \(T\) steps, \(\mathbf x_T \approx \mathcal N(\mathbf 0, \mathbf I)\) (pure noise)

Why We Need New Bayesian Machinery: Diffusion Models

2. Reverse process (learned): Gradually denoise to generate data \[p_\theta(\mathbf x_{t-1} | \mathbf x_t) = \mathcal N(\mathbf x_{t-1} | \mu_\theta(\mathbf x_t, t), \Sigma_\theta(\mathbf x_t, t))\]
   • Learn to predict the noise added at each step
   • Generate new data by starting from noise and denoising

The same intractable integral appears - we can’t compute \(p(\mathbf x_0)\) directly!

Diffusion Models: Sequential Decomposition Solution

Key insight: Replace one hard integral with a sequence of tractable ones

Instead of: \[ p(\mathbf x_0) = \int p(\mathbf x_0 | \mathbf z) p(\mathbf z) \, d\mathbf z \quad \text{(intractable)} \]

We have: \[ p_\theta(\mathbf x_0) = \int p(\mathbf x_T) \prod_{t=1}^T p_\theta(\mathbf x_{t-1} | \mathbf x_t) \, d\mathbf x_{1:T} \]

Why This Comes Before Deep Learning

Before we dive into neural networks, we need to understand why they work through a probabilistic lens

Deep learning didn’t appear from nowhere - it emerged from decades of Bayesian thinking:

• Weight decay = Bayesian regularization
• Dropout = Bayesian model averaging
• Batch normalization = controlling posterior geometry
• Data augmentation = expanding the prior
• Early stopping = implicit regularization

Regularization Controls Posterior Geometry

Neural networks have vastly more parameters than data points:

• ResNet-50: ~25 million parameters
• Typical ImageNet: ~1 million images
• We’re solving an underdetermined system!

Bayesian Language Everywhere in Deep Learning

Dropout → Approximate Bayesian model averaging:

• Each forward pass samples a different sub-network
• Predictions average over ~\(2^H\) models (where \(H\) = # hidden units)
• Training approximates marginalizing over model uncertainty: \(\int P(\mathcal D | \mathbf M) P(\mathbf M) d\mathbf M\)

Bayesian Language Everywhere in Deep Learning

Convolutional Neural Networks (CNNs) → Structured priors:

• Translation invariance: “features should work the same everywhere”
• Local connectivity: “nearby pixels are more related than distant ones”
• Weight sharing: “the same features appear in different locations”
• These architectural choices = strong inductive biases = priors on function class

Bayesian Language Everywhere in Deep Learning

Attention & Latent Representations → Probabilistic inference:

• Attention weights as posterior distributions over relevance
• VAE latent codes as approximate posteriors \(q_\phi(\mathbf z | \mathbf x)\)
• Diffusion models as sequential posterior inference
• Modern architectures are doing approximate Bayesian inference!

Wrap-Up

I think this all needed to come first.

• And it gave us a chance to review a lot of the stuff that you talked about in your earlier classes

Next class:

• Modern predictive methods, Curse of Dimensionality, and Why KNN and Trees fail