Usando a decomposição de valor singular para calcular a matriz de covariância de variância a partir do modelo de regressão linear

Eu tenho uma matriz de projeto de p regressores, n observações e estou tentando calcular a matriz de variância-covariância da amostra dos parâmetros. Eu estou tentando calcular diretamente usando svd.

Estou usando R, quando tomo svd da matriz de projeto, recebo três componentes: uma matriz que é , uma matriz que é (presumivelmente valores próprios) e uma matriz que é . Eu diagonalizei , tornando-o uma matriz com zeros nas diagonais externas.$U$$n \times p$$D$$1\times 3$$V$$3\times 3$$D$$3\times 3$

Supostamente, a fórmula da covariância é: , no entanto, a matriz não corresponde, nem é nem perto da função incorporada de R ,. Alguém tem algum conselho / referência? Admito que sou um pouco não qualificado nessa área.$V D^2 V'$vcov

Primeiro, lembre-se de que, sob premissas de normalidade multivariada do modelo de regressão linear, temos que

$$\hat{\beta} \sim \mathcal{N}( \beta, \sigma^2 (X^T X)^{-1} ) .$$

Agora, se , onde o lado direito é a SVD de X, então temos que . Portanto, $X = U D V^T$$X^T X = V D U^T U D V = V D^2 V^T$

$$(X^T X)^{-1} = V D^{-2} V^T .$$

Ainda estamos perdendo a estimativa da variação, que é

$$\hat{\sigma}^2 = \frac{1}{n - p} (y^T y - \hat{\beta}^T X^T y) .$$

Embora eu não tenho verificado, espero vcov retornos .$\hat{\sigma}^2 V D^{-2} V^T$

Nota: Você escreveu , que é , mas precisamos do inverso para a matriz de variância-covariância. Observe também que em , para fazer esse cálculo, você precisa fazer$V D^2 V^T$$X^T X$$R$

vcov.matrix <- var.est * (v %*% d^(-2) %*% t(v))

observando que para multiplicação de matrizes usamos em %*%vez de apenas *. var.estacima é a estimativa da variação do ruído.

(Além disso, fiz as suposições de que é de classificação completa e toda parte. Se esse não for o caso, você precisará fazer pequenas modificações no que foi dito acima.)$X$$n \geq p$