In Bayesian model selection (and parameter estimation) one is often faced with difficult, high-dimensional integrals. For example, suppose we have two models for our data, M1(x1, x2) and M2 (X1, x2), both of which are two-dimensional with model parameters x1 and x2. To find out which model does a better job at describing our dataset, within a Bayesian framework (and assuming no prior preference for either model) we would compute the Bayesian evidence factors, Z1 and Z2, where

Z_i = ∫ ζ (x₁, x₂; M_i) Π(x₁, x₂; M_i) dx₁dx₂.     (8)

Here ζ(x1, x2; Mi) is the likelihood function and Π(x1, x2; Mi) is the prior distribution of parameters. If Z1 > Z2 we will prefer model 1 over model 2, for example. Here, the semi-colon notation in ζ(x1, x2; Mi) means that the model Mi(xi, x2) has been assumed in the computation of the likelihood function. How this is done is beyond the scope of this class. Instead, we will assume the likelihood function is given. For simplicity we shall assume ir(xi, x2; Mi) = 1. (The prior must satisfy f it = 1, so strictly speaking we cannot set it = 1 in general.)

In this problem you will tackle a challenging case when the parameters are degenerate: that is, there is a special combination if x1 and x2 where the model makes the same prediction while simultaneously varying the values of both x1 and x2. In such situations, neural networks trying to learn the likelihood surface may also experience significant difficulties. Obviously, the problem gets even worst when the dimensionality of the parameter space grows beyond two.

In the presence of a model degeneracy, one may encounter a likelihood function such as

ζ₁(x1, x2) = exp (- (1 - x₁)² - 100 (x₂ - x₁²)²).

(a) Plot In (L1 (x1, x2)) over the region x1, X2 ∈ [-5, 5]. The model's degeneracy will show up as lines in the x1-x2 plane where the value of ln ζ₁ does not change.

(b) Modify the code you wrote to solve problem 2 to compute the value of

z₁ = _-5∫⁵_-5∫⁵ ζ₁(x₁, x₂) dx₁ dx₂

(c) Let N = 10^k for k = 2,3, 4, and compute In (Z₁) by continually increasing N. Continue increasing N until you can convincingly get at least 2 digits of precision. Quote some values of N and In (Z₁).

(d) By any means possible (using your code or some other way entirely), compute Z₂ for likelihood function ζ₂(x₁, x₂) = 1.

What fits the data better, model 1 or model 2?