# Covariance between a normal variable (x1) and a sum including it (x1 + x2)

I'm interested in the joint distribution of two variables, $$x_1$$ and $$x_2$$.

$$x_1 \sim Normal(0, \sigma^2_1);\\ \epsilon \sim Normal(0, \sigma^2_{\epsilon});\\ x_2 = x_1 + \epsilon;$$

as a bivariate normal distribution

$$\begin{bmatrix} x_1 \\ x_2 \end{bmatrix} \sim Normal( 0, \Sigma); \\$$

$$\Sigma = \begin{bmatrix} \sigma^2_1 & \rho\sigma_1\sigma_2 \\ \rho\sigma_1\sigma_2 & \sigma^2_2 \end{bmatrix}$$

$$\sigma^2_2 = \sigma^2_1 + \sigma^2_{\epsilon}$$

I can simulate from this distribution quite easily, but would like the closed-form solution for the correlation coefficient $$\rho$$. I can see that the $$\rho$$ should be proportional to $$\sigma^2_1$$: with higher $$\sigma^2_1$$, the signal present in $$x_1$$ is stronger and less influenced by the noise $$\epsilon$$.

Simulations show that the covariance is actually very close to the variance $$\sigma^2_1$$...

library(tidyverse)
vals = 2^seq(1, 8)
vals = seq(1, 10, 1)
d = expand.grid(var1=vals, var2=vals)

d$cov = pmap_dbl(d, function(var1, var2){ # Generate random variables and check covariance between x1 and x3 (= x1 + x2) n = 1000 x1 = rnorm(n, 0, sqrt(var1)) x2 = rnorm(n, 0, sqrt(var2)) x3 = x1 + x2 cov(x1, x3) }) ggplot(d, aes(var1, cov, color=factor(var2))) + geom_path() + geom_abline(intercept=0, slope=1, size=1, linetype='dashed') + labs(x='Var(x1)', y='Cov(x1, x2)', color='Var(x2)') + coord_equal() jet.colors <- c("#00007F", "blue", "#007FFF", "cyan", "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000") ggplot(d, aes(var1, var2, fill=cov)) + geom_tile() + scale_fill_gradientn(colors=jet.colors) + labs(x='Var(x1)', fill='Cov(x1, x2)', y='Var(x2)') + coord_equal()  ...but not quite right either, presumably as I need to incorporate $$\sigma^2_{\epsilon}$$ into the equation somewhere. Is there a known analytic solution to this? • You can calculate this extremely easily by applying definition of covariance and the correlation coefficient. It's even easier if you use the properties of covariance (which all are trivial consequences of the definition of covariance). – Matthew Gunn Oct 30 '18 at 11:19 • Do you mean by using the identitiy$ \sigma ^{2} (\sum _{i=1}^{n}a_{i}X_{i}) =\sum _{i=1}^{n}a_{i}^{2}\sigma ^{2}(X_{i})+ 2 \sum_{i, j : i\lt j} a_{i}a_{j} cov (X_{i},X_{j}) $? – Eoin Oct 30 '18 at 12:57 • You need the covariance between$x_1$and$x_2$. Then you can get the correlation coefficient given you already have variance of$x_1$and$x_2\$. – user158565 Oct 30 '18 at 14:07
• Thank you. I understand the relationship between correlation and covariance. My issue is that I have a closed-form solution for neither. – Eoin Oct 30 '18 at 21:58

Suppose $$x_1$$ and $$\epsilon$$ are independent, then

$$\mathrm{Cov}\left(\matrix{x_1\\ \epsilon}\right) = \left (\matrix{\sigma_1^2 & 0 \\0&\sigma_\epsilon^2}\right)$$ $$\left(\matrix{x_1\\ x_2}\right) = \left (\matrix{1 & 0 \\1 &1}\right)\left(\matrix{x_1\\ \epsilon}\right)$$

So $$\mathrm{Cov}\left(\matrix{x_1\\ x_2}\right) = \left (\matrix{1 & 0 \\1 &1}\right)\left (\matrix{\sigma_1^2 & 0 \\0&\sigma_\epsilon^2}\right)\left (\matrix{1 & 1 \\0 &1}\right) = \left (\matrix{\sigma_1^2 & \sigma_1^2 \\\sigma_1^2 &\sigma_1^2+\sigma_\epsilon^2}\right)$$

So the correlation coefficient $$ρ = \frac {\sigma_1^2}{\sqrt{\sigma_1^2 (\sigma_1^2+\sigma_\epsilon^2)}} = \frac {\sigma_1}{\sqrt{\sigma_1^2+\sigma_\epsilon^2}}$$

• Excellent, thank you. I've also realised why my simulation results didn't quite match this - stupid sampling error. I've posted the corrected simulations below. – Eoin Oct 31 '18 at 9:50

Thank you to @a_statistician for the correct analytical solution. I've since realised why my simulation results haven't agreed with this analytic result (cov(x1, x2) == var(x1)): the effect of sampling error is surprisingly strong here, even with large N, and so the correlation between $$x_1$$ and $$\epsilon$$ weren't necessarily 0, and the sample variances of $$x_1$$ and $$\epsilon$$ weren't necessarily their population values of $$\sigma^2_1$$ and $$\sigma^2_{\epsilon}$$.

Here's the corrected simulation code.

library(tidyverse)
orthogonalise = function(x1, x2){
X = data.frame(x1, x2) %>% as.matrix
pca = prcomp(X)
pX = X %*% pca$$rotation return(list(x1=pX[,1], x2=pX[,2])) } exact.orthogonal.samples = function(var1, var2){ n = 1000 x1 = rnorm(n, 0, 1) x2 = rnorm(n, 0, 1) # Remove correlation ox = orthogonalise(x1, x2) x1 = oxx1 x2 = ox$$x2
# Set variance
x1 = x1 * (sqrt(var1)/sd(x1))
x2 = x2 * (sqrt(var2)/sd(x2))
return(list(x1=x1, x2=x2))
}

vals = 2^seq(1, 8)
vals = seq(1, 10, 1)
d = expand.grid(var1=vals, var2=vals)
d$$cov = res = pmap_dbl(d, function(var1, var2){ # Generate random variables and check covariance between x1 and x2 (= x1 + e) X = exact.orthogonal.samples(var1, var2) x1 = Xx1 e = X$$x2
x2 = x1 + e
cov(x1, x2)
})

jet.colors <- c("#00007F", "blue", "#007FFF", "cyan", "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000")
ggplot(d, aes(var1, var2, fill=cov)) +
geom_tile() +