Generation of a random vector on an affine hyperplane

Question

I would like to design a proposal of the form:

$$ p(t=(t_i)|\hat{t}=(\hat{t}_i)) $$

where $t$ (and $\hat{t}$) lies in an affine hyperplane $T \subset R^n$:

$$ t \in T \Leftrightarrow \sum_i t_i=1 $$

Ideally I would like to approach something like: $$ (t_i)| (\hat{t}_i) \sim \mathcal{N}((\hat{t}_i),\sigma^2 I) $$ with the constraint that $\sum_i t_i=1$ (or anything in a similar flavor).

An intuitive strategy could result in the following algorithm:

1) sample individually $\tilde{t}_i \sim \mathcal{N}(\hat{t_i},\sigma)$

2) then renormalise $t_i = \frac{\tilde{t}_i}{\sum_j \tilde{t}_j}$

but it seems to me a bit difficult to characterise the resulting proposal (and particularly to compute $\frac{p(t=(t_i)|\hat{t}=(\hat{t}_i))}{p(\hat{t}=(\hat{t}_i))|t=(t_i))}$ as I am applying it within a metropolis hasting algorithm).

Is anybody has an alternative proposal or some idea on how to deal with $\frac{p(t=(t_i)|\hat{t}=(\hat{t}_i))}{p(\hat{t}=(\hat{t}_i))|t=(t_i))}$ ?

Note: I deal with large $n$ ($T \subset R^n$) and thus need a quite efficient strategy.

Answer 1

There are two general approaches: the implicit and the parametric. Stéphane Laurent describes the implicit method in a comment: generate a $N(\mu, \sigma^2\mathbb{I}_n)$ variate in $\mathbb{R}^n$ and project that orthogonally onto the affine subset. This works due to the spherical symmetry of the distribution in $\mathbb{R}^n$. A little computation is needed to choose $\mu$ correctly.

The implicit approach will be inefficient when $n$ is much larger than the dimension of the affine subspace. It also does not readily generalize to other multivariate normal distributions without a bit of work--and that work will be more extensive than the parametric method.

The parametric method chooses an origin $O$ (intended to be the mean) and orthonormal basis $(e_1, \ldots, e_d)$ for the affine subset $\mathbb{A}^d\in \mathbb{R}^n$. Generate a $N(0, \sigma^2\mathbb{I}_d)$ variate $y = (y_1, y_2, \ldots, y_d)$ and form

$$Y = y_1 e_1 + y_2 e_2 + \cdots + y_d e_d + O.$$

It should be obvious that this works and produces the desired result. Moreover, any desired covariance structure can be induced on this distribution by changing $\sigma^2\mathbb{I}_d$ to an arbitrary positive-definite symmetric matrix $\Sigma$; nothing else changes. Conceptually, all we have done is taken a $d$-dimensional multivariate normal distribution and moved it into $n$ dimensions via rotation and translation.

(The parametric method, then, easily generalizes to non-spherical and non-normal distributions whereas the implicit method does not, because the implicit method involves irrelevant random components distributed within the remaining $n-d$ dimensions that can make it difficult or impossible to achieve the desired distribution on $\mathbb{A}^d$.)

If (*) $\mathbb{A}^d$ is presented to you as the span of a set of vectors (and an origin), then the Gram-Schmidt process directly produces the $(e_i)$. Otherwise, if $\mathbb{A}^d$ is given as the solution of a set of linear equations (and an origin), then various well-known, widely-implemented algorithms will find a set of vectors spanning the kernel (null space) of those equations, putting us into the first situation (*).