Génération de deux nombres pseudo-aléatoires corrélés (avec exemple Python)

Lorsque vous effectuez des expériences numériques statistiques, vous souhaiterez peut-être générer des nombres aléatoires corrélés. Une recherche rapide révélera une méthode utilisant la décomposition choleskey. C'est une façon de gérer n'importe quel nombre de variables, mais je ne voulais que la variable $ 2 $, donc je vais l'expliquer d'une manière un peu plus confuse [^ 1]. [^ 1]: Après tout, le résultat est le même que pour la décomposition en choleskey pour la variable $ 2 $.

Aperçu

Supposons que vous ayez des nombres aléatoires indépendants $ X $ et $ Y $ avec une moyenne nulle et une variance égale. A ce moment, le nombre aléatoire $ Z $ qui est corrélé avec $ X $ par le coefficient de corrélation $ \ rho $ et a la même variance

Z = \rho X + \sqrt{1 - \rho^2} Y

Vous pouvez faire comme ça.

Dans cet article

Je vais écrire sur.

théorie

Pour les nombres aléatoires indépendants $ X $ et $ Y $, la moyenne est zéro et la variance est la même valeur, $ \ sigma ^ 2 $. Puisque $ X $ et $ Y $ sont indépendants, la covariance et le coefficient de corrélation sont nuls.

Aussi, pour le moment, le nombre aléatoire $ Z $ que vous voulez corréler avec $ X $ se superpose à $ X $ et $ Y $.

Z = a X + b Y

Je vais le définir comme. $ a $ et $ b $ sont des constantes.

La raison en est que ce qui suit est vrai pour les exemples extrêmes suivants.

a b \rho
> 0 0 1
=0 \neq 0 0
< 0 0 -1

En regardant ce tableau, j'ai le sentiment qu'un coefficient de corrélation libre peut être créé en modifiant l'équilibre entre $ a $ et $ b $ [^ 2]. Ainsi, l'expression de $ Z $ ci-dessus est prometteuse en tant que nombre aléatoire qui a une corrélation arbitraire avec $ X $.

Calculons maintenant le coefficient de corrélation $ \ rho $ entre $ X $ et $ Z $. Le coefficient de corrélation est

\rho = \frac{\mathrm{Cov}[X, Z]}{\sqrt{\mathrm{Var}[X]} \sqrt{\mathrm{Var}[Z]}}

Il peut être calculé avec. Ici, $ \ mathrm {Cov} $ et $ \ mathrm {Var} $ sont co-distribués et distribués, respectivement. Commençons par calculer à partir de la covariance.

\begin{align}
\mathrm{Cov}[X, Z] &= \mathrm{Cov}[X, a X + b Y] & (\because Z \Définition de)
\\
&= a \ \mathrm{Cov}[X, X] + b \ \mathrm{Cov}[X, Y] & (\because \mathrm{Cov} \Nature de)
\\
&= a \ \sigma^2 & (\because \mathrm{Cov}[X, X] = \mathrm{Var}[X] = \sigma^2, \mathrm{Cov}[X, Y] = 0)
\end{align}

C'est un calcul déroutant, mais si vous ne le comprenez pas, vous devriez pouvoir l'atteindre en revenant à la définition et en calculant calmement. Si c'est un problème, vous n'avez pas à suivre le calcul.

Aussi, concernant la distribution de $ Z $,

\begin{align}
\mathrm{Var}[Z] &= \mathrm{Var}[a X + b Y] & (\because Z \Définition de)
\\
&= a^2 \mathrm{Var}[X] + b^2 \mathrm{Var}[Y] & (\because \mathrm{Var} \Nature de)
\\
&= a^2 \sigma^2 + b^2 \sigma^2 &
\end{align}

Par conséquent, le coefficient de corrélation que je voulais calculer est

\begin{align}
\rho &= \frac{\mathrm{Cov}[X, Z]}{\sqrt{\mathrm{Var}[X]} \sqrt{\mathrm{Var}[Z]}}
\\
&= \frac{a \ \sigma^2}{\sqrt{\sigma^2} \sqrt{a^2 \sigma^2 + b^2 \sigma^2}}
\\
&= \frac{a}{\sqrt{a^2 + b^2}}
\end{align}

Et c'était une très belle cérémonie. Résoudre cette équation pour $ b $

b = \frac{\sqrt{1 - \rho^2}}{\rho} a

Il peut être obtenu. Si vous mettez $ a = c \ rho $ ($ c $ est une constante positive) juste pour nettoyer l'expression,

Z = c \left( \rho X + \sqrt{1 - \rho^2} Y \right)

Il peut être obtenu. Ce $ c $ est un paramètre pour contrôler la distribution de $ Z $, et en fait,

\begin{align}
\mathrm{Var}[Z] &= \mathrm{Var}\left[{c \left( \rho X + \sqrt{1 - \rho^2} Y \right)} \right]
\\
&= c^2 \left[ \rho^2 \sigma^2 + (1 - \rho) \sigma^2 \right]
\\
&= c^2 \sigma^2
\end{align}

Et a pour effet de multiplier la variance par $ c ^ 2 $ (c'est-à-dire multiplier l'écart type par $ c $). En définissant ceci sur 1 $, la formule d'ouverture

Z = \rho X + \sqrt{1 - \rho^2} Y

Peut être obtenu.

[^ 2]: La liberté est comprise entre $ -1 $ et $ 1 $.

Vérification

Vérifions maintenant que les nombres aléatoires réellement obtenus en utilisant Python 3 sont ce que nous voulons. Pour l'instant, utilisons des nombres aléatoires uniformes pour $ X $ et $ Y $.

import numpy as np
import numpy.random as rand
import matplotlib.pyplot as plt

size = int(1e3) # Size of the vector
rho_in = 0.6 # Correlation coefficient between two random variables

# Generate uniform random numbers with mean = 0
x = rand.rand(size) - 0.5
y = rand.rand(size) - 0.5

# Generate random numbers correlated with x
z = rho_in * x + (1 - rho_in ** 2) ** 0.5 * y

# Calculate stats
variance = np.var(z)
rho_out = np.corrcoef(x, z)[1][0]
print("variance: {}, correlation: {}".format(variance, rho_out))

# Plot results
fig, ax = plt.subplots()
ax.scatter(x, z, s=1)
ax.set_xlabel('X')
ax.set_ylabel('Z')
plt.show()

Le résultat ressemble à ceci.

variance: 0.08332907524577293, correlation: 0.5952102586280674

Oh? La valeur de la dispersion est à mi-chemin. Le coefficient de corrélation semble être plus proche de la cible (0,6 $). ..

Ça devrait être ça. La distribution de $ Z $ générée était la même que celle de $ X $. Maintenant $ X $ est un nombre aléatoire uniforme dans l'intervalle $ [-0.5, 0.5] $, et sa variance est [1/12](http://www.geocities.jp/jun930/etc/varianceofrandom.html "1" Puisque la distribution des nombres aléatoires est "), une valeur proche de 1/12 $ = 0,0833 ... $ sort.

Le tracé $ Z $ - $ X $ qui en résulte ressemble à ceci. Puisque le nombre aléatoire uniforme est ajouté au nombre aléatoire uniforme, il ressemble à un quadrilatère parallèle. uniform_uniform_rho0.6.png

Si vous n'aimez pas ça, allez à y

y = rand.randn(size) * (1 / 12) ** 0.5

Vous pouvez le changer en un nombre aléatoire normal comme celui-ci. L'écart type des nombres aléatoires uniformes est appliqué pour aligner les variances de $ X $ et $ Y $.

La sortie est ci-dessous. variance: 0.0823649369923887, correlation: 0.5983078612793685 uniform_normal_rho0.6.png

Cela a l'air plutôt bien, alors je l'ai utilisé.

Cependant, une mise en garde est que les nombres aléatoires $ Z $ créés par ces deux méthodes ne sont pas des nombres aléatoires uniformes. Si $ X $ et $ Y $ utilisent des nombres aléatoires uniformes, l'histogramme ressemblera à un trapèze comme celui ci-dessous. histo_uniform_uniform_rho0.6.png

C'est également le cas lorsque vous utilisez un nombre aléatoire normal pour $ Y $, ce qui donne au cado de montagne un aspect un peu arrondi.

C'est parce que l'ajout de deux variables stochastiques qui suivent des distributions uniformes différentes ne conduit pas toujours à une distribution uniforme (c'est vrai), et en termes statistiques, "[reproductibilité](https: /) /ja.wikipedia.org/wiki/%E5%86%8D%E7%94%9F%E6%80%A7 "Reproductibility-Wikipedia") n'est pas disponible. " Donc, de force encore Il y a aussi une tentative de revenir à un nombre aléatoire uniforme.

Le type de distribution qui a une reproductibilité est, par exemple, la distribution normale que tout le monde aime. Donc $ X $ et $ Y $ devraient être des nombres aléatoires réguliers! C'est la version suivante. Eh bien, la seule chose qui a changé est le «x» et le «y» supérieurs.

import numpy as np
import numpy.random as rand
import matplotlib.pyplot as plt

size = int(1e4) # Size of the vector
rho_in = 0.6 # Correlation coefficient between two random variables

# Generate normal random numbers
x = rand.randn(size)
y = rand.randn(size)

# Generate random numbers correlated with x
z = rho_in * x + (1 - rho_in ** 2) ** 0.5 * y

# Calculate stats
variance = np.var(z)
rho_out = np.corrcoef(x, z)[1][0]
print("variance: {}, correlation: {}".format(variance, rho_out))

# Plot results
fig, ax = plt.subplots()
ax.scatter(x, z, s=1)
ax.set_xlabel('X')
ax.set_ylabel('Z')
plt.show()

La sortie est comme ça. variance: 0.9884508169450733, correlation: 0.5936089719053868

J'ai essayé d'obtenir une corrélation de 0,6 $ en utilisant un nombre aléatoire normal avec une variance de 1 $, donc c'est assez proche! Bien sûr, si le coefficient de corrélation est compris entre $ [-1, 1] $, il peut être négatif ou $ 1 $ [^ 3].

Cliquez ici pour le tracé. normal_normal_rho0.6.png

Je publierai également un histogramme pour le moment. La ligne pointillée rouge est la distribution normale théoriquement prévue avec une moyenne de zéro et une variance de 1 $. Cela semble correspondre à la bonne sensation! histo_normal_normal_rho0.6.png

[^ 3]: Au fait, si vous spécifiez une valeur autre que $ [-1, 1] $ pour $ \ rho $, quelque chose d'étrange se produit. Si vous êtes intéressé, pensez ou essayez.

Résumé

Comment créer des nombres aléatoires corrélés? Et pourquoi? J'ai vu. Les gens qui n'y sont pas habitués peuvent être confus, mais je pense qu'il n'y a pas de perte à essayer de calculer la dispersion de nouveaux nombres aléatoires en raison de la nature de la dispersion (car c'est un calcul plus simple qu'il n'y paraît). Si vous souhaitez apprendre des statistiques à partir de maintenant, vous pouvez suivre les calculs de cet article.

On a également constaté qu'une simple combinaison de nombres aléatoires uniformes aboutit à une distribution étrange. J'ai utilisé $ X $: nombre aléatoire uniforme et $ Y $: version aléatoire normale parce que c'était pour une expérience légère qui ne se soucie pas beaucoup de la forme de la distribution, mais ceux qui se soucient de la distribution utilisent le nombre aléatoire normal. Cela peut être rapide après tout.

Je veux le lire ensemble

** Génération de n nombres pseudo-aléatoires corrélés (avec échantillon Python) --Qiita ** C'est la suite. Comment créer des nombres aléatoires $ n $? J'explique cela.

Recommended Posts

Génération de deux nombres pseudo-aléatoires corrélés (avec exemple Python)
Générer n nombres pseudo-aléatoires corrélés (avec l'exemple Python)
Générez des nombres de Fibonacci avec des fermetures, des itérateurs et des générateurs Python
Exemple de données créées avec python
Générer du XML (RSS) avec Python
Juger les nombres premiers avec python
Jouez des nombres manuscrits avec Python Partie 1
Tester avec des nombres aléatoires en Python
[Python] Joindre deux tables avec des pandas
[Python] Générer un mot de passe avec Slackbot
Jouez des nombres manuscrits avec python, partie 2 (identifier)
Générer des données de test japonais avec Python Faker
Exemple de conversion en ondelettes d'images en Python
Résumé de la comparaison des bibliothèques pour générer des PDF avec Python
Essayez de générer automatiquement des documents Python avec Sphinx
Remarque pour le formatage des nombres avec la fonction de format python
Exemple de notification Slack avec python lambda
Générez une instruction d'insertion à partir de CSV avec Python.
Exemple de programme qui génère un syslog avec la journalisation Python
FizzBuzz en Python3
Grattage avec Python
Statistiques avec python
Grattage avec Python
Python avec Go
Twilio avec Python
Exemple de fermeture Python
Intégrer avec Python
Jouez avec 2016-Python
AES256 avec python
Testé avec Python
python commence par ()
avec syntaxe (Python)
Bingo avec python
Zundokokiyoshi avec python
Excel avec Python
Micro-ordinateur avec Python
Cast avec python
Exemple de code spécifique pour travailler avec SQLite3 en Python
[Python] Obtenez les nombres dans l'image graphique avec OCR
[Python] Génère des nombres aléatoires qui suivent la distribution de Rayleigh