6

I have to implement a multivariate analysis on $n$ random variables with a sample of $m$ data points. I would like to get a matrix with the $\beta$ (as in $n$ $\beta$ vectors put together).

Is there an analytical formula for that or do I have to calculate each vector individually and put them back together ?

Also is there a formula that would give me directly a vector of the estimation error for each variable?

StasK
  • 29,235
  • 2
  • 80
  • 165
0x26res
  • 61
  • 1
  • 5
  • I have to admit that I don't understand your question. However, your question seems to be related to vectorisation of matrices, and operations on these. You might want to look at http://www.imm.dtu.dk/pubdb/views/edoc_download.php/3274/pdf/imm3274.pdf – queenbee May 15 '14 at 22:04
  • 1
    Can you say more on what you mean by $n \beta$ "put together"? Do you mean finding the sum of their products, as you would to find a fitted value in standard linear regression? To your second question, you'd like the vector to contain the standard error of the $\beta$ coefficients, is that correct? – Sean Easter May 19 '14 at 14:48

3 Answers3

3

I think I understand what you're asking, but correct me if I'm wrong. The analytical formula for $\beta$ is the same for the multivariate case as the univariate case: $$ \hat \beta = (X'X)^{-1}X'Y $$

You find this the same way as for the univariate case, by taking the first derivative of residual sum of squares. It is relatively straightforward to calculate using matrix calculus (which is covered in the matrix cookbook linked to by queenbee). You can test whether this solution works in R:

y <- cbind(rnorm(10), rnorm(10), rnorm(10))

x <- cbind(1, rnorm(10), rnorm(10), rnorm(10),
       rnorm(10), rnorm(10), rnorm(10))
colnames(x) <- paste("x", 1:6, sep = "")
colnames(y) <- paste("y", 1:3, sep = "")

fit <- lm(y ~ x - 1)
summary(fit)

anaSol <- solve((t(x) %*% x)) %*% t(x) %*% y
anaSol

coef(fit) - anaSol

Here's another reference, specifically related to multivariate analysis: http://socserv.mcmaster.ca/jfox/Books/Companion/appendix/Appendix-Multivariate-Linear-Models.pdf

Eric Brady
  • 728
  • 6
  • 9
1

If you have $q$ equations and $p$ independent variables (including a constant) that appear in every equation, the parameter estimates are given by the $p \times q$ matrix:

$$M=(X'IX)^{-1}X'IY$$

where

  • $Y$ is $n \times q$ matrix of dependent variables
  • X is $n \times p$ matrix of covariates
  • I is the identity matrix
dimitriy
  • 31,081
  • 5
  • 63
  • 138
0

Other answers nicely cover how to derive the $\beta$ coefficients. I'm not sure what you mean by $n, \beta$ "put together." But, if it means that you'd like to use the coefficients to derive the model's predicted values using the coefficients, it's simply the product $XB$, where $X$ is an $m \times n$ matrix of $m$ observations, each with $n$ independent variables, and $B$ a $n \times p$ matrix of regression coefficients. ($p$ here is the number of dependent variables.)

To your question about the standard error, for a single independent variable and single coefficient, the formula is:

$s.e.(\beta_j) = \sqrt{s^2 (X'X)^{-1}_{jj} }$

where $s^2$ is the sum of squared residuals, given by $\sum_i y_i -\hat y_i $, over $m - n$. (More here.) To broaden the formula to return a vector of standard errors corresponding to each coefficient:

$s.e.(\beta) = \sqrt{s^2 diag(X'X)^{-1} }$

where $diag$ returns the diagonal of the matrix. To further broaden that formula to return a matrix of standard error corresponding to each $\beta$ coefficient, replace $s^2$ with $S^2$, a row vector containing the $s^2$ for each independent variable:

$s.e.(B) = \sqrt{diag(X'X)^{-1} S^2 }$

Note that the vector returned by $diag$ should be $n \times 1$, and $S^2$ $1 \times p$, making their product $n \times p$, standard errors corresponding to the coefficients in $B$. (Square root again applied element-wise.)

Sean Easter
  • 8,359
  • 2
  • 29
  • 58