Given the normal target Np (0,(1 − ρ)I + ρJ): a. Write a Gibbs sampler using the...

Given the normal target Np (0,(1 − ρ)I + ρJ):

a. Write a Gibbs sampler using the conditional distributions provided in Example 7.4. Run your R code for p = 5 and Vρ = .25, and verify graphically that the marginals are all N (0, 1).

b. Compare your algorithm using T = 500 iterations with rmnorm described in Section 2.2.1 in terms of execution time.

c. Propose a constrained subset that is not a hypercube, and derive the corresponding Gibbs sampler.

Example 7.4
As an extension of Example 7.1, consider the multivariate normal

where I is the p × p identity matrix and J is a p × p matrix of ones. This is a model for equicorrelation, as corr(Xi, Xj) = ρ for every i and j. Using standard formulas for the conditional distributions of a multivariate normal random variable (see, for example, Johnson and Wichern, 1988), it is straightforward but tedious to verify that

where x(−i) = (x1, x2,…,xi−1, xi+1,…, p) and (−i) is the mean of this vector. The Gibbs sampler that generates from these univariate normals can then be easily derived, although it is useless for this problem (Exercise 7.5). It is, however, a short step to consider the setup where the components of the normal vector are restricted to a subset of R p. If this subset is a hypercube,