Kammregression umkehren: Finden Sie bei gegebener Antwortmatrix und gegebenen Regressionskoeffizienten geeignete Prädiktoren

Betrachten wir ein Standard OLS Regressionsproblem : Ich habe Matrizen und und ich möchte \ B finden, um L = \ | \ Y- \ X \ B \ | ^ 2 zu minimieren . Die Lösung ist gegeben durch \ hat \ B = \ argmin_ \ B \ {L \} = (\ X ^ \ top \ X) ^ + \ X ^ \ top \ Y.$\newcommand{\Y}{\mathbf Y}\newcommand{\X}{\mathbf X}\newcommand{\B}{\boldsymbol\beta}\DeclareMathOperator*{argmin}{argmin}$$\Y$$\X$$\B$

$$L=\|\Y-\X\B\|^2.$$ $$\hat\B=\argmin_\B\{L\} = (\X^\top\X)^+\X^\top \Y.$$

Ich kann auch ein "umgekehrtes" Problem aufwerfen: Wenn \ $\Y$ und $\B^*$ , finde $\hat\X$ , das $\hat\B\approx \B^*$ , dh minimiere $\|\argmin_\B\{L\}-\B^*\|^2$ . In Worten, ich habe die Antwortmatrix $\Y$ und den Koeffizientenvektor $\B^*$ und ich möchte die Prädiktormatrix finden, die Koeffizienten in der Nähe von $\B^*$ . Dies ist natürlich auch ein OLS-Regressionsproblem mit der Lösung

$$\hat\X = \argmin_\X\Big\{\|\argmin_\B\{L\}-\B^*\|^2\Big\} = \Y\B^\top(\B\B^\top)^{+}.$$

Klarstellungsaktualisierung: Wie in seiner Antwort unter @ GeoMatt22 erläutert, ist $\Y$ ein Vektor (dh, wenn es nur eine Antwortvariable gibt), dann hat $\hat \X$ den ersten Rang, und das umgekehrte Problem ist massiv unterbestimmt. In meinem Fall ist $\Y$ tatsächlich eine Matrix (dh es gibt viele Antwortvariablen, es ist eine multivariate Regression). So $\X$ ist $n\times p$ , $\Y$ ist $n\times q$ und $\B$ ist $p\times q$ .

---

Ich bin daran interessiert, ein "umgekehrtes" Problem für die Gratregression zu lösen. Nämlich, meine Verlustfunktion ist jetzt

$$L=\|\Y-\X\B\|^2+\mu\|\B\|^2$$ und die Lösung ist $$\hat\B=\argmin_\B\{L\}=(\X^\top \X+\mu\mathbf I)^{-1}\X^\top \Y.$$

Das "umgekehrte" Problem besteht darin,

$$\hat\X = \argmin_\X\Big\{\|\argmin_\B\{L\}-\B^*\|^2\Big\} = \;?$$

Wieder habe ich eine Antwortmatrix $\Y$ und einen Koeffizientenvektor $\B^*$ und ich möchte eine Prädiktormatrix finden, die Koeffizienten nahe bei \ B ^ * ergibt $\B^*$.

Tatsächlich gibt es zwei verwandte Formulierungen:

1. Find $\hat\X$ gegeben $\Y$ und $\B^*$ und $\mu$ .
2. Finden Sie $\hat\X$ und $\hat \mu$ gegeben $\Y$ und $\B^*$ .

Hat einer von ihnen eine direkte Lösung?

---

Hier ist ein kurzer Matlab-Auszug, um das Problem zu veranschaulichen:

% generate some data
n = 10; % number of samples
p = 20; % number of predictors
q = 30; % number of responses
Y = rand(n,q);
X = rand(n,p);
mu = 0;
I = eye(p);

% solve the forward problem: find beta given y,X,mu
betahat = pinv(X'*X + mu*I) * X'*Y;

% backward problem: find X given y,beta,mu
% this formula works correctly only when mu=0
Xhat =  Y*betahat'*pinv(betahat*betahat');

% verify if Xhat indeed yields betahat
betahathat = pinv(Xhat'*Xhat + mu*I)*Xhat'*Y;
max(abs(betahathat(:) - betahat(:)))

Dieser Code gibt Null aus, wenn mu=0nicht anders.

Nachdem sich die Frage auf eine genauere Formulierung des interessierenden Problems konzentriert hat, habe ich eine Lösung für Fall 1 (bekannter Kammparameter) gefunden. Dies sollte auch für Fall 2 hilfreich sein (keine exakte analytische Lösung, sondern eine einfache Formel und einige Einschränkungen).

Zusammenfassung: Keine der beiden inversen Problemformulierungen hat eine eindeutige Antwort. In Fall 2 , in dem der Kammparameter unbekannt ist, gibt es unendlich viele Lösungen für . In Fall 1, in dem angegeben ist, gibt es aufgrund der Mehrdeutigkeit im Singularwertspektrum eine endliche Anzahl von Lösungen für .$\mu\equiv\omega^2$$X_\omega$$\omega\in[0,\omega_\max]$$\omega$$X_\omega$

(Die Ableitung ist etwas langwierig, also TL, DR: Am Ende befindet sich ein funktionierender Matlab-Code.)

---

Unterbestimmter Fall ("OLS")

Das Vorwärtsproblem ist wobei , und .

$$\min_B\|XB-Y\|^2$$ $X\in\mathbb{R}^{n\times p}$$B\in\mathbb{R}^{p\times q}$$Y\in\mathbb{R}^{n\times q}$

Basierend auf der aktualisierten Frage nehmen wir an, dass , so dass bei und . Wie in der Frage wird die "Standard" (Minimum ) wobei die Pseudoinverse von .$n<p<q$$B$$X$$Y$$L_2$

$$B=X^+Y$$ $X^+$$X$

Aus der Singularwertzerlegung ( SVD ) von , gegeben durch * die Pseudoinverse berechnet werden als ** (* Die ersten Ausdrücke verwenden die vollständige SVD, während die zweiten Ausdrücke die reduzierte SVD verwenden. ** Der Einfachheit halber gehe ich davon aus, dass den vollen Rang hat, dh existiert.)$X$

$$X=USV^T=US_0V_0^T$$ $$X^+=VS^+U^T=V_0S_0^{-1}U^T$$ $X$$S_0^{-1}$

So ist die Vorwärts-Problem hat Lösung Für die Zukunft stelle ich fest , dass , wo ist der Vektor der Singularwerte.

$$B\equiv X^+Y=\left(V_0S_0^{-1}U^T\right)Y$$ $S_0=\mathrm{diag}(\sigma_0)$$\sigma_0>0$

Bei dem inversen Problem werden wir gegeben und . Wir wissen, dass aus dem obigen Prozess stammt, aber wir kennen . Die Aufgabe besteht dann darin, das entsprechende zu bestimmen .$Y$$B$$B$$X$$X$

Wie in der aktualisierten Frage erwähnt, in diesem Fall können wir erholen im Wesentlichen mit dem gleichen Ansatz, dh jetzt der Pseudo - Inverse der Verwendung .$X$

$$X_0=YB^+$$ $B$

---

Überbestimmter Fall (Ridge Estimator)

Im "OLS" -Fall wurde das unterbestimmte Problem durch die Wahl der Minimum-Norm-Lösung gelöst , dh unsere "einzigartige" Lösung wurde implizit geregelt .

Anstatt die minimale Normlösung zu wählen , führen wir hier einen Parameter zu steuern, "wie klein" die Norm sein soll, dh wir verwenden die Kammregression .$\omega$

In diesem Fall haben wir eine Reihe von Vorwärtsproblemen für , , die durch Sammeln der verschiedenen Vektoren der linken und rechten Seite in diese Sammlung von Probleme können auf das folgende "OLS" -Problem bei dem wir die erweiterten Matrizen $\beta_k$$k=1,\ldots,q$

$$\min_\beta\|X\beta-y_k\|^2+\omega^2\|\beta\|^2$$ $$B_{\omega}=[\beta_1,\ldots,\beta_k] \quad,\quad Y=[y_1,\ldots,y_k]$$ $$\min_B\|\mathsf{X}_\omega B-\mathsf{Y}\|^2$$ $$\mathsf{X}_\omega=\begin{bmatrix}X \\ \omega I\end{bmatrix} \quad , \quad \mathsf{Y}=\begin{bmatrix}Y \\ 0 \end{bmatrix}$$

In diesem überbestimmten Fall wird die Lösung immer noch durch das pseudoinverse aber das pseudoinverse wird jetzt geändert, was zu * wobei die neue Matrix "Singularitätsspektrum" eine (inverse) Diagonale ** (* Die etwas Berechnung, die erforderlich ist, um dies abzuleiten, wurde der Kürze halber weggelassen. Sie ähnelt der Darstellung hier für den Fall . ** Hier die Einträge des Vektor werden als Vektor , bei dem alle Operationen sind.)

$$B_\omega = \mathsf{X}^+\mathsf{Y}$$ $$B_\omega = \left(V_0S_\omega^{-2}U^T\right) Y$$ $$\sigma_\omega^2 = \frac{\sigma_0^2+\omega^2}{\sigma_0}$$ $p\leq n$$\sigma_\omega$$\sigma_0$

Jetzt können wir in diesem Problem noch formal eine "Basislösung" als aber dies ist keine echte Lösung mehr.

$$X_\omega=YB_\omega^+$$

Die Analogie gilt jedoch immer noch, dass diese "Lösung" SVD mit den oben angegebenen singulären Werten .

$$X_\omega=US_\omega^2V_0^T$$ $\sigma_\omega^2$

Wir können also eine quadratische Gleichung ableiten, die die gewünschten Singularwerte mit den wiederherstellbaren Singularwerten und dem Regularisierungsparameter . Die Lösung wird dann $\sigma_0$$\sigma_\omega^2$$\omega$

$$\sigma_0=\bar{\sigma} \pm \Delta\sigma \quad , \quad \bar{\sigma} = \tfrac{1}{2}\sigma_\omega^2 \quad , \quad \Delta\sigma = \sqrt{\left(\bar{\sigma}+\omega\right)\left(\bar{\sigma}-\omega\right)}$$

---

Die folgende Matlab-Demo (online über Octave getestet ) zeigt, dass diese Lösungsmethode sowohl in der Praxis als auch in der Theorie funktioniert. Die letzte Zeile zeigt, dass alle singulären Werte von in der Rekonstruktion , aber ich habe nicht vollständig herausgefunden, welche Wurzel ich nehmen soll ( = vs. ). Bei dies immer die Wurzel . Diese Regel scheint für „klein“ zu halten , während für „große“ der Wurzel scheint zu übernehmen. (Demo unten ist derzeit auf "Groß" eingestellt.)$X$$\bar{\sigma}\pm\Delta\sigma$sgn$+$$-$$\omega=0$$+$$\omega$$\omega$$-$

% Matlab demo of "Reverse Ridge Regression"
n = 3; p = 5; q = 8; w = 1*sqrt(1e+1); sgn = -1;
Y = rand(n,q); X = rand(n,p);
I = eye(p); Z = zeros(p,q);
err = @(a,b)norm(a(:)-b(:),Inf);

B = pinv([X;w*I])*[Y;Z];
Xhat0 = Y*pinv(B);
dBres0 = err( pinv([Xhat0;w*I])*[Y;Z] , B )

[Uw,Sw2,Vw0] = svd(Xhat0, 'econ');

sw2 = diag(Sw2); s0mid = sw2/2;
ds0 = sqrt(max( 0 , s0mid.^2 - w^2 ));
s0 = s0mid + sgn * ds0;
Xhat = Uw*diag(s0)*Vw0';

dBres = err( pinv([Xhat;w*I])*[Y;Z] , B )
dXerr = err( Xhat , X )
sigX = svd(X)', sigHat = [s0mid+ds0,s0mid-ds0]' % all there, but which sign?

Ich kann nicht sagen, wie robust diese Lösung ist, da inverse Probleme im Allgemeinen schlecht gestellt sind und analytische Lösungen sehr fragil sein können. Allerdings scheinen flüchtige Experimente, die mit Gaußschem Rauschen belasten (dh es hat den vollen Rang gegenüber dem reduzierten Rang ), darauf hinzudeuten, dass sich die Methode angemessen verhält.$B$$p$$n$

Für Problem 2 (dh unbekannt) gibt das Obige mindestens eine Obergrenze für . Damit die quadratische Diskriminante nicht negativ ist, muss $\omega$$\omega$

$$\omega \leq \omega_{\max} = \bar{\sigma}_n = \min[\tfrac{1}{2}\sigma_\omega^2]$$

---

Für die Mehrdeutigkeit des quadratischen Wurzelzeichens zeigt der folgende Codeausschnitt, dass unabhängig vom Vorzeichen jedes dieselbe Vorwärts- Grat-Lösung ergibt, auch wenn von abweicht .$\hat{X}$$B$$\sigma_0$$\mathrm{SVD}[X]$

Xrnd=Uw*diag(s0mid+sign(randn(n,1)).*ds0)*Vw0'; % random signs
dBrnd=err(pinv([Xrnd;w*I])*[Y;Z],B) % B is always consistent ...
dXrnd=err(Xrnd,X) % ... even when X is not