# How to combine a known and unknown covariance in a linear model

I am faced with a problem where I have a situation that can be fit with a linear multivariate regression type model. But, there are two types residual error that will occur from the lack of fit.

Right now, my $\mathbf{Y}$ data is uncorrupted by noise source $\mathbf{R}_{2}$, but will be in the "future". Therefore, I need to incorporate the extra uncertainty it will cause into my test statistic.

Is this a multivariate random effects model? How do I deal with these two sources of variance?

So far, what I have done is as follows. I can write:

$\mathbf{Y} = \mathbf{X}\mathbf{B} + \mathbf{R}_{1} + \mathbf{R}_{2}$,

where the $\mathbf{R}$'s are a matrices of the residuals, such that of course,

$\mathbf{Y} - \mathbf{X}\mathbf{B} = \mathbf{R}_{1} + \mathbf{R}_{2} = 0$ if I have a perfectly fitting model.

So, the distributional assumption on each $\mathbf{R}$ is $\mathbf{R} \sim N(\mathbf{0},\mathbf{I} \otimes \mathbf{\Sigma_{R}})$.

Except, $\mathbf{\Sigma_{R_{1}}}$ is unknown while $\mathbf{\Sigma_{R_{2}}}$ is known. I can get an unbiased estimate of $\mathbf{\Sigma_{R_{1}}}$ by computing

$\mathbf{\hat{\Sigma}_{R_{1}}} = \frac{1}{n-p} \mathbf{Y} [\mathbf{I} - \mathbf{X}(\mathbf{X}^{'}\mathbf{X})^{-1}\mathbf{X}^{'}]\mathbf{Y}^{'}$.

This looks odd, but it is because $\mathbf{Y}$ has "not yet" been corrupted by the noise due to the known $\mathbf{R}_{2}$. However, when I compute test statistics based on this model (eg. null hypothesis: $\mathbf{B}_{1,2} = 0$ and so on), I want this "future" noise to be taken into consideration, where its covariance structure is known. aka, I have at this moment the matrix $\mathbf{\Sigma_{R_{2}}}$ on hand.

So, how can I combine an known and unknown covariance matrix, one I know, and one I can easily estimate?

• By "..be estimated in the standard way" do you mean that the covariance of $n$ can be estimated separately from the data? If not, then exactly what do you mean by "the standard way"? – whuber Jul 28 '13 at 22:13
• I will rewrite the question. Its actually a multivariate model, where the covariace I'm referring to is the analog of the s^2 in a multiple regression. I just thought it would be simpler to ask in the univariate case and I could generalize. – bill_e Jul 28 '13 at 22:44
• My question is much more specific now – bill_e Jul 28 '13 at 23:09
• @PeterRabbit, Can you define a clear cut between future and past? e.g. For your samples, $y_i$, you know when $i > t$, the noise source $R_2$ will play a role. – zhanxw Aug 5 '13 at 3:50
• No, what I mean by 'future' isn't anything like that. The scenario might be like: a month from now when some process is implemented, that noise source will play a role. There's no reason the statistical analysis can't be done now, since the covariance structure of that noise is known exactly. I'm just unsure of how to proceed. – bill_e Aug 5 '13 at 8:06

There are two answers to your question. The first has to do with a Monte Carlo approach. It is along the lines of your answer to my comment-question. For this, you just need to know the following: don't generate just one sample of $R_2$, generate many. I believe you can fill in the rest.
What follows is the second answer:
In the future, your random disturbance will be $U=R_1 + R_2$. I presume that $R_1$ and $R_2$ are independent. Then $\Sigma_U = \Sigma_{R_{1}} + \Sigma_{R_{2}}$. Given your current sample, you can estimate $\hat \Sigma_{R_{1}}$ separately, and you are very lucky to know the true $\Sigma_{R_{2}}$. Then a perfectly valid estimate of $\Sigma_U$ is $\hat \Sigma_U= \hat \Sigma_{R_{1}}+\Sigma_{R_{2}}$. What you don't have is estimates $\hat B^*$ that would be derived using some $Y^*=Y+R_2$ (if it was available). So it appears somehow "internally inconsistent" to use $\hat \Sigma_U$ together with $\hat B =\left(X'X\right)^{-1}X'Y$ for hypothesis testing and statistical inference in general. Ideally you would want to have available estimates $\hat B^* =\left(X'X\right)^{-1}X'Y^*$, to use with $\hat \Sigma_U$... Ideally yes, -but is it erroneous to do statistical inference by using $\hat B$ and $\hat \Sigma_U$? The answer is, thankfully, no.

First let's consider what exactly would be the difference between $\hat B$ and $\hat B^*$:

$$\hat B^* =\left(X'X\right)^{-1}X'Y^* = \left(X'X\right)^{-1}X'\left(Y+R_2\right)= \hat B+ \left(X'X\right)^{-1} X' R_2$$ which can be written

$$\hat B^* = \hat B+ \left(\frac1nX'X\right)^{-1} \left(\frac 1nX' R_2\right)$$ where $n$ is the sample size. If your sample size is respectably large, and the usual uncorrelatedness assumption between regressors and error hopefully hold, then the second term will be close to zero and hence we expect $\hat B^* \approx \hat B$ (note that the regressor matrix is the same). But this is the practical approach to the matter. The very helpful fact is, that if we can incorporate more information in our inference, then we can (and should) do it even if this new information won't affect all our estimates. This is fundamental: through statistics we approximate reality, and the approximation error accumulates from different sources - we do not have to wipe out all approximation error, or do nothing at all. We can partially improve by eliminating one source of it. Think of the following cases:
A) We suspect heteroskedasticity. One way to proceed is to replace the estimated covariance matrix with a "heteroskedasticity-robust" covariance matrix. But this is not the matrix the OLS estimator gave us. Still, we can use it.
B) Assume no new future disturbance will appear. You estimate your current sample, and you obtain $\hat B$ and $\hat \Sigma_{R_{1}}$. Then somebody comes along and reveals to you the true $\Sigma_{R_{1}}$. You would use $\Sigma_{R_{1}}$ instead of $\hat \Sigma_{R_{1}}$ in your hypothesis testing, because then your t-statistic would have an exact standard normal distribution: you are better off.

Your question describes essentially an equivalent situation to these cases, as regards "what coefficient estimates and what covariance matrix". So you can go ahead and use $\hat \Sigma_U= \hat \Sigma_{R_{1}}+\Sigma_{R_{2}}$ together with $\hat B$.

• Awesome answer, just what I was looking for. Thank you so much for the careful thought. As a brief follow-up, I made estimates for $\hat{B}$ and $\hat{\Sigma_{U}}$, and then computed the 'theoretical' test statistic for a few hypotheses. I compared these to a monte carlo simulations, and found the mean of the monte carlo results to pretty much equal the 'theoretical' results. So, awesome! thank you again – bill_e Aug 7 '13 at 7:55
• Glad to be of help. – Alecos Papadopoulos Aug 7 '13 at 13:37
• I had a super quick follow-up question.. if the expectation value of $\mathbf{R}_{2}$ is the zero matrix.. then $(X'X)^{-1}X'R_{2} = 0$, right? – bill_e Aug 12 '13 at 23:10
• By all means, since $X$ is exogenous to $R_2$. – Alecos Papadopoulos Aug 13 '13 at 0:08