[R] qr() and Gram-Schmidt
Peter Dalgaard
P.Dalgaard at biostat.ku.dk
Mon Nov 3 11:38:46 CET 2008
cruz wrote:
> Hi,
>
> Why the qr() produces a negative Q compared with Gram-Schmidt? (note
> example below, except Q[2,3])
This is a recurrent question in various guises (related to sign issues
in factor analysis and PCA). Probably the easiest answer is "Why not?".
Notice that at each step of Gram-Schmidt, the sign is really arbitrary:
It gives you _an_ orthonormal basis, not the only one. You can't expect
a _different_ method for constructing an orthonormal basis to make the
_same_ arbitrary choices as G-S!
> Here is an example, I calculate the Q by Gram-Schmidt process and
> compare the output with qr.Q()
>
>
> a <- c(1,0,1)
> b <- c(1,0,0)
> c <- c(2,1,0)
> x <- matrix(c(a,b,c),3,3)
>
> ##########################
> # Gram-Schmidt
> ##########################
>
> A <- matrix(a,3,1)
> q1 <- (1/sqrt(sum(A^2)))*A
> B <- b - (q1%*%b)%*%q1
> q2 <- (1/sqrt(sum(B^2)))*B
> C <- c - (q1%*%c)%*%q1 - (q2%*%c)%*%q2
> q3 <- (1/sqrt(sum(C^2)))*C
> Orthonormal.basis <- matrix(c(q1,q2,q3),3,3)
>> Orthonormal.basis
> [,1] [,2] [,3]
> [1,] 0.7071068 0.7071068 0
> [2,] 0.0000000 0.0000000 1
> [3,] 0.7071068 -0.7071068 0
>
>
> ##########################
> # QR Factorisation X = QR
> ##########################
>
> x.qr <- qr(x)
> Q <- qr.Q(x.qr)
> R <- qr.R(x.qr)
> X <- qr.X(x.qr)
>> Q
> [,1] [,2] [,3]
> [1,] -0.7071068 -0.7071068 0
> [2,] 0.0000000 0.0000000 1
> [3,] -0.7071068 0.7071068 0
>
>
> Thanks,
> cruz
>
> ______________________________________________
> R-help at r-project.org 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.
--
O__ ---- Peter Dalgaard Øster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk) FAX: (+45) 35327907
More information about the R-help
mailing list