I have a symmetric matrix B (17x17), and a (17x17) square matrix A. If do
the following matrix multiplication I SHOULD get a symmetric matrix, however
i don't. The computation required is:
C = t(A)%*%B%*%A
here are some checks for symmetry
-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
On Behalf Of statfan
Sent: Saturday, October 22, 2011 10:45 PM
To: r-help at r-project.org
Subject: [R] symmetric matrix multiplication
I have a symmetric matrix B (17x17), and a (17x17) square matrix A. If do
the following matrix multiplication I SHOULD get a symmetric matrix,
however
i don't. The computation required is:
C = t(A)%*%B%*%A
here are some checks for symmetry
(max(abs(B - t(B))))
[1] 0
C = t(A)%*%B%*%A
(max(abs(C - t(C))))
[1] 3.552714e-15
Any help on the matter would be very much appreciated.
Welcome to the world of floating-point calculation on finite precision computers. You need to read R FAQ 7.31. Your maximum difference is for all intents and purposes equal to zero.
Hope this is helpful,
Dan
Daniel Nordlund
Bothell, WA USA
-----Original Message-----
From: r-help-bounces at r-project.org
[mailto:r-help-bounces at r-project.org]
On Behalf Of statfan
Sent: Saturday, October 22, 2011 10:45 PM
To: r-help at r-project.org
Subject: [R] symmetric matrix multiplication
I have a symmetric matrix B (17x17), and a (17x17) square matrix A.
If do
the following matrix multiplication I SHOULD get a symmetric matrix,
however
i don't. The computation required is:
C = t(A)%*%B%*%A
here are some checks for symmetry
(max(abs(B - t(B))))
[1] 0
C = t(A)%*%B%*%A
(max(abs(C - t(C))))
[1] 3.552714e-15
Any help on the matter would be very much appreciated.
Welcome to the world of floating-point calculation on
finite precision computers. You need to read R FAQ 7.31.
Your maximum difference is for all intents and purposes
equal to zero.
Hope this is helpful,
Dan
Daniel Nordlund
Bothell, WA USA
In addition to Dan's comment, let me point out that you
can convert your very nearly symmetric matrix C to an
exactly (even by R's finite-precision standards) symmetric
matrix by using (C + t(C))/2. The result will differ from
the original matrix C by similar "for all intents and purposes
zero" amounts. Here is an example, using 4x4 matrices:
##[1]: The symmetric matrix B:
B <- matrix( c(
1.1, 1.2, 1.3, 1.4,
1.2, 2.2, 2.3, 2.4,
1.3, 2.3, 3.3, 3.4,
1.4, 2.4, 3.4, 4.4), byrow=TRUE, nrow=4 )
##[2]: The non-symmetric matrix B:
A <- matrix( c(
1.1, 1.2, 1.3, 1.4,
2.1, 2.2, 2.3, 2.4,
3.1, 3.2, 3.3, 3.4,
4.1, 4.2, 4.3, 4.4), byrow=TRUE, nrow=4 )
##[3]: An allegedly symmetric matrix C1 (constructed
## like your C):
C1 <- t(A)%*%B%*%A
##[4]: But it isn't exactly symmetric:
max(abs(C1 - t(C1)))
# [1] 5.684342e-14
##[5]: So construct an exactly symmetric version:
C2 <- (C1 + t(C1))/2
##[6]: Check that it is exactly symmetric:
max(abs(C2 - t(C2)))
# [1] 0
##[7]: And check how close it is to the original C1:
max(abs(C2 - C1))
# 1] 5.684342e-14
Hoping this helps!
Ted.
--------------------------------------------------------------------
E-Mail: (Ted Harding) <ted.harding at wlandres.net>
Fax-to-email: +44 (0)870 094 0861
Date: 23-Oct-11 Time: 08:43:18
------------------------------ XFMail ------------------------------
Thank you Dan and Ted for these helpful comments. I will implement this
simple force symmetry code you suggested and make sure I familiarize with
this floating-point calculation problem so I can recognize such issues in
the future.
--
View this message in context: http://r.789695.n4.nabble.com/symmetric-matrix-multiplication-tp3929739p3930832.html
Sent from the R help mailing list archive at Nabble.com.