B %*% t(B) = R , then solve for B

An embedded and charset-unspecified text was scrubbed...
Name: not available
URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20110412/4feaac31/attachment.pl>
There are easier solutions. Suppose you have a matrix A, such as:

### Use the info from lm() help
ctl <- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt <- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group <- gl(2,10,20, labels=c("Ctl","Trt"))
weight <- c(ctl, trt)
lm.D9 <- lm(weight ~ group)

A <- model.matrix(lm.D9)
L <- chol(crossprod(A))
all.equal(crossprod(L), crossprod(A))
[1] TRUE

I think this answers your remaining two questions below as ewll.
-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
Behalf Of Shawn Koppenhoefer
Sent: Tuesday, April 12, 2011 11:11 AM
To: r-help at r-project.org
Subject: Re: [R] B %*% t(B) = R , then solve for B
Importance: High

Hello.

I found a solution that may interest others.

Recall that my problem was how to use R to decompose a matrix into the
product of a matrix and its transpose. or, symbolically:

     For known matrix M (3x3 matrix) and unknown matrix F and its
transpose t(F)
     where  F * t(F) = M
     determine F

The solution using R seems to be :

     U=eigen(M)$vectors
     D=diag(x=eigen(M)$values)
     F=U %*% sqrt(D)

But now I have two new questions:
1. How can I find a solution where F is a triangular matrix.
2. How can I find solutions to non-square matrices?

/shawn

p.s. Here's a numerical example that demonstrates the above.

 > M
           [,1]       [,2]       [,3]
[1,] 0.6098601  0.2557882  0.1857773
[2,] 0.2557882  0.5127065 -0.1384238
[3,] 0.1857773 -0.1384238  0.9351089

 > U=eigen(M)$vectors
 > D=diag(x=eigen(M)$values)
 > F=U %*% sqrt(D)

 > *F %*% t(F)*
           [,1]       [,2]       [,3]
[1,] 0.6098601  0.2557882  0.1857773
[2,] 0.2557882  0.5127065 -0.1384238
[3,] 0.1857773 -0.1384238  0.9351089

	[[alternative HTML version deleted]]

______________________________________________
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.
For non-square matrices, see qr and svd.

Spencer
There are easier solutions. Suppose you have a matrix A, such as:

### Use the info from lm() help
ctl<- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt<- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group<- gl(2,10,20, labels=c("Ctl","Trt"))
weight<- c(ctl, trt)
lm.D9<- lm(weight ~ group)

A<- model.matrix(lm.D9)

L<- chol(crossprod(A))
all.equal(crossprod(L), crossprod(A))
[1] TRUE

I think this answers your remaining two questions below as ewll.

-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
Behalf Of Shawn Koppenhoefer
Sent: Tuesday, April 12, 2011 11:11 AM
To: r-help at r-project.org
Subject: Re: [R] B %*% t(B) = R , then solve for B
Importance: High

Hello.

I found a solution that may interest others.

Recall that my problem was how to use R to decompose a matrix into the
product of a matrix and its transpose. or, symbolically:

      For known matrix M (3x3 matrix) and unknown matrix F and its
transpose t(F)
      where  F * t(F) = M
      determine F

The solution using R seems to be :

      U=eigen(M)$vectors
      D=diag(x=eigen(M)$values)
      F=U %*% sqrt(D)

But now I have two new questions:
1. How can I find a solution where F is a triangular matrix.
2. How can I find solutions to non-square matrices?

/shawn

p.s. Here's a numerical example that demonstrates the above.

  >  M
            [,1]       [,2]       [,3]
[1,] 0.6098601  0.2557882  0.1857773
[2,] 0.2557882  0.5127065 -0.1384238
[3,] 0.1857773 -0.1384238  0.9351089

  >  U=eigen(M)$vectors
  >  D=diag(x=eigen(M)$values)
  >  F=U %*% sqrt(D)

  >  *F %*% t(F)*
            [,1]       [,2]       [,3]
[1,] 0.6098601  0.2557882  0.1857773
[2,] 0.2557882  0.5127065 -0.1384238
[3,] 0.1857773 -0.1384238  0.9351089

	[[alternative HTML version deleted]]

______________________________________________
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.
______________________________________________
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.

Spencer Graves, PE, PhD
President and Chief Operating Officer
Structure Inspection and Monitoring, Inc.
751 Emerson Ct.
San Jos?, CA 95126
ph:  408-655-4567
Thanks to both Doran and Spencer for your helpful answers.

@Doran: I still have to parse your answer, as I'm not experienced enough 
to see immediately how it relates to my M. Sorry for being daft! I don't 
immediately see how doing a linear regression on two vectors. I'm 
reading about design matrices now to try to understand your answer 
(?model.matrix). And when I look at A, I don't see anything that leads 
me to finding a diagonal matrix. :-(

@Spencer: Yup. SVD (Single Value Decomposition) gave me a solution 
(non-triangular unfortunately). From it, I got my F and t(F) however F 
is not triangular.
I didn't know about about QR which gives me orthogonal and upper 
triangular matrix. So this is closer!

In the meantime, I've learned about the "LU Decomposition"... which 
sounded promising:

library(Matrix)
lum=lu(M)
L=expand(lum)$L
U=expand(lum)$U
L%*%U

I get lower L and upper U triangular matrices,
However it is not the case that U is the transpose of L :-(
which is what I'm after (certainly Doran's solution gives me the answer 
once I figure it out).

3 x 3 Matrix of class "dgeMatrix"
           [,1]       [,2]       [,3]
[1,] 0.6098601  0.2557882  0.1857773
[2,] 0.2557882  0.5127065 -0.1384238
[3,] 0.1857773 -0.1384238  0.9351089
 > L
3 x 3 Matrix of class "dtrMatrix" (unitriangular)
      [,1]       [,2]       [,3]
[1,]  1.0000000          .          .
[2,]  0.4194211  1.0000000          .
[3,]  0.3046228 -0.5336215  1.0000000
 > U
3 x 3 Matrix of class "dtrMatrix"
      [,1]       [,2]       [,3]
[1,]  0.6098601  0.2557882  0.1857773
[2,]          .  0.4054235 -0.2163427
[3,]          .          .  0.7630718
 >

See?
U is not t(L) :-(

remember that I'm trying to get this solution:

       [,1]   [,2]  [,3]
[1,] 0.781  0.000 0.000
[2,] 0.328  0.637 0.000
[3,] 0.238 -0.341 0.873

/shawn