Problem 660428 · medium · Level 06 Heuristics & Optimization

A Prior Is a Penalty

MAP estimate · normal prior · conjugate prior · regularisation · ridge regression

A retail chain changes its advertising spend in a few stores and records, for each store, the change in spend x (thousands of pounds) and the change in weekly sales y (thousands of pounds). The model is y = w * x + noise, with independent normal noise of known standard deviation sigma, so the log-likelihood of a slope w is -sum((y - w * x)**2) / (2 * sigma**2) plus a constant.

Past campaigns suggest that effects are usually small: the analyst's prior belief is that w is normal with mean prior_mean and standard deviation prior_sd.

Write shrunk_slope(xs, ys, sigma, prior_mean, prior_sd) that returns a tuple of four values:

  1. the least-squares (maximum likelihood) slope, or None if every x is 0,
  2. the slope where the posterior density is highest (the MAP estimate),
  3. the weight lam for which the MAP slope is also the w that minimises sum((y - w * x)**2) + lam * (w - prior_mean)**2,
  4. the standard deviation of the posterior distribution of w (it is normal).

The setup provides ad_trials(n, effect, noise, seed), which returns (xs, ys) for n simulated stores with true slope effect.

Examples

Input:  xs = [1, 2, -1, 0.5], ys = [0.9, 2.6, -0.4, 0.2], sigma = 1.0, prior_mean = 0.0, prior_sd = 0.5
Output: (1.056, 0.6439024390243903, 4.0, 0.31234752377721214)
Explanation: four stores suggest a slope of 1.056, but the prior says slopes that
large are unusual, and the MAP estimate is pulled down to 0.644.

Input:  xs = [0, 0], ys = [1.0, -3.0], sigma = 2.0, prior_mean = 0.3, prior_sd = 1.0
Output: (None, 0.3, 4.0, 1.0)
Explanation: with no change in spend the data say nothing about w, and the
posterior is the prior.

Constraints

  • 1 <= len(xs) == len(ys) <= 10**5
  • sigma > 0, prior_sd > 0
  • floats are compared with a tolerance of 1e-6

Goals

  • Combine a normal likelihood for a slope with a normal prior and find the posterior's peak
  • Show that the MAP estimate minimises squared error plus a penalty, and find the penalty's weight
  • See how the prior's strength pulls a noisy estimate towards the prior mean
Starting Python…