I'm trying to obtain an estimator $f(x)=y$ where $x \in \mathbb{R}^{D_1}$ and $y \in \mathbb{R}^{D_2}$, both are column vectors. So my training set $X$ and $Y$ are data matrices of size $D_1 \times N$ and $D_2 \times N$, respectively, where $N$ is the number of samples, and $D$'s are the input (feature) and output dimensions.
So I want to learn $\beta$ that gives $\beta x \sim y$ in a least-squares fashion. I was doing this in MATLAB simply by beta_hat = Y * pinv(X); and it seems like working without a problem. Though I want to ask, is this correct?
My question:
Now I want to implement this without pinv because I want to add regularization to it, so I came up with this solution (this is without regularization) :
$\hat \beta = Y (X^TX)^{-1}X^T$ is this correct? It also works but MATLAB complains about this :
Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 2.565271e-20.
And even crashes sometimes. So I think I'm making a mistake somewhere, but where?
Thanks in advance,
Edit
Here is what my MATLAB code looks like : (I know there are non-initialized variables like N, but I just cropped them out, they are working as expected) :
Ntr = round(N * 0.7); % Assign first 70% of the samples as training set
Trains = [1:Ntr]; Tests = [Ntr+1:N];
XData = zeros(FeatureSize, N);
YData = zeros(OutputSize, N);
for n=1:N
% Collect the independent data (into the columns of X)
XData(:,n) = getFeature(sample(n));
% Collect output variable for Train samples :
if find(Trains==n)
YData(:,n) = getLabel(sample(n));
end
end % for each sample
% Learn model:
if strcmp(RegressionType, 'ordinary')
C = YData(:,Trains) * pinv(XData(:,Trains));
elseif strcmp(RegressionType, 'ordinary_myImplementation')
X = XData(:,Trains);
Y = YData(:,Trains);
C = Y * inv(X'*X)*X'; % this is where the error happens. Isn't this the same with pinv(X) ?
elseif strcmp(RegressionType, 'ridge')
X = XData(:,Trains);
Y = YData(:,Trains);
C = Y * inv(X'*X + alpha*eye(Ntr,Ntr)) * X';
else, error('Unknown regression type'); end
% Apply model on Test samples :
YData(:,Tests) = C * XData(:,Tests);
Edit 2
After @Matthew Drury's suggestion, I replaced the line C = Y * inv(X'*X)*X'; to C = linsolve(X',Y')';
But now I'm getting this error:
Warning: Rank deficient, rank = 17, tol = 2.729816e-12.
Is this normal?
pinvif possible and use themldivideoperator. You do not show the way you try to solve this system. Especially if you use ridge regression, ie. you solve $(X^TX + \lambda I)^{-1} (X^T Y)$ it is even harder to get a singular matrix (you amp all the eigenvalues by $\lambda$ essentially. Please give some reproduction code so we can see your estimation procedure. I strongly suspect that the issue stems from there. – usεr11852 Oct 16 '15 at 17:49solve(A, y)which returns the solution(s) to the linear equation $Ax = y$. Always favor this. Here's a nice blog on the subject: http://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/ – Matthew Drury Oct 16 '15 at 18:38linsolve: http://www.mathworks.com/help/matlab/ref/linsolve.html – Matthew Drury Oct 16 '15 at 18:44X\Y, you have overcomplicated yourself. You do not tell us what your $N$ or $D1$, $D2$ so we cannot say if you might have sample sizes. If you have $D1$ or $D2$ larger than 17 this will be rank deficient (based on the warning message). – usεr11852 Oct 16 '15 at 19:20Ax = B, whereas I want to solve forxA = B. – jeff Oct 16 '15 at 19:41