david.clayton@cimr.cam.ac.uk
2005-Oct-27 14:50 UTC
[Rd] Column names in qr() and chol() (PR#8258)
I am using 2.2.0 If the QR decomposition of an N*M matrix is such that the pivoting order is not 1:M, Q%*%R does not result in the original matrix but in a matrix with the columns permuted. This is clearly intentional, and probably to be expected if pivoting is used --- chol() behaves in the same manner (it would perhaps be nice if the qr help page made that clear in the same way that the chol() help page does). The small bug is that column names are not permuted in the same way, so that, after permuting the columns of Q%*%R so that it is equal to the input matrix, it is mislabelled. The same small bug is present in chol() when pivoting is switched on. It's a minor matter but, perhaps, easily fixed. I think one line which reordered the column names of the $qr list element of the object returned by qr() would fix things (and similarly for chol) Example:> A <- 1:10 > B <- A^2 > C <- B-A > D <- rep(1, 10) > X <- cbind(A,B,C,D) > qrX <- qr(X) > qrX$pivot[1] 1 2 4 3> oo <- order(qrX$pivot) > Q <- qr.Q(qrX) > R <- qr.R(qrX) > Q%*%RA B C D [1,] 1 1 1 2.875732e-14 [2,] 2 4 1 2.000000e+00 [3,] 3 9 1 6.000000e+00 [4,] 4 16 1 1.200000e+01 [5,] 5 25 1 2.000000e+01 [6,] 6 36 1 3.000000e+01 [7,] 7 49 1 4.200000e+01 [8,] 8 64 1 5.600000e+01 [9,] 9 81 1 7.200000e+01 [10,] 10 100 1 9.000000e+01> (Q%*%R)[,oo]A B D C [1,] 1 1 2.875732e-14 1 [2,] 2 4 2.000000e+00 1 [3,] 3 9 6.000000e+00 1 [4,] 4 16 1.200000e+01 1 [5,] 5 25 2.000000e+01 1 [6,] 6 36 3.000000e+01 1 [7,] 7 49 4.200000e+01 1 [8,] 8 64 5.600000e+01 1 [9,] 9 81 7.200000e+01 1 [10,] 10 100 9.000000e+01 1 # Note: column names of last two columns have been interchanged> XtX <- t(X)%*%X > U <- chol(XtX, pivot=T) > oo <- order(attr(U, "pivot")) > XtXA B C D A 385 3025 2640 55 B 3025 25333 22308 385 C 2640 22308 19668 330 D 55 385 330 10> t(U[,oo])%*% U[,oo]D A B C D 385 3025 2640 55 A 3025 25333 22308 385 B 2640 22308 19668 330 C 55 385 330 10 # Matrix has been successfully reproduced but mislabelled David Clayton