Décorrélation de variables aléatoires - correction¶

On construit des variables corrélées gaussiennes et on cherche à construire des variables décorrélées en utilisant le calcul matriciel. (correction)

Création d’un jeu de données¶

Q1¶

La première étape consiste à construire des variables aléatoires normales corrélées dans une matrice $N \times 3$ . On cherche à construire cette matrice au format numpy. Le programme suivant est un moyen de construire un tel ensemble à l’aide de combinaisons linéaires. Complétez les lignes contenant des .....

[1]:

import random
import numpy as np


def combinaison():
    x = random.gauss(0, 1)  # génère un nombre aléatoire
    y = random.gauss(0, 1)  # selon une loi normale
    z = random.gauss(0, 1)  # de moyenne null et de variance 1
    x2 = x
    y2 = 3 * x + y
    z2 = -2 * x + y + 0.2 * z
    return [x2, y2, z2]


li = [combinaison() for i in range(0, 100)]
mat = np.matrix(li)
mat[:5]

[1]:

matrix([[ 0.17179957, -0.80900898, -1.53083278],
        [-0.48698277, -0.95015277,  1.66842305],
        [ 0.53400158,  0.42468545, -2.41318128],
        [ 2.97269198,  9.62935172, -5.0672538 ],
        [ 0.88952227,  2.09758362, -2.36428305]])

Q2¶

[2]:

npm = mat
t = npm.transpose()
a = t @ npm
a /= npm.shape[0]
a

[2]:

matrix([[ 1.00968737,  3.07225686, -2.0111879 ],
        [ 3.07225686, 10.26264038, -5.16277675],
        [-2.0111879 , -5.16277675,  5.05453658]])

a est la matrice de covariance.

Corrélation de matrices¶

Q3¶

[3]:

cov = a

[4]:

var = np.array([cov[i, i] ** (-0.5) for i in range(cov.shape[0])])
var.resize((3, 1))
varvar = var @ var.transpose()
varvar

[4]:

array([[0.99040557, 0.31065402, 0.4426555 ],
       [0.31065402, 0.09744081, 0.13884485],
       [0.4426555 , 0.13884485, 0.19784207]])

[5]:

cor = np.multiply(cov, varvar)
cor

[5]:

matrix([[ 1.        ,  0.95440895, -0.89026339],
        [ 0.95440895,  1.        , -0.71682496],
        [-0.89026339, -0.71682496,  1.        ]])

Q4¶

[6]:

def correlation(npm):
    t = npm.transpose()
    a = t @ npm
    a /= npm.shape[0]
    var = np.array([cov[i, i] ** (-0.5) for i in range(cov.shape[0])])
    var.resize((3, 1))
    varvar = var @ var.transpose()
    return np.multiply(cov, varvar)


correlation(npm)

[6]:

matrix([[ 1.        ,  0.95440895, -0.89026339],
        [ 0.95440895,  1.        , -0.71682496],
        [-0.89026339, -0.71682496,  1.        ]])

Calcul de la racine carrée¶

Q6¶

Le module numpy propose une fonction qui retourne la matrice $P$ et le vecteur des valeurs propres $L$ :

L,P = np.linalg.eig(a)

Vérifier que $P'P=I$ . Est-ce rigoureusement égal à la matrice identité ?

[7]:

L, P = np.linalg.eig(a)

[8]:

P.transpose() @ P

[8]:

matrix([[ 1.00000000e+00,  2.91433544e-16, -5.55111512e-17],
        [ 2.91433544e-16,  1.00000000e+00, -2.77555756e-17],
        [-5.55111512e-17, -2.77555756e-17,  1.00000000e+00]])

C’est presque l’identité aux erreurs de calcul près.

Q7¶

np.diag(l) construit une matrice diagonale à partir d’un vecteur.

[9]:

np.diag(L)

[9]:

array([[1.44438859e+01, 0.00000000e+00, 0.00000000e+00],
       [0.00000000e+00, 1.60670907e-03, 0.00000000e+00],
       [0.00000000e+00, 0.00000000e+00, 1.88137173e+00]])

Q8¶

Ecrire une fonction qui calcule la racine carrée de la matrice $\frac{1}{n}M'M$ (on rappelle que $M$ est la matrice npm). Voir aussi Racine carrée d’une matrice.

[10]:

def square_root_matrix(M):
    L, P = np.linalg.eig(M)
    L = L**0.5
    root = P @ np.diag(L) @ P.transpose()
    return root


root = square_root_matrix(cov)
root

[10]:

matrix([[ 0.30467222,  0.77731521, -0.55914512],
        [ 0.77731521,  2.95469126, -0.96344225],
        [-0.55914512, -0.96344225,  1.9528626 ]])

On vérifie qu’on ne s’est pas trompé.

[11]:

root @ root - cov

[11]:

matrix([[ 1.77635684e-15,  1.77635684e-15, -1.33226763e-15],
        [ 1.33226763e-15, -1.06581410e-14,  4.44089210e-15],
        [-1.33226763e-15,  4.44089210e-15, -8.88178420e-16]])

Décorrélation¶

[12]:

np.linalg.inv(cov)

[12]:

matrix([[ 577.59780655, -117.8510365 ,  109.45002008],
        [-117.8510365 ,   24.24634204,  -22.12707857],
        [ 109.45002008,  -22.12707857,   21.14682284]])

Q9¶

Chaque ligne de la matrice $M$ représente un vecteur de trois variables corrélées. La matrice de covariance est $V=\frac{1}{n}M'M$ . Calculer la matrice de covariance de la matrice $N=M V^{-\frac{1}{2}}$ (mathématiquement).

Q10¶

Vérifier numériquement.

Simulation de variables corrélées¶

Q11¶

A partir du résultat précédent, proposer une méthode pour simuler un vecteur de variables corrélées selon une matrice de covariance $V$ à partir d’un vecteur de lois normales indépendantes.

Q12¶

Proposer une fonction qui crée cet échantillon :

[13]:

def simultation(N, cov):
    # simule un échantillon de variables corrélées
    # N : nombre de variables
    # cov : matrice de covariance
    # ...
    return M

Q13¶

Vérifier que votre échantillon a une matrice de corrélations proche de celle choisie pour simuler l’échantillon.

Notebook on github