From ligges at statistik.tu-dortmund.de  Mon Aug 31 10:21:34 2009
From: ligges at statistik.tu-dortmund.de (Uwe Ligges)
Date: Mon, 31 Aug 2009 10:21:34 +0200
Subject: [Rd] Problem in matrix definition?
In-Reply-To: <895394.79925.qm@web59711.mail.ac4.yahoo.com>
References: <895394.79925.qm@web59711.mail.ac4.yahoo.com>
Message-ID: <4A9B880E.9000405@statistik.tu-dortmund.de>

Try to learn how to debug.
The following copied from my R session might give you some hints:



 > options(error=recover)
 > geninv(a)
Error in L[k:n, 1:(r - 1)] %*% (t(L[k, 1:(r - 1)])) :
   non-conformable arguments

Enter a frame number, or 0 to exit

1: geninv(a)

Selection: 1
Called from: eval(expr, envir, enclos)
Browse[1]> (t(L[k, 1:(r - 1)]))
          [,1]   [,2]
[1,] 75.68261 29.277
Browse[1]> L[k:n, 1:(r - 1)]
           [,1]    [,2]
[1,]  75.68261 29.2770
[2,] 103.71320 43.9155
[3,] 131.74380 58.5540



Best,
Uwe Ligges




Fabio Mathias Corr?a wrote:
> I'm implementing a function to compute the moore-penrose inverse, using a code from the article: Fast Computation of Moore-Penrose Inverse Matrices. Neural Information Processing - Letters and Reviews. Vol.8, No.2, August 2005
> 
> However, the R presents an error message when I use the geninv.
> 
> The odd thing is that the error occurs for some arrays, however they have the same size. And the R indicates the lack of compatibility between the matrix!
> 
> Below is an example:
> 	
> Creating the function geninv
> 
> geninv <- function(x)
> {
>  m <- dim(x)[1]
>  n <- dim(x)[2]
>  tr <- 0
>  if(m < n) {
>    a <- tcrossprod(x)
>    n <- m
>    tr <- 1
>  }
>  else  a <- crossprod(x)
>  dA <- diag(a)
>  tol=min(dA[dA>0])*1e-9
>  L = a*0
>  r = 0
>  for(k in 1:n){
>   r = r+1
>   L[k:n,r] = a[k:n,k]-(L[k:n,1:(r-1)]%*%(t(L[k,1:(r-1)])))
>   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]
>   }
>   else r <- r-1
>  }
>  L <- L[,1:r]
>  M <- solve(crossprod(L))
>  if (tr == 1) Y <- t(x)%*%L%*%M%*%M%*%t(L)
>  else Y <- L%*%M%*%M%*%t(L)%*%t(x)
>  return(Y)
>  }
>  
> # Perfect result! This result is identical of the function ginv!
> 
> library(MASS)
> mp <- 10
> np <- 5
> a <- matrix(c(1:mp*np),mp,np)
> dim(a) # 10,5
> geninv(a)
> ginv(a)
> 
> # Problem
> a <- matrix(c(1:50),mp,np) # The difference is the vector (1:50)
> dim(a) # 10,5
> 
> geninv(a)
> Error in L[k:n, 1:(r - 1)] %*% (t(L[k, 1:(r - 1)])) : 
>   arguments are not compatible
> 
> 
> The problem this in matrix definition?
> 
> Thanks very much!
> 
>                F?bio Mathias Corr?a
> Estat?stica e Experimenta??o Agropecu?ria/UFLA
> 
> 
> 
>       ____________________________________________________________________________________
> Veja quais s?o os assuntos do momento no Yahoo! +Buscados
> http://br.maisbuscados.yahoo.com
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel



[Rd] Problem in matrix definition?
Uwe Ligges