# Covariance matrix of multivariate multiple regression coefficients

I would like to perform a regression analysis on a dataset comprising one independent variable (X) and two dependent variables (Y1 and Y2) which may be affected by correlated errors. R's stats::lm function handles this situation nicely, producing a 4x4 covariance matrix of the regression parameters (slopes and intercepts). I noticed that the correlation coefficient of the intercepts equals the correlation coefficient of the residuals. Is this true for all other fitting models, including Generalised Linear Models? Am I correct in thinking that:

1. The regression coefficients for the dependent variable may be obtained by fitting each of them separately?
2. The covariance/correlation matrix of the intercepts may always be obtained from the residuals?

Edit:

Here's an example:

library(MASS)

a1 <- 5  # intercept 1
a2 <- 10 # intercept 2
b1 <- 1  # slope 1
b2 <- 2  # slope 2
x <- 1:10 # independent variable
e <- mvrnorm(n=length(x),mu=c(0,0),Sigma=matrix(c(1,0.5,0.5,1),nrow=2)) # residuals

# generating some synthetic data:
y1 <- a1 + b1 * x + e[,1] # dependent variable 1
y2 <- a2 + b2 * x + e[,2] # dependent variable 2
y <- cbind(y1,y2)

# linear regression
fit <- lm(y ~ x)


Which produces the following output:

> fit

Call:
lm(formula = y ~ x)

Coefficients:
y1     y2
(Intercept)  4.603  9.591
x            1.095  2.066

> vcov(fit)
y1:(Intercept)         y1:x y2:(Intercept)         y2:x
y1:(Intercept)     0.39622492 -0.056603560     0.35146659 -0.050209513
y1:x              -0.05660356  0.010291556    -0.05020951  0.009129002
y2:(Intercept)     0.35146659 -0.050209513     0.44698258 -0.063854654
y2:x              -0.05020951  0.009129002    -0.06385465  0.011609937


You can see that stats::lm has correctly retrieved both the regression coefficients and the covariance matrix of the regression coefficients. However, stats::glm does not accept multiple dependent variables. So I was wondering if I could calculate their covariance structure from the residuals as well.

• The answer to #1 is yes. The #2 is much more interesting (and I don't know the answer). – amoeba Apr 26 '16 at 22:24
• I don't follow this. What "4x4 covariance matrix" are you referring to, eg? Is that from the vcov()? What do you mean that the "correlation coefficient of the intercepts equals the correlation coefficient of the residuals"? Can you paste in your code & output, or otherwise clarify what you are referring to? – gung Apr 26 '16 at 22:24
• @amoeba, I wouldn't necessarily say that #1 is appropriate. – gung Apr 26 '16 at 22:25
• @gung I have added an example to my post. – PVM Apr 27 '16 at 8:10