Utilizo la descomposición de Cholesky para simular variables aleatorias correlacionadas dada una matriz de correlaciones. El caso es que el resultado nunca reproduce la estructura de correlación tal y como se da. He aquí un pequeño ejemplo en Python para ilustrar la situación.
import numpy as np
n_obs = 10000
means = [1, 2, 3]
sds = [1, 2, 3] # standard deviations
# generating random independent variables
observations = np.vstack([np.random.normal(loc=mean, scale=sd, size=n_obs)
for mean, sd in zip(means, sds)]) # observations, a row per variable
cor_matrix = np.array([[1.0, 0.6, 0.9],
[0.6, 1.0, 0.5],
[0.9, 0.5, 1.0]])
L = np.linalg.cholesky(cor_matrix)
print(np.corrcoef(L.dot(observations)))
Esto imprime:
[[ 1. 0.34450587 0.57515737]
[ 0.34450587 1. 0.1488504 ]
[ 0.57515737 0.1488504 1. ]]
Como puede ver, la matriz de correlaciones estimada post hoc difiere drásticamente de la anterior. Hay algún error en mi código o existe alguna alternativa al uso de la descomposición Cholesky?
Editar
Le pido perdón por este lío. No creía que hubiera un error en el código y/o en la forma de aplicar la descomposición de Cholesky debido a algún malentendido del material que había estudiado antes. De hecho estaba seguro de que el método en sí no pretendía ser preciso y me había parecido bien hasta la situación que me hizo publicar esta pregunta. Gracias por señalarme el error de concepto que tenía. He editado el título para reflejar mejor la situación real propuesta por @Silverfish.