[PYTHON] Was ist eine reduzierte Rangkammregression?

Referenzlink

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

einfach! : bow_tone1:

"** Reduzierte Rangregression " und " Ridge-Regression **" werden kombiniert, um gleichzeitig ** erklärende variable Co-Linearität ** und ** objektivvariable niedrigdimensionale Struktur ** zu erreichen. Eine Methode der linearen Modellregression.

Was ist eine reduzierte Rangregression?

Verallgemeinertes lineares Modell

Betrachten Sie als Problemeinstellung die Vorhersage (Rückgabe) der Zielvariablen $ Y_ {N \ times Q} $ durch die erklärende Variable $ X_ {N \ times P} $.

Hier wird häufig das ** verallgemeinerte lineare Modell ** verwendet, das durch die folgende Gleichung ausgedrückt wird. $ Y = XB + E $ Wobei $ B_ {P \ times Q} $ die Koeffizientenmatrix und $ E_ {N \ times Q} $ die Zufallsfehlermatrix ist.

Die Lösung durch diese Methode der kleinsten Quadrate (OLS) ist wie folgt. $ \ hat {B} _ {OLS} = (X ^ TX) ^ {-1} X ^ TY $

Dieses Mal konzentrieren wir uns auf ** zwei Probleme **, wenn wir das verallgemeinerte lineare Modell mit der Methode der kleinsten Quadrate lösen.

Die Leistung sinkt, wenn die wahre Dimension der Zielvariablen $ Y_ {N \ times Q} $ kleiner als $ Q $ ist. (** $ Y $ Rangabfall **)
Die Leistung sinkt, wenn die Zeilen der erklärenden Variablen $ X_ {N \ times P} $ stark korreliert sind. (** $ X $ Co-Linearität **)

Dimensionsreduzierung

Somit ist eine Dimensionsreduzierung effektiv, wenn die hochdimensionalen Daten eine potentiell niedrigdimensionale Struktur aufweisen. Eine der bekanntesten ist ** Principal Analysis (PCA) **. Ein Modell, das nach Reduzierung der Dimensionen durch Hauptkomponentenanalyse eine Regression durchführt, wird als ** Hauptkomponentenregression (PCR) ** bezeichnet.

Referenz: Grundlegende Komponentenanalyse verstehen, [Grundlegende Komponentenanalyse mit Python verstehen](https://qiita.com/maskot1977/items/ 082557fcda78c4cdb41f)

Ein Modell, das die Zielvariable durch Auswahl eines niedrigdimensionalen Faktors (Faktors) aus den erklärenden Variablen, wie z. B. der Hauptkomponentenregression, rekursiv macht, wird als ** lineares Faktormodell ** bezeichnet. Zum Beispiel unabhängige Komponentenregression, partielle Regression kleinster Quadrate, kanonische Korrelationsanalyse usw.

Eine andere Methode zur Dimensionsreduzierung ist ** Reduzierte Rangregression **. Eine reduzierte Rangregression ermöglicht es der Regression, potenzielle niedrigdimensionale Strukturen anzunehmen, indem die Fehlerfunktion minimiert und gleichzeitig der Rang der Koeffizientenmatrix $ B $ begrenzt wird.

Referenz: Reduzierte Rangregression

Regulierung

Für das zweite Problem wird häufig die Kongruenz ** der erklärenden Variablen ** $ X $ **, ** Regularisierung ** durchgeführt. Bekannte Techniken umfassen ** Ridge Regression ** und ** LASSO Regression **. Die LASSO-Regression wird hauptsächlich zur Schätzung der Spärlichkeit verwendet und auch als Methode zur Auswahl von Merkmalen. Die Ridge-Regression wird häufig verwendet, um das Problem der schlechten Einstellungen aufgrund der Co-Linearität der erklärenden Variablen anzugehen.

Erweiterung auf Reduced Rank Ridge Regression

Führen Sie zwei Einschränkungen in die Fehlerfunktion ein

Basierend auf dem Obigen werden wir, um die beiden Probleme zu lösen, die folgenden zwei Einschränkungen zum quadratischen Fehler hinzufügen.

Regularisierungsterm der Gratregression
Rangbeschränkung von $ \ hat {B} $

Daher können Sie $ \ hat {B} (\ lambda, r) $ erhalten, indem Sie den Fehler wie folgt minimieren.

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

$ R \ leq \ min \ lbrace P, Q \ rbrace $ und $ \ Vert \ cdot \ Vert_F ^ 2 $ ist jedoch die Frobenius-Norm (Matrixnorm).

Dies im Rahmen einer reduzierten Rangregression zu betrachten

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)

Auf diese Weise wird die Fehlerfunktion wie folgt ausgedrückt.

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

Verwenden Sie außerdem den geschätzten Wert der Gratregression $ \ hat {Y_R} ^ * = X ^ * \ hat {B_R} ^ * $ aus der normalen Orthogonalität.

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

Da der erste Term nicht von $ B $ abhängt, kann die Fehlerfunktion wie folgt ausgedrückt werden:

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

Ableitung einer Lösung durch ein lineares Problem

Nehmen wir nun an, dass die Singularitätszerlegung wie folgt gegeben ist:

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

Wobei $ \ sigma_i $ der Singularwert ist, $ u_i, v_i $ der linke und rechte Singularwertvektor und $ \ tau $ der Rang von $ \ hat {Y} _R ^ * $ ($ Q $ durch Regularisierung) Es wird).

Dies ist Eckert Youngs Theorem in Frobenius Norm (Referenz Stellen Sie sich das als% 86% E6% 9E% 90)) vor, und die optimale Annäherung an den Rang $ r $ lautet:

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

Verwenden wir dies, um die optimale Koeffizientenmatrix $ \ hat {B} (\ lambda, r) $ darzustellen.

\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)

Auf diese Weise erhalten wir $ \ hat {B} (\ lambda, r) = \ hat {B} _R ^ * P_r $ als Lösung.

Darüber hinaus kann es unter Verwendung der Lösung der Kammregression im kleinsten Quadrat wie folgt ausgedrückt werden.

\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

Hier ist $ \ hat {Y} _ {\ lambda} $ die geschätzte Lösung der Gratregression in $ \ lambda $ und "** die optimale Lösung durch Projizieren der geschätzten Lösung der Gratregression in den Raum der $ r $ -Dimension. Kann als ** "interpretiert werden. Daher ist für $ r = Q $ eine einfache Ridge-Regressionslösung die optimale Lösung.

Parametereinstellung

Durch Einstellen der Parameter $ \ lambda, r $ wird die Leistung der reduzierten Rangkammregression bestimmt. In dem Referenzpapier wurde eine K-fache Kreuzvalidierung durchgeführt, um die optimalen Parameter zu bestimmen.

Probieren Sie das Paket aus

Führen Sie Demo im RRRR-Paket des Referenzlinks aus. Sah.

Daten

Verwenden wir die bereitgestellten Demodaten.

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()

Aus der Korrelationsmatrix können wir erkennen, dass die Zielvariable $ Y $ eine niedrigdimensionale Struktur hat.

Gegenüberstellung

Es wird eine Hyperparameter-Kreuzvalidierung durchgeführt und eine Optimierung für den Einschränkungsrang und die Stärke der Regularisierung durchgeführt.

 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

Lassen Sie uns das Ergebnis der Kreuzungsüberprüfung visualisieren. Sie können sehen, dass das Ranking der Punktzahl im Teil des optimalen Rankings hoch ist (um 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()

Inferenz

Schlussfolgerungen basierend auf den optimalen Hyperparametern, die durch Kreuzungsüberprüfung erhalten wurden.

 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

Die MSE (mittlerer quadratischer Fehler) ist kleiner und $ \ hat {Y} $ scheint gut vorhergesagt zu sein.

Zusammenfassung

Eine reduzierte Rangkammregression hat eine hohe Fähigkeit, niedrigdimensionale Zielvariablen anzupassen, und kann in großem Umfang auf reale Datenstrukturen angewendet werden. Außerdem ist das Paket gemäß scicit-learn implementiert, sodass es einfach zu bedienen zu sein scheint.