Finding the marginal posterior distribution of future prediction, $y_{n+1}$

Question

Assume the following bivariate regression model:

$y_i = \beta x_i + u_i$ where $u_i$ is i.i.d $N(0, \sigma^2 = 9)$ for $i = 1, 2, ..., n$.

Assume a noninformative prior of the form:

$p(\beta) \propto constant$, then it can be shown that the posterior pdf for $\beta$ is:

$p(\beta|\mathbf{y}) = (18\pi)^{-\frac{1}{2}}\left(\sum_{i=1}^n x_i^2\right)^{\frac{1}{2}} \exp\left[-\frac{1}{18}\sum_{i=1}^n x_i^2 (\beta - \hat{\beta})^2\right]$

where $\displaystyle{\hat{\beta} = \frac{\sum_{i=1}^n y_ix_i}{\sum_{i=1}^n x_i^2}}$

Now consider the value of $y$ with a given future value of $x$, $x_{n+1}$:

$y_{n+1} = \beta x_{n+1} + u_{n+1}$ where u$_{n+1}$ is i.i.d $N(0, \sigma^2 = 9)$ , show that:

$p(y_{n+1}|x_{n+1},\mathbf{y}) = \int_{\beta} p(y_{n+1}|x_{n+1}, \beta, \mathbf{y}) p(\beta|\mathbf{y})d\beta$ is a normal density with:

$E[y_{n+1}|x_{n+1},\mathbf{y}] = \hat{\beta}x_{n+1}$

and

$\displaystyle{var[y_{n+1}|x_{n+1},\mathbf{y}] = \frac{9[x_{n+1}^2 + \sum_{i=1}^n x_i^2]}{\sum_{i=1}^n x_i^2}}$

Now my approach to this question is like this:

We can derive $p(y_{n+1}|x_{n+1}, \beta, \mathbf{y}) = (18\pi)^{-\frac{1}{2}} \exp\left[-\frac{1}{18}(y_{n+1} - \beta x_{n+1})^2\right]$

So $p(y_{n+1}|x_{n+1}, \beta, \mathbf{y}) \propto \exp\left[-\frac{1}{18}(y_{n+1} - \beta x_{n+1})^2\right]$

and $p(\beta|\mathbf{y}) \propto \exp\left[-\frac{1}{18}\sum_{i=1}^n x_i^2 (\beta - \hat{\beta})^2\right]$

Hence: $p(y_{n+1}|x_{n+1},\mathbf{y}) \propto \int_{\beta} \exp\left[-\frac{1}{18}(y_{n+1} - \beta x_{n+1})^2\right] \exp\left[-\frac{1}{18}\sum_{i=1}^n x_i^2 (\beta - \hat{\beta})^2\right] d\beta$

Now I am stuck on how I should manipulate the integrand to a known distribution so I can recover the integrating constants and hence find an expression for the integral.

Any help would be appreciated.

Answer 1

Expanding on your comment, you can do it like this:

$$ \text{exp}\left[-\frac{1}{18}\left(y^2_{n+1}-2y_{n+1}\beta x_{n+1}+\beta^2x^2_{n+1}+\beta^2\sum x^2_i-2\beta\hat{\beta}\sum x^2_i+\hat{\beta}^2\sum x_i^2\right)\right]\\ \propto \text{exp}\left[-\frac{1}{18}\left(\beta^2(x^2_{n+1}+\sum x^2_i)-2\beta(\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1})\right)\right]\\ =\text{exp}\left[-\frac{1}{18}\left(\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)^2-2\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)\right)\right]\\ =\text{exp}\left[-\frac{1}{18}\left(\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)^2-2\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)+\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)^2-\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)^2\right)\right]\\ \propto\text{exp}\left[-\frac{1}{18}\left(\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)^2-2\left(\beta\sqrt{(x^2_{n+1}+\sum x^2_i)}\right)\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)+\left(\frac{\hat{\beta}\sum x^2_i+y_{n+1}x_{n+1}}{\sqrt{(x^2_{n+1}+\sum x^2_i)}}\right)^2\right)\right]\\ $$ Note that the proportionality here is used in terms of the integral. You need to keep these terms in the density function (as it’s the density for $y$), but not in the integrand.

Let’s call the terms next to the $\beta$ $a$ and the other term $\mu$. The variance is $\sigma^2$. Then the expression is: $$ \text{exp}[-\frac{1}{2\sigma^2}(\beta a - \mu)]= \text{exp}[-\frac{a^2}{2\sigma^2}(\beta -\frac{ \mu}{a})]=)]= \text{exp}[-\frac{1}{2(\sigma a^{-1})^2}(\beta -\frac{ \mu}{a})] $$ Thus, you can rewrite the integrand as the density of a normal distribution (where $\beta$ is the random variable) with mean $\frac{\mu}{a}$ and variance $\sigma^2a^{-2}$.