An errors-in-variables model

Problem and solution

Suppose we wish to estimate the regression coefficient for

$$Y_0\mid X_0={\alpha }_0+{\beta }_0{X}_0+{\sigma }_0{ϵ}_0,$$ where $E[ϵ\mid X_0]=0$.

However, we do not observe $Y_0$ and $X_0$. Instead, we observe $$Y_1=Y_0+S_Y{ϵ}_Y$$ for some $\delta$ with $E[\delta \mid Y_0]=0$ for some random variable $S$ and $$X_1=X_0+S_X{ϵ}_X,$$ for some $\eta$ with $E[\delta \mid Y_0]=0$.

As is well known, the regression coefficient ${\beta }_0=\frac{\text{Cov}(X_1,Y_1)}{\text{Var}{X}_1}$. However, $\text{Cov}(Y_1,X_1)$ equals $\text{Cov}(Y_0,X_0)$. About $\text{Var}{X}_1$, we can employ the Law of Total Variance $$\text{Var}{X}_1=E\text{Var}(X\mid S_X)+\text{Var}E(X_0\mid S_X).$$ The term $\text{Var}E(X_0\mid S_X)$ vanishes, as $E(X_0\mid S_X)$ is constant. Moreover, $$\text{Var}{X}_1=\text{Var}{X}_0+\text{Var}{S}_X,$$ hence $$\text{Var}{X}_0=\text{Var}{X}_1-\text{Var}{S}_X.$$

Define the regression model $$Y_1={\alpha }_1+{\beta }_1{X}_1+(S_Y\delta +{\sigma }_0\eta )$$

It follows that $${\beta }_0=\frac{\text{Cov}(Y_1,X_1)}{\text{Var}{X}_1-\text{Var}{S}_X}={\beta }_1\frac{\text{Var}{X}_1}{\text{Var}{X}_1-\text{Var}{S}_X},$$ Moreover, $$\begin{array}{rcl}{\alpha }_0& =& E{Y}_0-{\beta }_{0}E{X}_0,\\ & =& E{Y}_1-{\beta }_1\frac{\text{Var}{X}_1}{\text{Var}{X}_1-\text{Var}{S}_X}E{X}_1.\end{array}$$ If $E{X}_1$ has been normalized to $0$, then ${\alpha }_0=E{Y}_1={\alpha }_1$.

Notice that $Y_1$ has known errors. This makes it – perhaps – possible to estimate ${\beta }_1$ with additional precision, using something similar to weighted least squares. The weights would be $\sqrt{S_Y^2+{\sigma }^2}$. However, as ${\sigma }^2$ is unknown, the resulting regression would not truly be weighted least squares.

Verification

Let’s simulate a bunch of values from the model.

n = 1000000
s_x = sqrt(3)*rexp(n)
s_y = 3*rexp(n)
s_0 = 1

x_0 = rnorm(n, 1, 2)
x_1 = x_0 + s_x * rnorm(n)
y_0 = 0.8 + 0.5 * x_0 + s_0*rnorm(n)
y_1 = y_0 + s_y * rnorm(n)

The calculated coefficients are

beta0_hat = cov(y_1, x_1)/(var(x_1) - 2*var(s_x))
alpha0_hat = mean(y_1) - beta0_hat * mean(x_1)
c(alpha0_hat, beta0_hat)

[1] 0.7964205 0.4951465

But the naive regression $Y_1\sim {\alpha }_1+{\beta }_1{X}_1$ yields

lm(y_1 ~ x_1)

Call:
lm(formula = y_1 ~ x_1)

Coefficients:
(Intercept)          x_1  
     1.0927       0.1988  

On the other hand, the correct (but unobserved) regression yields

lm(y_0~x_0)

Call:
lm(formula = y_0 ~ x_0)

Coefficients:
(Intercept)          x_0  
     0.7987       0.4999  

Inference and literature

To do inference on this method, use the delta method and large-sample theory (together with the studentized bootstrap), or perhaps the bias-corrected accelerated bootstrap (BCa). The delta method should be fairly easy to derive using the formulation of the “covariance of the covariance” foundin e.g. Magnus and Neudecker’s Matrix differential calculus.

There is a sizable literature on error-in-variable models, and inference for this simple model has probably been worked out, but a very rudimentary search yielded nothing for me. I think it’s uncommon to know the variances of the $X$ errors. Moreover, the problem can probably be cast in the language of structrual equation models. But I’m unsure if software (such as lavaan) will help, because you don’t know the item variances in a typical application of structural equations models.

A final option is to assume bivariate normality and use maximum likelihood. This is also likely to be possible using an R package, but I’m not sure the estimates would be consistent. Probably you’d have to use a sandwich matrix for correct standard errors.

To make things easy on yourself, if you’re faced with a problem of this kind, I would suggest just going with the BCa + the equations above. The equations are trivial to compute and BCa will be fairly simple as well; it might be possible to calculate using packages such as bootstrap. Do something else only if the reviewers demand it.