1

I am trying to follow some simple hmatrix examples:

import qualified Numeric.LinearAlgebra as La

w = 4 La.|> [2, 0, -3, 0 :: Double]
m = (3 La.>< 4) [1::Double ..]
x = La.singularValues m
x' = sqrt . La.eigenvalues $ m La.<> La.trans m

Load into ghci:

> x'
fromList [25.436835633480243 :+ 0.0,1.7226122475210675 :+ 0.0,2.3006847884268434e-7 :+ 0.0]

Results in R:

> x = matrix(1:12, 3, 4, byrow=T)
> x
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    5    6    7    8
[3,]    9   10   11   12
> sqrt(eigen(x %*% t(x))$values)
[1] 2.543684e+01 1.722612e+00 1.159372e-07

The first two eigenvalues agree, but they differ quite a bit in the third. Why?

Note that hmatrix uses complex numbers.

qed
  • 22,298
  • 21
  • 125
  • 196
  • 3
    1e-7 is a small number. They agree to 6 d.p. – cchalmers Mar 27 '15 at 17:03
  • 2
    Since hmatrix is just a wrapper over LAPACK, BLAS and GSL, the disagreement between the singular values is due to the algorithm and precision use for the singular-value decomposition (svd). Have you seen that R provides 2 different function for it : svd and La.svd ? – felipez Mar 27 '15 at 17:17
  • wow, indeed, there is a La.svd in R. – qed Mar 27 '15 at 19:58

0 Answers0