[Rd] Problem in matrix definition?

[(Ted Harding)]

Does the following help to clarify things?
>From your original two different definitions of the matrix 'a':

  mp <- 10
  np <- 5
  a <- matrix(c(1:mp*np),mp,np)
  a
  #       [,1] [,2] [,3] [,4] [,5]
  #  [1,]    5    5    5    5    5
  #  [2,]   10   10   10   10   10
  #  [3,]   15   15   15   15   15
  #  [4,]   20   20   20   20   20
  #  [5,]   25   25   25   25   25
  #  [6,]   30   30   30   30   30
  #  [7,]   35   35   35   35   35
  #  [8,]   40   40   40   40   40
  #  [9,]   45   45   45   45   45
  # [10,]   50   50   50   50   50

  a <- matrix(c(1:50),mp,np)
  a
  #       [,1] [,2] [,3] [,4] [,5]
  #  [1,]    1   11   21   31   41
  #  [2,]    2   12   22   32   42
  #  [3,]    3   13   23   33   43
  #  [4,]    4   14   24   34   44
  #  [5,]    5   15   25   35   45
  #  [6,]    6   16   26   36   46
  #  [7,]    7   17   27   37   47
  #  [8,]    8   18   28   38   48
  #  [9,]    9   19   29   39   49
  # [10,]   10   20   30   40   50

Check what is said about precedence in '?Syntax'. The effect of
"c(1:mp*np)" is (1:mp)*np, i.e. first create (1:mp), and then
multiuply each element by np. So, with mp=10, np=5, you first get
c(1,2,3,4,5,6,7,8,9,10) from (1:10), and then
c(5,10,15,20,25,30,35,40,45,50) after multiplying this by 5.

This vector with 10 elements is then recycled 5 times over
when you ask for "a <- matrix(c(1:mp*np),mp,np)", i.e. the
equivalent of

  a <- matrix(c(5,10,15,20,25,30,35,40,45,50),10,5)

On the other hand, a <- matrix(c(1:50),10,5) does what you expected!

Hoping this helps,
Ted.

On 31-Aug-09 08:21:34, Uwe Ligges wrote:
> 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
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel

--------------------------------------------------------------------
E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk>
Fax-to-email: +44 (0)870 094 0861
Date: 31-Aug-09                                       Time: 10:00:04
------------------------------ XFMail ------------------------------