LU bug in Matrix package

There is a bug in Matrix package, please check it, thanks!

You didn't say what R code you ran to get that output or why you think it 
is wrong.

Let us experiment to see if we can guess what the problem is from the 
limited information provided

x<-t(Matrix(1:25,5))
x

5 x 5 Matrix of class "dgeMatrix"
      [,1] [,2] [,3] [,4] [,5]
[1,]    1    2    3    4    5
[2,]    6    7    8    9   10
[3,]   11   12   13   14   15
[4,]   16   17   18   19   20
[5,]   21   22   23   24   25

lux<-lu(x)

Warning message:
Exact singularity detected during LU decomposition.

Now check that the decomposition is correct

with(expand(a), P%*%L%*%U)

5 x 5 Matrix of class "dgeMatrix"
      [,1] [,2] [,3] [,4] [,5]
[1,]    1    2    3    4    5
[2,]    6    7    8    9   10
[3,]   11   12   13   14   15
[4,]   16   17   18   19   20
[5,]   21   22   23   24   25

Ok, it is correct.
Now let's Look at the matrices

expand(a)

$L
5 x 5 Matrix of class "dtrMatrix"
      [,1]       [,2]       [,3]       [,4]       [,5]
[1,] 1.00000000          .          .          .          .
[2,] 0.04761905 1.00000000          .          .          .
[3,] 0.52380952 0.50000000 1.00000000          .          .
[4,] 0.28571429 0.75000000 0.35400130 1.00000000          .
[5,] 0.76190476 0.25000000 0.50000000 0.00000000 1.00000000

$U
5 x 5 Matrix of class "dtrMatrix"
      [,1]          [,2]          [,3]          [,4]          [,5]
[1,]  2.100000e+01  2.200000e+01  2.300000e+01  2.400000e+01  2.500000e+01
[2,]             .  9.523810e-01  1.904762e+00  2.857143e+00  3.809524e+00
[3,]             .             . -1.919092e-15 -3.711332e-15 -5.614541e-15
[4,]             .             .             .  3.781500e-16  6.288326e-16
[5,]             .             .             .             .  0.000000e+00

$P
5 x 5 sparse Matrix of class "pMatrix"

[1,] . 1 . . .
[2,] . . . 1 .
[3,] . . 1 . .
[4,] . . . . 1
[5,] 1 . . . .


Perhaps the output was a trimmed version of this? The L and U matrices 
agree, but you didn't show the permutation matrix.

Did you just miss the fact that the decomposition is P%*%L%*%U? If you 
aren't used to S4 classes it is easy to miss the fact that class?LU 
gives help on the structure of the "LU" class, and if you try to guess 
what the components are without the documentation it is easy to be 
confused.


 	-thomas

Matlab result:



x =

    1     2     3     4     5

    6     7     8     9    10

   11    12    13    14    15

   16    17    18    19    20

   21    22    23    24    25

lu(x)

ans =



  21.0000   22.0000   23.0000   24.0000   25.0000

   0.0476    0.9524    1.9048    2.8571    3.8095

   0.7619    0.2500   -0.0000   -0.0000   -0.0000

   0.5238    0.5000    0.4839         0    0.0000

   0.2857    0.7500   -0.0645         0    0.0000





Gsl result:



21 22 23 24 25

0.047619 0.952381 1.90476 2.85714 3.80952

0.761905 0.25 -3.44169e-015 -6.88338e-015 -1.03251e-014

0.52381 0.5 -0.0322581 -1.77636e-015 -1.66533e-015

0.285714 0.75 -0.0645161 -0 2.22045e-016





R result:

$L

5 x 5 Matrix of class "dtrMatrix"

    [,1]       [,2]       [,3]       [,4]       [,5]

[1,] 1.00000000          .          .          .          .

[2,] 0.04761905 1.00000000          .          .          .

[3,] 0.52380952 0.50000000 1.00000000          .          .

[4,] 0.28571429 0.75000000 0.35400130 1.00000000          .

[5,] 0.76190476 0.25000000 0.50000000 0.00000000 1.00000000



$U

5 x 5 Matrix of class "dtrMatrix"

    [,1]          [,2]          [,3]          [,4]          [,5]

[1,]  2.100000e+01  2.200000e+01  2.300000e+01  2.400000e+01  2.500000e+01

[2,]             .  9.523810e-01  1.904762e+00  2.857143e+00  3.809524e+00

[3,]             .             . -1.919092e-15 -3.711332e-15 -5.614541e-15

[4,]             .             .             .  3.781500e-16  6.288326e-16

[5,]             .             .             .             .  0.000000e+00












	[[alternative HTML version deleted]]

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Thomas Lumley			Assoc. Professor, Biostatistics
tlumley at u.washington.edu	University of Washington, Seattle