I'm not getting that difference:
> .Machine$double.eps
[1] 2.220446e-16
# so I don't think that explains the different behavior.
> WB<-as.matrix(read.table(file.choose()))
> tcp1 <- tcrossprod(WB)
> tcp2 <- crossprod(t(WB))
> tcp1
[,1] [,2] [,3]
[1,] 1916061939 2281366606 678696067
[2,] 2281366606 3098975504 1092911209
[3,] 678696067 1092911209 452399849
WBtWB
[,1] [,2] [,3]
[1,] 1916061939 2281366606 678696067
[2,] 2281366606 3098975504 1092911209
[3,] 678696067 1092911209 452399849
> tcp2
[,1] [,2] [,3]
[1,] 1916061939 2281366606 678696067
[2,] 2281366606 3098975504 1092911209
[3,] 678696067 1092911209 452399849
> solve(tcp1)
Error in solve.default(tcp1) :
system is computationally singular: reciprocal condition number =
9.60696e-17
> solve(tcp2)
Error in solve.default(tcp2) :
system is computationally singular: reciprocal condition number =
9.60696e-17
That is somewhat interesting since yesterday my machine solved my
input version of your:
WBtWB
[,1] [,2] [,3]
[1,] 1916061939 2281366606 678696067
[2,] 2281366606 3098975504 1092911209
[3,] 678696067 1092911209 452399849
I assume that despite those matrices being displayed as the same they
are represented differently in the machine.
Berry wrote:
I suppose that it is possible that the difference between what you
report
and what I see lies in the numerical libraries (LINPACK/LAPACK) that R
calls upon.
That would seem to be a possibility. You are using an out-of-date
version which may limit people's interest in investigating the problem.
--
David Winsemius
On Jan 15, 2009, at 12:47 AM, dos Reis, Marlon wrote:
Hi,
I attached the files I'm using, it may help.
I'm using Windows XP
sessionInfo()
R version 2.6.0 (2007-10-03)
i386-pc-mingw32
locale: LC_COLLATE=English_New Zealand.1252;LC_CTYPE=English_New
Zealand.1252;LC_MONETARY=English_New
Zealand.1252;LC_NUMERIC=C;LC_TIME=English_New Zealand.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
Try for example:
WB<-as.matrix(read.table("WB.dat"))
tcp1 <- tcrossprod(WB)
tcp2 <- crossprod(t(WB))
solve(tcp1)
[,1] [,2] [,3]
[1,] -41692.80 58330.89 -78368.17
[2,] 58330.89 -81608.66 109642.09
[3,] -78368.17 109642.09 -147305.32
solve(tcp2)
Error in solve.default(tcp2) :
system is computationally singular: reciprocal condition number =
2.17737e-17
Marlon.
-----Original Message-----
From: Charles C. Berry [mailto:cbe...@tajo.ucsd.edu]
Sent: Thursday, 15 January 2009 5:16 p.m.
To: dos Reis, Marlon
Cc: r-help@r-project.org
Subject: Re: [R] Precision in R
Marlon,
Are you using a current version of R? sessionInfo()?
It would help if you had something we could _fully_ reproduce.
Taking the _printed_ values you have below (WBtWB) and adding or
subtracting what you have printed as the difference of the two
visually
equal matrices ( say Delta ) , I am able to run
solve( dat3 )
solve( WBtWB + Delta )
solve( WBtWB - Delta )
solve( WBtWB + 2*Delta )
solve( WBtWB - 2*Delta )
and get the results to agree to 3 significant digits. And perturbing
things even more I still get solve() to return a value:
for ( i in 1:1000 ) solve(WBtWB - tcrossprod(rnorm(3)))
for ( i in 1:1000 ) solve(WBtWB + tcrossprod(rnorm(3)))
And I cannot get condition numbers anything like what you report:
range(replicate( 10000, 1/kappa(dat3-tcrossprod(matrix(rnorm(9),
3)))))
[1] 5.917764e-11 3.350445e-09
So I am very curious that you got the results that you print below.
I suppose that it is possible that the difference between what you
report
and what I see lies in the numerical libraries (LINPACK/LAPACK) that R
calls upon.
This was done on a windows XP PC. Here is my sessionInfo()
sessionInfo()
R version 2.8.1 Patched (2008-12-22 r47296)
i386-pc-mingw32
locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United
States.1252;LC_MONETARY=English_United
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
HTH,
Chuck
On Thu, 15 Jan 2009, dos Reis, Marlon wrote:
Dear All,
I'm preparing a simple algorithm for matrix multiplication for a
specific purpose, but I'm getting some unexpected results.
If anyone could give a clue, I would really appreciate.
Basically what I want to do is a simple matrix multiplication:
(WB) %*% t(WB).
The WB is in the disk so I compared to approaches:
- Load 'WB' using 'read.table' (put it in WB.tmp) and then to the
simple matrix multiplication
WB.tmp%*%t(WB.tmp)
- Scan each row of WB and do the cross products 'sum(WB.i*WB.i)'
and 'sum(WB.i*WB.j)', which proper arrangement leads to WBtWB.
Comparing these two matrices, I get the very similar values, however
when I tried their inverse, WBtWB leads to a singular system.
I've tried different tests and my conclusion is that my precision
problem is related to cross products 'sum(WB.i*WB.i)' and
'sum(WB.i*WB.j)'.
Does it makes sense?
Thanks,
Marlon
WB.tmp%*%t(WB.tmp)
WB.i WB.i WB.i
WB.i 1916061939 2281366606 678696067
WB.i 2281366606 3098975504 1092911209
WB.i 678696067 1092911209 452399849
WBtWB
[,1] [,2] [,3]
[1,] 1916061939 2281366606 678696067
[2,] 2281366606 3098975504 1092911209
[3,] 678696067 1092911209 452399849
WBtWB-WB.tmp%*%t(WB.tmp)
WB.i WB.i WB.i
WB.i 2.861023e-06 4.768372e-07 4.768372e-07
WB.i 4.768372e-07 3.814697e-06 -2.622604e-06
WB.i 4.768372e-07 -2.622604e-06 5.960464e-08
solve(WB.tmp%*%t(WB.tmp))
WB.i WB.i WB.i
WB.i -41692.80 58330.89 -78368.17
WB.i 58330.89 -81608.66 109642.09
WB.i -78368.17 109642.09 -147305.32
solve(WBtWB)
Error in solve.default(WBtWB) :
system is computationally singular: reciprocal condition number =
2.17737e-17
WB.tmp<-NULL
WBtWB<-matrix(NA,n,n)
for (i in 1:n)
{
setwd(Home.dir)
WB.i<-scan("WB.dat", skip = (i-1), nlines = 1)
WB.tmp<-rbind(WB.tmp,WB.i)
WBtWB[i,i]<-sum(WB.i*WB.i)
if (i<n)
{
for (j in (i+1):n)
{
setwd(Home.dir)
WB.j<-scan("WB.dat", skip = (j-1), nlines = 1)
WBtWB[i,j]<-sum(WB.i*WB.j)
WBtWB[j,i]<-sum(WB.i*WB.j)
}
}
}
=
======================================================================
Attention: The information contained in this message
and...{{dropped:15}}
______________________________________________
R-help@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.
Charles C. Berry (858) 534-2098
Dept of Family/Preventive
Medicine
E mailto:cbe...@tajo.ucsd.edu UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego
92093-0901
=
======================================================================
Attention: The information contained in this message and/or
attachments
from AgResearch Limited is intended only for the persons or entities
to which it is addressed and may contain confidential and/or
privileged
material. Any review, retransmission, dissemination or other use of,
or
taking of any action in reliance upon, this information by persons or
entities other than the intended recipients is prohibited by
AgResearch
Limited. If you have received this message in error, please notify the
sender immediately.
=
======================================================================
<WB.dat>
______________________________________________
R-help@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.
