[PYTHON] Qu'est-ce que la régression de crête de rang réduit?

Lien de référence

Mukherjee A, Zhu J. Reduced Rank Ridge Regression and Its Kernel Extensions. Stat Anal Data Min. 2011;4(6):612–622. doi:10.1002/sam.10138
RRRR - Python implementation of Reduced Rank Ridge Regression

simplement! : bow_tone1:

«** Régression de rang réduit » et « Régression de crête **» sont combinés pour correspondre à la ** co-linéarité de la variable explicative ** et à la ** structure à faible dimension de la variable objective ** en même temps. Une méthode de régression de modèle linéaire.

Qu'est-ce que la régression de rang réduit?

Modèle linéaire généralisé

Comme paramètre de problème, envisagez de prédire (renvoyer) la variable objectif $ Y_ {N \ fois Q} $ par la variable explicative $ X_ {N \ fois P} $.

On utilise souvent ici le ** modèle linéaire généralisé ** exprimé par l'équation suivante. $ Y = XB + E $ Où $ B_ {P \ fois Q} $ est la matrice des coefficients et $ E_ {N \ fois Q} $ est la matrice des erreurs aléatoires.

La solution par cette méthode des moindres carrés (OLS) est la suivante. $ \ hat {B} _ {OLS} = (X ^ TX) ^ {-1} X ^ TY $

Cette fois, nous nous concentrerons sur ** deux problèmes ** lors de la résolution du modèle linéaire généralisé par la méthode des moindres carrés.

Les performances diminuent lorsque la vraie dimension de la variable objectif $ Y_ {N \ fois Q} $ est inférieure à $ Q $. (** Baisse de rang $ Y $ **)
Les performances diminuent lorsque les lignes de la variable explicative $ X_ {N \ fois P} $ sont fortement corrélées. (** $ X $ colinéarité **)

Réduction de dimension

Ainsi, la réduction de dimension est efficace lorsqu'il existe une structure potentielle de faible dimension dans les données de forte dimension. L'une des plus connues est ** l'analyse principale (ACP) **. Un modèle qui effectue une régression après avoir réduit les dimensions par analyse en composantes principales est appelé ** régression en composantes principales (PCR) **.

Référence: Comprendre l'analyse des composants principaux, [Comprendre l'analyse des composants principaux avec Python](https://qiita.com/maskot1977/items/ 082557fcda78c4cdb41f)

Un modèle qui récurent la variable objective en sélectionnant un facteur de faible dimension (facteur) parmi les variables explicatives, comme la régression en composantes principales, est appelé un ** modèle à facteurs linéaires **. Par exemple, régression en composantes indépendantes, régression partielle des moindres carrés, analyse de corrélation canonique, etc.

Une autre méthode de réduction de dimension est la ** régression de rang réduit **. La régression de rang réduit permet à la régression de supposer des structures potentielles de faible dimension en minimisant la fonction d'erreur tout en limitant le rang de la matrice de coefficients $ B $.

Référence: Régression de rang réduit

Régularisation

Pour le deuxième problème, la congruence ** de la variable explicative ** $ X $ **, ** régularisation ** est souvent effectuée. Les techniques bien connues incluent la ** régression de crête ** et la ** régression de LASSO **. La régression LASSO est principalement utilisée pour estimer la parcimonie et est également utilisée comme méthode de sélection des entités. La régression Ridge est largement utilisée pour résoudre le problème des mauvais réglages dus à la colinéarité des variables explicatives.

Extension à la régression de crête de rang réduit

Introduisez deux contraintes dans la fonction d'erreur

Sur la base de ce qui précède, afin de traiter les deux problèmes, nous ajouterons les deux contraintes suivantes à l'erreur quadratique.

Terme de régularisation de la régression des crêtes
Contrainte de rang de $ \ hat {B} $

Par conséquent, vous pouvez obtenir $ \ hat {B} (\ lambda, r) $ en minimisant l'erreur comme suit.

\underset{ \lbrace B:rank(B) \leq r \rbrace}{argmin} \Vert Y-XB \Vert_F^2 + \lambda \Vert B \Vert_F^2

Cependant, $ r \ leq \ min \ lbrace P, Q \ rbrace $ et $ \ Vert \ cdot \ Vert_F ^ 2 $ est la norme de Frobenius (norme de la matrice).

Penser cela dans le cadre de la régression de rang réduit

X_{(N+P)\times P}^* = 
\left(  \begin{array}{c}    X \\    \sqrt{\lambda}I  \end{array}\right), \ 
Y_{(N+P)\times Q}^* = 
\left(  \begin{array}{c}    Y \\    0  \end{array}\right)

Ce faisant, la fonction d'erreur est exprimée comme suit.

\underset{ \lbrace B:rank(B) \leq r \rbrace}{argmin} \Vert Y^*-X^*B \Vert_F^2

De plus, en utilisant la valeur estimée de la régression de crête $ \ hat {Y_R} ^ * = X ^ * \ hat {B_R} ^ * $, à partir de l'orthogonalité normale,

\Vert Y^*-X^*B \Vert_F^2 = \Vert Y^*-\hat{Y}_R^* \Vert_F^2 + \Vert \hat{Y}_R^*-X^*B \Vert_F^2

Puisque le premier terme ne dépend pas de $ B $, la fonction d'erreur peut être exprimée comme suit:

\underset{ \lbrace B:rank(B) \leq r \rbrace}{argmin} \Vert \hat{Y}_R^*-X^*B \Vert_F^2

Dérivation d'une solution par un problème linéaire

Supposons maintenant que la décomposition de singularité soit donnée comme suit:

\hat{Y}_R^* = \sum_{i=1}^{\tau} \sigma_i u_i v_i^T

Ici, $ \ sigma_i $ est une valeur singulière, $ u_i, v_i $ est un vecteur de valeur singulière gauche et droite, et $ \ tau $ est le rang de $ \ hat {Y} _R ^ * $ Il devient).

Il s'agit du théorème d'Eckert Young dans la norme Frobenius ([Reference](https://ja.wikipedia.org/wiki/%E4%B8%BB%E6%88%90%E5%88%86%E5%88] Considérez-le comme% 86% E6% 9E% 90)), et l'approximation optimale du rang $ r $ est:

\hat{Y}_r^* = \sum_{i=1}^{r} \sigma_i u_i v_i^T

Utilisons ceci pour représenter la matrice de coefficients optimale $ \ hat {B} (\ lambda, r) $.

\hat{Y}_r^* = \sum_{i=1}^{r} \sigma_i u_i v_i^T = \left(\sum_{i=1}^{\tau} \sigma_i u_i v_i^T\right) \left(\sum_{j=1}^{r} u_j v_j^T\right) = \hat{Y}_r^* P_r = X^* \hat{B}_R^* P_r = X^* \hat{B}(\lambda,r)

Alors, la solution est $ \ hat {B} (\ lambda, r) = \ hat {B} _R ^ * P_r $.

De plus, elle peut être exprimée comme suit en utilisant la solution des moindres carrés de régression des crêtes.

\hat{B}(\lambda,r) = \hat{B}_R^* P_r = (X^T X + \lambda I)^{-1}X^T Y P_r \\
\hat{Y}(\lambda,r) = X \hat{B}(\lambda,r) = X (X^T X + \lambda I)^{-1}X^T Y P_r = \hat{Y}_{\lambda} P_r

Ici, $ \ hat {Y} _ {\ lambda} $ est la solution estimée de la régression de crête dans $ \ lambda $, et "** La solution optimale en projetant la solution estimée de la régression de crête dans l'espace de la dimension $ r $. Cela peut être interprété comme «obtenir **». Par conséquent, pour $ r = Q $, une solution de régression de crête simple est la solution optimale.

Réglage des paramètres

Le réglage des paramètres $ \ lambda, r $ déterminera les performances de la régression de crête de rang réduit. Dans l'article de référence, une validation croisée du pli K a été effectuée pour déterminer les paramètres optimaux.

Essayez le package

Exécutez Demo dans RRRR Package du lien de référence. Vu.

Les données

Utilisons les données de démonstration fournies.

import numpy as np
from sklearn.model_selection import GridSearchCV
from sklearn.model_selection import PredefinedSplit
import reduced_rank_regressor as RRR

import matplotlib.pyplot as plt
%matplotlib inline

N_PARAMETERS_GRID_SEARCH = 20 
Data_path = "data/"

 Load Data
trainX = np.loadtxt(Data_path+"trainX.txt")
testX = np.loadtxt(Data_path+"testX.txt")
validX = np.loadtxt(Data_path+"validX.txt")

trainY = np.loadtxt(Data_path+"trainY.txt")
testY = np.loadtxt(Data_path+"testY.txt")
validY = np.loadtxt(Data_path+"validY.txt")

 Inspection of Data
f,ax = plt.subplots(1,2,figsize=(10,10))
ax[0].imshow(np.corrcoef(trainX), cmap='jet')
ax[0].set_title('trainX')
ax[1].imshow(np.corrcoef(trainY), cmap='jet')
ax[1].set_title('trainY')
plt.show()

A partir de la matrice de corrélation, nous pouvons voir que la variable objective $ Y $ a une structure de faible dimension.

Vérification croisée

Une validation croisée des hyperparamètres est effectuée et une optimisation est effectuée pour le rang de contrainte et la force de la régularisation.

 Cross-validation setup. Define search spaces
rank_grid = np.linspace(1,min(min(trainX.shape),min(trainY.shape)), num=N_PARAMETERS_GRID_SEARCH)
rank_grid = rank_grid.astype(int)
reg_grid = np.power(10,np.linspace(-20,20, num=N_PARAMETERS_GRID_SEARCH+1))
parameters_grid_search  = {'reg':reg_grid, 'rank':rank_grid}

valid_test_fold = np.concatenate((np.zeros((trainX.shape[0],))-1,np.zeros((validX.shape[0],))))
ps_for_valid = PredefinedSplit(test_fold=valid_test_fold)

 Model initialisation
rrr = RRR.ReducedRankRegressor()#rank, reg)
grid_search = GridSearchCV(rrr, parameters_grid_search, cv=ps_for_valid,
                           scoring='neg_mean_squared_error')

print("fitting...")
grid_search.fit(np.concatenate((trainX,validX)), np.concatenate((trainY,validY)))

 Display the best combination of values found
print(grid_search.best_params_)
means = grid_search.cv_results_['mean_test_score']
means = np.array(means).reshape(N_PARAMETERS_GRID_SEARCH, N_PARAMETERS_GRID_SEARCH+1)
print(grid_search.best_score_)

[output]

fitting...
{'rank': 316, 'reg': 1e-20}
-75126.47521541138

Visualisons le résultat de la vérification d'intersection. Vous pouvez voir que le classement du score est élevé dans la partie du classement optimal (environ 316).

 Show CV results
scores = [x for x in grid_search.cv_results_['rank_test_score']]
scores = np.array(scores).reshape(N_PARAMETERS_GRID_SEARCH, N_PARAMETERS_GRID_SEARCH+1)

f,ax = plt.subplots(1,1,figsize=(10,10))
cbar = ax.imshow(scores, cmap='jet')
ax.set_title('Test score Ranking')
ax.set_xlabel('Regression')
ax.set_ylabel('Rank')
ax.set_xticks(list(range(len(reg_grid))))
ax.set_yticks(list(range(len(rank_grid))))
ax.set_xticklabels([str("%0.*e"%(0,x)) for x in reg_grid],rotation=30)
ax.set_yticklabels([str(x) for x in rank_grid])
f.colorbar(cbar)
plt.show()

inférence

Diminution basée sur les hyperparamètres optimaux obtenus par vérification d'intersection.

 Train a model with the best set of hyper-parameters found
rrr.rank                = int(grid_search.best_params_['rank'])
rrr.reg                 = grid_search.best_params_['reg']
rrr.fit(trainX, trainY)


 Testing 
Yhat = rrr.predict(testX).real
MSE = (np.power((testY - Yhat),2)/np.prod(testY.shape)).mean()
print("MSE =",MSE)

f,ax = plt.subplots(1,2,figsize=(10,10))
ax[0].imshow(np.corrcoef(testY), cmap='jet')
ax[0].set_title('testY')
ax[1].imshow(np.corrcoef(Yhat), cmap='jet')
ax[1].set_title('Yhat')
plt.show()

[output]

MSE = 0.1508996355795968

La MSE (erreur quadratique moyenne) est plus petite et $ \ hat {Y} $ semble bien prédite.

Résumé

La régression par crête de rang réduit a une grande capacité à s'adapter à des variables objectives de faible dimension et a le potentiel d'être largement applicable aux structures de données du monde réel. De plus, le package est implémenté selon scicit-learn, il semble donc facile à utiliser.