-----Original Message-----
From: r-devel-bounces at r-project.org
[mailto:r-devel-bounces at r-project.org] On Behalf Of Mathieu Ribatet
Sent: Monday, August 31, 2009 1:57 PM
To: Fabio Mathias Corr?a
Cc: r-devel at r-project.org; Uwe Ligges
Subject: Re: [Rd] Problem in matrix definition?
Dear Fabio,
The problem is that L[k,1:(r-1)] is not anymore a matrix but a vector.
Hence when you do t(L[k,1:(r-1)]) you get a matrix with only
one row while I think you expected one column instead. You
can see this feature with a simpler example like the following one:
x <- runif(10)
is.vector(x)
dim(t(x))##1,10
dim(t(t(x)))##10,1
Getting back to your code, you probably want that
L[k,1:(r-1)] is a matrix with r-1 columns. This is possible
by substituting
t(L[k,1:(r-1)]) for t(t(L[k,1:(r-1)])) or
matrix(L[k,1:(r-1)], nrow = 1)
Maybe somebody else will find a more elegant fix for your problem.
However, I tried and it did solve your issue.
Cheers,
Mathieu
Le lundi 31 ao?t 2009 ? 19:09 +0200, Fabio Mathias Corr?a a ?crit :
The problem is that arrays are the same size. The only
difference is how they were built!
The way to build the matrix should not influence in the outcome!
Below is the commented code in MatLab!
Neural Information Processing - Letters and Reviews Vol.8, No.2,
August 2005
function Y = geninv(G)
% Returns the Moore-Penrose inverse of the argument
% Transpose if m < n
[m,n]=size(G); transpose=false;
if m<n
transpose=true;
A=G*G';
n=m;
else
A=G'*G;
end
% Full rank Cholesky factorization of A
dA=diag(A); tol= min(dA(dA>0))*1e-9;
L=zeros(size(A));
r=0;
for k=1:n
r=r+1;
L(k:n,r)=A(k:n,k)-L(k:n,1:(r-1))*L(k,1:(r-1))';
% Note: for r=1, the substracted vector is zero
if L(k,r)>tol
L(k,r)=sqrt(L(k,r));
if k<n
L((k+1):n,r)=L((k+1):n,r)/L(k,r);
end
else
r=r-1;
end
end
L=L(:,1:r);
% Computation of the generalized inverse of G
M=inv(L'*L);
if transpose
Y=G'*L*M*M*L';
else
Y=L*M*M*L'*G';
end
Thanks!
F?bio Mathias Corr?a
Estat?stica e Experimenta??o Agropecu?ria/UFLA
Brazil