Explicar como o `eigen` ajuda a inverter uma matriz

13

Minha pergunta se refere a uma técnica de computação explorada em geoR:::.negloglik.GRFou geoR:::solve.geoR.

Em uma configuração linear modelo misto: onde e são os efeitos fixos e aleatórios, respectivamente. Além disso,β b Σ = cov ( Y )

Y=Xβ+Zb+e
βbΣ=cov(Y)

Ao estimar os efeitos, é necessário calcular que normalmente pode ser feito usando algo como , mas às vezes é quase invertível, portanto, use o truque( X Σ - 1 X )

(XΣ-1X)-1XΣ-1Y
solve(XtS_invX,XtS_invY)(XΣ-1X)geoR
t.ei=eigen(XtS_invX)
crossprod(t(t.ei$vec)/sqrt(t.ei$val))%*%XtS_invY

(pode ser visto em geoR:::.negloglik.GRFe geoR:::.solve.geoR) que equivale a decomposição onde e, portanto,

(XΣ-1X)=ΛDΛ-1
Λ=Λ-1
(XΣ1X)1=(D1/2Λ1)(D1/2Λ1)

Duas questões:

  1. Como essa decomposição de eigen ajuda a inverter ?(XΣ1X)
  2. Existem outras alternativas viáveis ​​(robustas e estáveis)? (por exemplo, qr.solveou chol2inv?)
qoheleth
fonte

Respostas:

15

1) A composição eigend realmente não ajuda muito. É certamente mais numericamente estável do que uma fatoração de Cholesky, o que é útil se sua matriz estiver mal condicionada / quase singular / tiver um número de condição alto. Assim, você pode usar a composição de eigend e ela fornecerá uma solução para o seu problema. Mas há pouca garantia de que será a solução CERTA . Honestamente, quando você inverte explicitamente , o dano já está feito. A formação de X T Σ - 1 X apenas torna as coisas piores. A composição do eigend ajudará você a vencer a batalha, mas a guerra certamente está perdida.ΣXTΣ-1X

2) Sem saber as especificidades do seu problema, é isso que eu faria. Em primeiro lugar, realizar uma fatoração de Cholesky sobre de modo que Σ = L L T . Em seguida, realizar um fatoração QR em G - 1 X de modo a que L - 1 X = Q R . Por favor, certifique-se de computar L - 1 X utilizando a substituição para a frente - NÃO explicitamente invertido L . Então você obtém: X T Σ - 1 X = X T ( LΣΣ=eueuTeu-1Xeu-1X=QReu-1Xeu A partir daqui, você pode resolver qualquer lado direito que desejar. Mas, novamente, por favor, não explicitamente invertidoR(ouRTR). Use substituições para frente e para trás, conforme necessário.

XTΣ-1X=XT(eueuT)-1X=XTeu-Teu-1X=(eu-1X)T(eu-1X)=(QR)TQR=RTQTQT=RTR
RRTR

BTW, estou curioso sobre o lado direito da sua equação. Você escreveu que é . Tem certeza de que não é X T Σ - 1 Y ? Porque, se fosse, você poderia usar um truque semelhante no lado direito: X T Σ - 1 Y = X T ( L L T ) - 1 YXTΣYXTΣ-1Y E então você pode entregar o golpe de misericórdia ao resolverβ: X T Σ - 1 X β = X T Σ - 1 Y R T R β = R T Q T L - 1 Y R β = Q T L - 1 Y β = R - 1 Q T L

XTΣ-1Y=XT(eueuT)-1Y=XTeu-Teu-1Y=(eu-1X)Teu-1Y=(QR)Teu-1Y=RTQTeu-1Y
β Claro, você nunca inverteria explicitamenteRpara a etapa final, certo? Isso é apenas uma substituição para trás. :-)
XTΣ-1Xβ=XTΣ-1YRTRβ=RTQTeu-1YRβ=QTeu-1Yβ=R-1QTeu-1Y
R
Bill Woessner
fonte
Obrigado. Esta é uma resposta útil. Só para ser explícito, a sua alternativa chol / qr vai ajudar a vencer a guerra? ou apenas ganhar o jogo melhor do que o eigen faz?
Qoheleth 30/10/2015
XΣXTΣ-1X