Tracer la sélection de modèles PCA vs FA

Introduction

Dans ce laboratoire, nous explorerons deux modèles probabilistes - la PCA probabiliste et l'analyse factorielle - et comparerons leur efficacité dans la sélection de modèles et l'estimation de la covariance. Nous effectuerons une validation croisée sur des données de rang faible corrompues par un bruit homoscédaastique ou hétéroscédaastique. De plus, nous comparerons la vraisemblance du modèle aux vraisemblances obtenues à partir d'estimateurs de covariance à rétrécissement.

Conseils sur la machine virtuelle

Une fois le démarrage de la machine virtuelle terminé, cliquez dans le coin supérieur gauche pour basculer vers l'onglet Carnet de notes pour accéder au carnet Jupyter pour pratiquer.

Parfois, vous devrez peut-être attendre quelques secondes pour que le carnet Jupyter ait fini de charger. La validation des opérations ne peut pas être automatisée en raison des limitations du carnet Jupyter.

Si vous rencontrez des problèmes pendant l'apprentissage, n'hésitez pas à demander à Labby. Donnez des commentaires après la session, et nous résoudrons rapidement le problème pour vous.

Créer les données

Nous allons créer un ensemble de données simulées composé de 500 échantillons, 25 caractéristiques et un rang de 5. Nous ajouterons également un bruit homoscédaastique et hétéroscédaastique à l'ensemble de données.

import numpy as np
from scipy import linalg

n_samples, n_features, rank = 500, 25, 5
sigma = 1.0
rng = np.random.RandomState(42)
U, _, _ = linalg.svd(rng.randn(n_features, n_features))
X = np.dot(rng.randn(n_samples, rank), U[:, :rank].T)

## Ajout de bruit homoscédaastique
X_homo = X + sigma * rng.randn(n_samples, n_features)

## Ajout de bruit hétéroscédaastique
sigmas = sigma * rng.rand(n_features) + sigma / 2.0
X_hetero = X + rng.randn(n_samples, n_features) * sigmas

Ajuster les modèles

Nous allons ajuster les modèles de PCA probabiliste et d'analyse factorielle à l'ensemble de données, et utiliser la validation croisée pour évaluer leur performance. Nous allons également calculer les scores pour les estimateurs de covariance à rétrécissement, et comparer les résultats.

import matplotlib.pyplot as plt
from sklearn.decomposition import PCA, FactorAnalysis
from sklearn.covariance import ShrunkCovariance, LedoitWolf
from sklearn.model_selection import cross_val_score, GridSearchCV

n_components = np.arange(0, n_features, 5)  ## options pour n_components

def compute_scores(X):
    pca = PCA(svd_solver="full")
    fa = FactorAnalysis()

    pca_scores, fa_scores = [], []
    for n in n_components:
        pca.n_components = n
        fa.n_components = n
        pca_scores.append(np.mean(cross_val_score(pca, X)))
        fa_scores.append(np.mean(cross_val_score(fa, X)))

    return pca_scores, fa_scores

def shrunk_cov_score(X):
    shrinkages = np.logspace(-2, 0, 30)
    cv = GridSearchCV(ShrunkCovariance(), {"shrinkage": shrinkages})
    return np.mean(cross_val_score(cv.fit(X).best_estimator_, X))

def lw_score(X):
    return np.mean(cross_val_score(LedoitWolf(), X))

for X, title in [(X_homo, "Bruit homoscédaastique"), (X_hetero, "Bruit hétéroscédaastique")]:
    pca_scores, fa_scores = compute_scores(X)
    n_components_pca = n_components[np.argmax(pca_scores)]
    n_components_fa = n_components[np.argmax(fa_scores)]

    pca = PCA(svd_solver="full", n_components="mle")
    pca.fit(X)
    n_components_pca_mle = pca.n_components_

    print("meilleur n_components par validation croisée PCA = %d" % n_components_pca)
    print("meilleur n_components par validation croisée FactorAnalysis = %d" % n_components_fa)
    print("meilleur n_components par estimation maximale de vraisemblance PCA = %d" % n_components_pca_mle)

    plt.figure()
    plt.plot(n_components, pca_scores, "b", label="Scores PCA")
    plt.plot(n_components, fa_scores, "r", label="Scores FA")
    plt.axvline(rank, color="g", label="VRAI: %d" % rank, linestyle="-")
    plt.axvline(
        n_components_pca,
        color="b",
        label="Validation croisée PCA: %d" % n_components_pca,
        linestyle="--",
    )
    plt.axvline(
        n_components_fa,
        color="r",
        label="Validation croisée FactorAnalysis: %d" % n_components_fa,
        linestyle="--",
    )
    plt.axvline(
        n_components_pca_mle,
        color="k",
        label="Estimation maximale de vraisemblance PCA: %d" % n_components_pca_mle,
        linestyle="--",
    )

    ## comparer avec d'autres estimateurs de covariance
    plt.axhline(
        shrunk_cov_score(X),
        color="violet",
        label="Estimation maximale de vraisemblance de covariance rétrécie",
        linestyle="-.",
    )
    plt.axhline(
        lw_score(X),
        color="orange",
        label="Estimation maximale de vraisemblance de LedoitWolf %d" % n_components_pca_mle,
        linestyle="-.",
    )

    plt.xlabel("nombre de composantes")
    plt.ylabel("Scores de validation croisée")
    plt.legend(loc="bas droite")
    plt.title(title)

plt.show()

Sommaire

Dans ce laboratoire, nous avons exploré l'efficacité des modèles de PCA probabiliste et d'analyse factorielle dans la sélection de modèles et l'estimation de la covariance. Nous avons créé un ensemble de données simulées avec un bruit homoscédaastique et hétéroscédaastique, et comparé les performances des modèles en utilisant la validation croisée. Nous avons également comparé la vraisemblance du modèle aux vraisemblances obtenues à partir d'estimateurs de covariance à rétrécissement. Les résultats ont montré que la PCA et l'AF étaient efficaces pour retrouver la taille du sous-espace de rang faible en présence de bruit homoscédaastique. Cependant, la PCA a échoué et a surestimé le rang en présence de bruit hétéroscédaastique. Dans des circonstances appropriées, les données réservées étaient plus probables pour les modèles de rang faible que pour les modèles à rétrécissement.