[Rcpp-devel] symmatu does not give symmetric matrix

Yixuan Qiu yixuan.qiu at cos.name
Wed Jul 27 15:33:57 CEST 2016


Hi Takoua,

Are you sure the Cholesky decomposition was complaining about symmetry and
not something else?
In fact in your R example, the result *IS* a symmetric matrix. The two
numbers are just displayed differently (the first column displays 6 digits,
and the second column displays 7 digits), but the values in memory are the
same. You can try

## A better way is m = tcrossprod(x)
m = x %*% t(x)
format(m, digits = 10)
m[1, 2] == m[2, 1]

to verify this.


Best,
Yixuan


2016-07-21 6:46 GMT-04:00 Jendoubi, Takoua <t.jendoubi14 at imperial.ac.uk>:

> Dear all,
>
>
>
> I am using RcppArmadillo to deal with some matrix computations.
> Specifically I need to cholesky factorization of some symmetric matrices.
>
> I am generating random vectors using Rcpp and using them to construct
> symmetric matrices.
>
> I have an error stating that my matrix is not symmetric although it
> definitely should be. Here is the example I am working with (in R):
>
>
>
> >x
>
> [1] -1.6683320 -0.8597148
>
> >x%*%t(x)
>
>          [,1]      [,2]
>
> [1,] 2.783332 1.4342896
>
> [2,] 1.434290 0.7391095
>
>
>
> Apparently, it is a rounding-off error. Is there any way to ensure that
> x%*%t(x) gives an exactly symmetric matrix to use for cholesky
> factorization?
>
> I tried symmatu and symmatl but the error is still there.
>
>
>
> Thanks,
>
> Best,
>
>
>
> Takoua Jendoubi
>
> PhD student
>
> School of Public Health & Department of Surgery and Cancer
>
> Main address :
>
> Imperial College London | Room 171 | St Mary’s Campus | Norfolk Place | London
> | W2 1PG |
>
>
>
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel at lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
>



-- 
Yixuan Qiu <yixuan.qiu at cos.name>
Department of Statistics,
Purdue University
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20160727/2132e503/attachment.html>


More information about the Rcpp-devel mailing list