Regresión lineal multivariada bayesiana

Enfoque bayesiano para la regresión lineal multivariada

En estadística , la regresión lineal multivariada bayesiana es un enfoque bayesiano de la regresión lineal multivariada , es decir, la regresión lineal donde el resultado previsto es un vector de variables aleatorias correlacionadas en lugar de una única variable aleatoria escalar. Un tratamiento más general de este enfoque se puede encontrar en el artículo Estimador MMSE .

Detalles

Considere un problema de regresión en el que la variable dependiente que se va a predecir no es un único escalar de valor real sino un vector de longitud m de números reales correlacionados. Como en la configuración de regresión estándar, hay n observaciones, donde cada observación i consta de k −1 variables explicativas , agrupadas en un vector de longitud k (donde se ha agregado una variable ficticia con un valor de 1 para permitir un coeficiente de intersección ). Esto puede verse como un conjunto de m problemas de regresión relacionados para cada observación i : ${\displaystyle \mathbf {x} _{i}}$

$${\displaystyle {\begin{aligned}y_{i,1}&=\mathbf {x} _{i}^{\mathsf {T}}{\boldsymbol {\beta }}_{1}+\epsilon _{i,1}\\&\;\;\vdots \\y_{i,m}&=\mathbf {x} _{i}^{\mathsf {T}}{\boldsymbol {\beta }}_{m}+\epsilon _{i,m}\end{aligned}}}$$

vector de fila ${\displaystyle \{\epsilon _{i,1},\ldots ,\epsilon _{i,m}\}}$ ${\displaystyle \mathbf {y} _{i}^{\mathsf {T}}}$

$${\displaystyle \mathbf {y} _{i}^{\mathsf {T}}=\mathbf {x} _{i}^{\mathsf {T}}\mathbf {B} +{\boldsymbol {\epsilon }}_{i}^{\mathsf {T}}.}$$

La matriz de coeficientes B es una matriz donde los vectores de coeficientes para cada problema de regresión se apilan horizontalmente: ${\displaystyle k\times m}$ ${\displaystyle {\boldsymbol {\beta }}_{1},\ldots ,{\boldsymbol {\beta }}_{m}}$

$${\displaystyle \mathbf {B} ={\begin{bmatrix}{\begin{pmatrix}\\{\boldsymbol {\beta }}_{1}\\\\\end{pmatrix}}\cdots {\begin{pmatrix}\\{\boldsymbol {\beta }}_{m}\\\\\end{pmatrix}}\end{bmatrix}}={\begin{bmatrix}{\begin{pmatrix}\beta _{1,1}\\\vdots \\\beta _{k,1}\end{pmatrix}}\cdots {\begin{pmatrix}\beta _{1,m}\\\vdots \\\beta _{k,m}\end{pmatrix}}\end{bmatrix}}.}$$

El vector de ruido para cada observación i es conjuntamente normal, de modo que los resultados de una observación determinada están correlacionados: ${\displaystyle {\boldsymbol {\epsilon }}_{i}}$

$${\displaystyle {\boldsymbol {\epsilon }}_{i}\sim N(0,{\boldsymbol {\Sigma }}_{\epsilon }).}$$

Podemos escribir todo el problema de regresión en forma matricial como:

$${\displaystyle \mathbf {Y} =\mathbf {X} \mathbf {B} +\mathbf {E} ,}$$

YEmatriz de diseño Xde regresión lineal estándar: ${\displaystyle n\times m}$ ${\displaystyle n\times k}$

$${\displaystyle \mathbf {X} ={\begin{bmatrix}\mathbf {x} _{1}^{\mathsf {T}}\\\mathbf {x} _{2}^{\mathsf {T}}\\\vdots \\\mathbf {x} _{n}^{\mathsf {T}}\end{bmatrix}}={\begin{bmatrix}x_{1,1}&\cdots &x_{1,k}\\x_{2,1}&\cdots &x_{2,k}\\\vdots &\ddots &\vdots \\x_{n,1}&\cdots &x_{n,k}\end{bmatrix}}.}$$

La solución clásica de mínimos cuadrados lineales frecuentistas es simplemente estimar la matriz de coeficientes de regresión utilizando la pseudoinversa de Moore-Penrose : ${\displaystyle {\hat {\mathbf {B} }}}$

$${\displaystyle {\hat {\mathbf {B} }}=(\mathbf {X} ^{\mathsf {T}}\mathbf {X} )^{-1}\mathbf {X} ^{\mathsf {T}}\mathbf {Y} .}$$

Para obtener la solución bayesiana, necesitamos especificar la probabilidad condicional y luego encontrar el conjugado previo apropiado. Al igual que con el caso univariante de regresión bayesiana lineal , encontraremos que podemos especificar un conjugado condicional natural previo (que depende de la escala).

Escribamos nuestra probabilidad condicional como ^[1]

$${\displaystyle \rho (\mathbf {E} |{\boldsymbol {\Sigma }}_{\epsilon })\propto |{\boldsymbol {\Sigma }}_{\epsilon }|^{-n/2}\exp \left(-{\tfrac {1}{2}}\operatorname {tr} \left(\mathbf {E} ^{\mathsf {T}}\mathbf {E} {\boldsymbol {\Sigma }}_{\epsilon }^{-1}\right)\right),}$$

${\displaystyle \mathbf {E} }$ ${\displaystyle \mathbf {Y} ,\mathbf {X} ,}$ ${\displaystyle \mathbf {B} }$

$${\displaystyle \rho (\mathbf {Y} |\mathbf {X} ,\mathbf {B} ,{\boldsymbol {\Sigma }}_{\epsilon })\propto |{\boldsymbol {\Sigma }}_{\epsilon }|^{-n/2}\exp(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {Y} -\mathbf {X} \mathbf {B} )^{\mathsf {T}}(\mathbf {Y} -\mathbf {X} \mathbf {B} ){\boldsymbol {\Sigma }}_{\epsilon }^{-1})),}$$

Buscamos un prior conjugado natural: una densidad conjunta que tenga la misma forma funcional que la verosimilitud. Dado que la probabilidad es cuadrática en , reescribimos la probabilidad para que sea normal en (la desviación de la estimación muestral clásica). ${\displaystyle \rho (\mathbf {B} ,\Sigma _{\epsilon })}$ ${\displaystyle \mathbf {B} }$ ${\displaystyle (\mathbf {B} -{\hat {\mathbf {B} }})}$

Usando la misma técnica que con la regresión lineal bayesiana , descomponemos el término exponencial usando una forma matricial de la técnica de suma de cuadrados. Aquí, sin embargo, también necesitaremos utilizar el cálculo diferencial matricial ( producto de Kronecker y transformaciones de vectorización ).

Primero, apliquemos la suma de cuadrados para obtener una nueva expresión para la probabilidad:

$${\displaystyle \rho (\mathbf {Y} |\mathbf {X} ,\mathbf {B} ,{\boldsymbol {\Sigma }}_{\epsilon })\propto |{\boldsymbol {\Sigma }}_{\epsilon }|^{-(n-k)/2}\exp(-\operatorname {tr} ({\tfrac {1}{2}}\mathbf {S} ^{\mathsf {T}}\mathbf {S} {\boldsymbol {\Sigma }}_{\epsilon }^{-1}))|{\boldsymbol {\Sigma }}_{\epsilon }|^{-k/2}\exp(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {B} -{\hat {\mathbf {B} }})^{\mathsf {T}}\mathbf {X} ^{\mathsf {T}}\mathbf {X} (\mathbf {B} -{\hat {\mathbf {B} }}){\boldsymbol {\Sigma }}_{\epsilon }^{-1})),}$$

$${\displaystyle \mathbf {S} =\mathbf {Y} -\mathbf {X} {\hat {\mathbf {B} }}}$$

Nos gustaría desarrollar una forma condicional para los anteriores:

$${\displaystyle \rho (\mathbf {B} ,{\boldsymbol {\Sigma }}_{\epsilon })=\rho ({\boldsymbol {\Sigma }}_{\epsilon })\rho (\mathbf {B} |{\boldsymbol {\Sigma }}_{\epsilon }),}$$

distribución de Wishart inversadistribución normalde vectorización ${\displaystyle \rho ({\boldsymbol {\Sigma }}_{\epsilon })}$ ${\displaystyle \rho (\mathbf {B} |{\boldsymbol {\Sigma }}_{\epsilon })}$ ${\displaystyle \mathbf {B} }$ ${\displaystyle \mathbf {B} ,{\hat {\mathbf {B} }}}$ ${\displaystyle {\boldsymbol {\beta }}=\operatorname {vec} (\mathbf {B} ),{\hat {\boldsymbol {\beta }}}=\operatorname {vec} ({\hat {\mathbf {B} }})}$

Escribir

$${\displaystyle \operatorname {tr} ((\mathbf {B} -{\hat {\mathbf {B} }})^{\mathsf {T}}\mathbf {X} ^{\mathsf {T}}\mathbf {X} (\mathbf {B} -{\hat {\mathbf {B} }}){\boldsymbol {\Sigma }}_{\epsilon }^{-1})=\operatorname {vec} (\mathbf {B} -{\hat {\mathbf {B} }})^{\mathsf {T}}\operatorname {vec} (\mathbf {X} ^{\mathsf {T}}\mathbf {X} (\mathbf {B} -{\hat {\mathbf {B} }}){\boldsymbol {\Sigma }}_{\epsilon }^{-1})}$$

Dejar

$${\displaystyle \operatorname {vec} (\mathbf {X} ^{\mathsf {T}}\mathbf {X} (\mathbf {B} -{\hat {\mathbf {B} }}){\boldsymbol {\Sigma }}_{\epsilon }^{-1})=({\boldsymbol {\Sigma }}_{\epsilon }^{-1}\otimes \mathbf {X} ^{\mathsf {T}}\mathbf {X} )\operatorname {vec} (\mathbf {B} -{\hat {\mathbf {B} }}),}$$

producto de KroneckerABproducto externo ${\displaystyle \mathbf {A} \otimes \mathbf {B} }$ ${\displaystyle m\times n}$ ${\displaystyle p\times q}$ ${\displaystyle mp\times nq}$

Entonces

$${\displaystyle {\begin{aligned}&\operatorname {vec} (\mathbf {B} -{\hat {\mathbf {B} }})^{\mathsf {T}}({\boldsymbol {\Sigma }}_{\epsilon }^{-1}\otimes \mathbf {X} ^{\mathsf {T}}\mathbf {X} )\operatorname {vec} (\mathbf {B} -{\hat {\mathbf {B} }})\\&=({\boldsymbol {\beta }}-{\hat {\boldsymbol {\beta }}})^{\mathsf {T}}({\boldsymbol {\Sigma }}_{\epsilon }^{-1}\otimes \mathbf {X} ^{\mathsf {T}}\mathbf {X} )({\boldsymbol {\beta }}-{\hat {\boldsymbol {\beta }}})\end{aligned}}}$$

${\displaystyle ({\boldsymbol {\beta }}-{\hat {\boldsymbol {\beta }}})}$

Con la probabilidad en una forma más manejable, ahora podemos encontrar un conjugado previo natural (condicional).

Distribución previa conjugada

El conjugado natural previo usando la variable vectorizada tiene la forma: ^[1] ${\displaystyle {\boldsymbol {\beta }}}$

$${\displaystyle \rho ({\boldsymbol {\beta }},{\boldsymbol {\Sigma }}_{\epsilon })=\rho ({\boldsymbol {\Sigma }}_{\epsilon })\rho ({\boldsymbol {\beta }}|{\boldsymbol {\Sigma }}_{\epsilon }),}$$

$${\displaystyle \rho ({\boldsymbol {\Sigma }}_{\epsilon })\sim {\mathcal {W}}^{-1}(\mathbf {V} _{0},{\boldsymbol {\nu }}_{0})}$$

$${\displaystyle \rho ({\boldsymbol {\beta }}|{\boldsymbol {\Sigma }}_{\epsilon })\sim N({\boldsymbol {\beta }}_{0},{\boldsymbol {\Sigma }}_{\epsilon }\otimes {\boldsymbol {\Lambda }}_{0}^{-1}).}$$

Distribución posterior

Usando lo anterior y la probabilidad, la distribución posterior se puede expresar como: ^[1]

$${\displaystyle {\begin{aligned}\rho ({\boldsymbol {\beta }},{\boldsymbol {\Sigma }}_{\epsilon }|\mathbf {Y} ,\mathbf {X} )\propto {}&|{\boldsymbol {\Sigma }}_{\epsilon }|^{-({\boldsymbol {\nu }}_{0}+m+1)/2}\exp {(-{\tfrac {1}{2}}\operatorname {tr} (\mathbf {V} _{0}{\boldsymbol {\Sigma }}_{\epsilon }^{-1}))}\\&\times |{\boldsymbol {\Sigma }}_{\epsilon }|^{-k/2}\exp {(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {B} -\mathbf {B} _{0})^{\mathsf {T}}{\boldsymbol {\Lambda }}_{0}(\mathbf {B} -\mathbf {B} _{0}){\boldsymbol {\Sigma }}_{\epsilon }^{-1}))}\\&\times |{\boldsymbol {\Sigma }}_{\epsilon }|^{-n/2}\exp {(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {Y} -\mathbf {XB} )^{\mathsf {T}}(\mathbf {Y} -\mathbf {XB} ){\boldsymbol {\Sigma }}_{\epsilon }^{-1}))},\end{aligned}}}$$

${\displaystyle \operatorname {vec} (\mathbf {B} _{0})={\boldsymbol {\beta }}_{0}}$ ${\displaystyle \mathbf {B} }$ ${\displaystyle {\boldsymbol {\Lambda }}_{0}=\mathbf {U} ^{\mathsf {T}}\mathbf {U} }$

$${\displaystyle {\begin{aligned}&\left(\mathbf {B} -\mathbf {B} _{0}\right)^{\mathsf {T}}{\boldsymbol {\Lambda }}_{0}\left(\mathbf {B} -\mathbf {B} _{0}\right)+\left(\mathbf {Y} -\mathbf {XB} \right)^{\mathsf {T}}\left(\mathbf {Y} -\mathbf {XB} \right)\\={}&\left({\begin{bmatrix}\mathbf {Y} \\\mathbf {U} \mathbf {B} _{0}\end{bmatrix}}-{\begin{bmatrix}\mathbf {X} \\\mathbf {U} \end{bmatrix}}\mathbf {B} \right)^{\mathsf {T}}\left({\begin{bmatrix}\mathbf {Y} \\\mathbf {U} \mathbf {B} _{0}\end{bmatrix}}-{\begin{bmatrix}\mathbf {X} \\\mathbf {U} \end{bmatrix}}\mathbf {B} \right)\\={}&\left({\begin{bmatrix}\mathbf {Y} \\\mathbf {U} \mathbf {B} _{0}\end{bmatrix}}-{\begin{bmatrix}\mathbf {X} \\\mathbf {U} \end{bmatrix}}\mathbf {B} _{n}\right)^{\mathsf {T}}\left({\begin{bmatrix}\mathbf {Y} \\\mathbf {U} \mathbf {B} _{0}\end{bmatrix}}-{\begin{bmatrix}\mathbf {X} \\\mathbf {U} \end{bmatrix}}\mathbf {B} _{n}\right)+\left(\mathbf {B} -\mathbf {B} _{n}\right)^{\mathsf {T}}\left(\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0}\right)\left(\mathbf {B} -\mathbf {B} _{n}\right)\\={}&\left(\mathbf {Y} -\mathbf {X} \mathbf {B} _{n}\right)^{\mathsf {T}}\left(\mathbf {Y} -\mathbf {X} \mathbf {B} _{n}\right)+\left(\mathbf {B} _{0}-\mathbf {B} _{n}\right)^{\mathsf {T}}{\boldsymbol {\Lambda }}_{0}\left(\mathbf {B} _{0}-\mathbf {B} _{n}\right)+\left(\mathbf {B} -\mathbf {B} _{n}\right)^{\mathsf {T}}\left(\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0}\right)\left(\mathbf {B} -\mathbf {B} _{n}\right),\end{aligned}}}$$

$${\displaystyle \mathbf {B} _{n}=\left(\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0}\right)^{-1}\left(\mathbf {X} ^{\mathsf {T}}\mathbf {X} {\hat {\mathbf {B} }}+{\boldsymbol {\Lambda }}_{0}\mathbf {B} _{0}\right)=\left(\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0}\right)^{-1}\left(\mathbf {X} ^{\mathsf {T}}\mathbf {Y} +{\boldsymbol {\Lambda }}_{0}\mathbf {B} _{0}\right).}$$

Esto ahora nos permite escribir la parte posterior de una forma más útil:

$${\displaystyle {\begin{aligned}\rho ({\boldsymbol {\beta }},{\boldsymbol {\Sigma }}_{\epsilon }|\mathbf {Y} ,\mathbf {X} )\propto {}&|{\boldsymbol {\Sigma }}_{\epsilon }|^{-({\boldsymbol {\nu }}_{0}+m+n+1)/2}\exp {(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {V} _{0}+(\mathbf {Y} -\mathbf {XB_{n}} )^{\mathsf {T}}(\mathbf {Y} -\mathbf {XB_{n}} )+(\mathbf {B} _{n}-\mathbf {B} _{0})^{\mathsf {T}}{\boldsymbol {\Lambda }}_{0}(\mathbf {B} _{n}-\mathbf {B} _{0})){\boldsymbol {\Sigma }}_{\epsilon }^{-1}))}\\&\times |{\boldsymbol {\Sigma }}_{\epsilon }|^{-k/2}\exp {(-{\tfrac {1}{2}}\operatorname {tr} ((\mathbf {B} -\mathbf {B} _{n})^{\mathsf {T}}(\mathbf {X} ^{T}\mathbf {X} +{\boldsymbol {\Lambda }}_{0})(\mathbf {B} -\mathbf {B} _{n}){\boldsymbol {\Sigma }}_{\epsilon }^{-1}))}.\end{aligned}}}$$

Esto toma la forma de una distribución de Wishart inversa multiplicada por una distribución normal Matrix :

$${\displaystyle \rho ({\boldsymbol {\Sigma }}_{\epsilon }|\mathbf {Y} ,\mathbf {X} )\sim {\mathcal {W}}^{-1}(\mathbf {V} _{n},{\boldsymbol {\nu }}_{n})}$$

$${\displaystyle \rho (\mathbf {B} |\mathbf {Y} ,\mathbf {X} ,{\boldsymbol {\Sigma }}_{\epsilon })\sim {\mathcal {MN}}_{k,m}(\mathbf {B} _{n},{\boldsymbol {\Lambda }}_{n}^{-1},{\boldsymbol {\Sigma }}_{\epsilon }).}$$

Los parámetros de este posterior vienen dados por:

$${\displaystyle \mathbf {V} _{n}=\mathbf {V} _{0}+(\mathbf {Y} -\mathbf {XB_{n}} )^{\mathsf {T}}(\mathbf {Y} -\mathbf {XB_{n}} )+(\mathbf {B} _{n}-\mathbf {B} _{0})^{\mathsf {T}}{\boldsymbol {\Lambda }}_{0}(\mathbf {B} _{n}-\mathbf {B} _{0})}$$

$${\displaystyle {\boldsymbol {\nu }}_{n}={\boldsymbol {\nu }}_{0}+n}$$

$${\displaystyle \mathbf {B} _{n}=(\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0})^{-1}(\mathbf {X} ^{\mathsf {T}}\mathbf {Y} +{\boldsymbol {\Lambda }}_{0}\mathbf {B} _{0})}$$

$${\displaystyle {\boldsymbol {\Lambda }}_{n}=\mathbf {X} ^{\mathsf {T}}\mathbf {X} +{\boldsymbol {\Lambda }}_{0}}$$

Ver también

• Regresión lineal bayesiana
• Distribución normal de la matriz

Referencias

1. Peter E. Rossi, Greg M. Allenby, Rob McCulloch. Estadística y Marketing Bayesianos . John Wiley e hijos, 2012, pág. 32.

• Caja, GEP ; Tiao, GC (1973). "8". Inferencia bayesiana en análisis estadístico . Wiley. ISBN 0-471-57428-7.
• Geisser, S. (1965). "Estimación bayesiana en análisis multivariado". Los anales de la estadística matemática . 36 (1): 150-159. JSTOR  2238083.
• Tiao, GC; Zellner, A. (1964). "Sobre la estimación bayesiana de la regresión multivariada". Revista de la Real Sociedad de Estadística. Serie B (Metodológica) . 26 (2): 277–285. JSTOR  2984424.