¿Cómo generar números aleatorios correlacionados (medias, variaciones y grado de correlación dados)?

Lo siento si esto parece demasiado básico, pero supongo que solo estoy buscando confirmar la comprensión aquí. Tengo la sensación de que tendría que hacer esto en dos pasos, y he comenzado a tratar de asimilar las matrices de correlación, pero está empezando a parecer realmente complicado. Estoy buscando una explicación concisa (idealmente con sugerencias para una solución de pseudocódigo) de una forma buena, idealmente rápida, de generar números aleatorios correlacionados.

Dadas dos variables pseudoaleatorias de altura y peso con medias y variaciones conocidas, y una correlación dada, creo que básicamente estoy tratando de entender cómo debería ser este segundo paso:

   height = gaussianPdf(height.mean, height.variance)
   weight = gaussianPdf(correlated_mean(height.mean, correlation_coefficient), 
                        correlated_variance(height.variance, 
                        correlation_coefficient))

• ¿Cómo calculo la media y la varianza correlacionadas? Pero quiero confirmar que ese es realmente el problema relevante aquí.
• ¿Necesito recurrir a la manipulación matricial? ¿O tengo algo más muy mal en mi enfoque básico de este problema?

Para responder a su pregunta sobre "una forma buena, idealmente rápida, de generar números aleatorios correlacionados": Dada una matriz varianza-covarianza deseada que es, por definición, positiva definida, la descomposición de Cholesky es: C = L L T ; L siendo matriz triangular inferior.$C$$C$$LL^T$$L$

Si ahora usa esta matriz para proyectar un vector variable aleatorio no correlacionado X , la proyección resultante Y = L X será la de las variables aleatorias correlacionadas.$L$$X$$Y = LX$

Puede encontrar una explicación concisa de por qué sucede esto aquí .

+1 a @ user11852 y @ jem77bfp, estas son buenas respuestas. Permítanme abordar esto desde una perspectiva diferente, no porque piense que es necesariamente mejor en la práctica , sino porque creo que es instructivo. Aquí hay algunos hechos relevantes que ya conocemos:

1. es la pendiente de la línea de regresión cuando X e Y estánestandarizados, es decir, N ( 0 , 1 ) , $r$$X$$Y$$\mathcal N(0,1)$

2. $r^2$$Y$$X$

   
   
   (también, de las reglas para las variaciones ):

3. 
   $$\text{Var}[aX]=a^2\text{Var}[X]$$
4. 
   $$\text{Var}[X+\varepsilon]=\text{Var}[X]+\text{Var}[\varepsilon]$$

$r$$\rho$$X$$\mathcal N(0,1)$$a$$v_e$$Y \sim\mathcal N(0,a^2+v_e)$$a^2+v_e=1$$|a|$ $\le 1$$a=r$$r$$a$$1-r^2$$x_i$$e_i$$v_e$$y_i$

Si desea hacer esto en R, el siguiente código podría funcionar para usted:

correlatedValue = function(x, r){
  r2 = r**2
  ve = 1-r2
  SD = sqrt(ve)
  e  = rnorm(length(x), mean=0, sd=SD)
  y  = r*x + e
  return(y)
}

set.seed(5)
x = rnorm(10000)
y = correlatedValue(x=x, r=.5)

cor(x,y)
[1] 0.4945964

$0$

Nuevamente, esto, en su forma más simple, solo le permite generar un par de variables correlacionadas (esto podría ampliarse, pero se vuelve feo rápido), y ciertamente no es la forma más conveniente de hacer el trabajo. En R, desearía usar ? Mvrnorm en el paquete MASS , tanto porque es más fácil como porque puede generar muchas variables con una matriz de correlación de población dada. No obstante, creo que vale la pena haber recorrido este proceso para ver cómo se desarrollan algunos principios básicos de una manera simple.

En general, esto no es algo simple de hacer, pero creo que hay paquetes para la generación de variables normales multivariadas (al menos en R, ver mvrnormen el MASSpaquete), donde simplemente ingresas una matriz de covarianza y un vector medio.

$(X_1,X_2)$$F(x_1,x_2)$$F$$x_2$

$F(x_1,x_2)$$x_1=\hat{x}_1$

$\xi_2$$[0,1]$$\xi_1$$F(x_2 | X_1=\hat{x}_1)$$\hat{x}_2=(F(x_2 | X_1=\hat{x}_1))^{-1}(\xi)$$\hat x_2$$F(\hat x_2 | X_1=\hat{x}_1) = \xi$

Si no comprende el significado de conectar una variable uniforme en una función de distribución de probabilidad inversa, intente hacer un bosquejo del caso univariado y luego recuerde cuál es la interpretación geométrica de la función inversa.

Si está listo para renunciar a la eficiencia, puede usar un alogoritmo desechable. Su ventaja es que permite cualquier tipo de distribución (no solo gaussiana).

$\{x_i\}_{i=1}^N$$\{y_i\}_{i=1}^N$$C$

$c_{old}=corr(\{x_i\},\{y_i\})$

$n_1$$n_2: 1 \leq n_{1,2} \leq N$

$x_{n_1}$$x_{n_2}$

$c_{new}=corr( \{x_i\},\{y_i\})$

$|C-c_{new}| < |C-c_{old}|$

$|C-c| < \epsilon$

${x_i}$

¡Buena suerte!