summaryrefslogtreecommitdiff
path: root/examples/pca2.hs
diff options
context:
space:
mode:
Diffstat (limited to 'examples/pca2.hs')
-rw-r--r--examples/pca2.hs17
1 files changed, 10 insertions, 7 deletions
diff --git a/examples/pca2.hs b/examples/pca2.hs
index e7ea95f..892d382 100644
--- a/examples/pca2.hs
+++ b/examples/pca2.hs
@@ -1,5 +1,7 @@
1-- Improved PCA, including illustrative graphics 1-- Improved PCA, including illustrative graphics
2 2
3{-# LANGUAGE FlexibleContexts #-}
4
3import Numeric.LinearAlgebra 5import Numeric.LinearAlgebra
4import Graphics.Plot 6import Graphics.Plot
5import System.Directory(doesFileExist) 7import System.Directory(doesFileExist)
@@ -10,27 +12,27 @@ type Vec = Vector Double
10type Mat = Matrix Double 12type Mat = Matrix Double
11 13
12-- Vector with the mean value of the columns of a matrix 14-- Vector with the mean value of the columns of a matrix
13mean a = constant (recip . fromIntegral . rows $ a) (rows a) <> a 15mean a = konst (recip . fromIntegral . rows $ a) (rows a) <# a
14 16
15-- covariance matrix of a list of observations stored as rows 17-- covariance matrix of a list of observations stored as rows
16cov x = (trans xc <> xc) / fromIntegral (rows x - 1) 18cov x = (mTm xc) -- / fromIntegral (rows x - 1)
17 where xc = x - asRow (mean x) 19 where xc = x - asRow (mean x)
18 20
19 21
20type Stat = (Vec, [Double], Mat) 22type Stat = (Vec, [Double], Mat)
21-- 1st and 2nd order statistics of a dataset (mean, eigenvalues and eigenvectors of cov) 23-- 1st and 2nd order statistics of a dataset (mean, eigenvalues and eigenvectors of cov)
22stat :: Mat -> Stat 24stat :: Mat -> Stat
23stat x = (m, toList s, trans v) where 25stat x = (m, toList s, tr v) where
24 m = mean x 26 m = mean x
25 (s,v) = eigSH' (cov x) 27 (s,v) = eigSH (cov x)
26 28
27-- creates the compression and decompression functions from the desired reconstruction 29-- creates the compression and decompression functions from the desired reconstruction
28-- quality and the statistics of a data set 30-- quality and the statistics of a data set
29pca :: Double -> Stat -> (Vec -> Vec , Vec -> Vec) 31pca :: Double -> Stat -> (Vec -> Vec , Vec -> Vec)
30pca prec (m,s,v) = (encode,decode) 32pca prec (m,s,v) = (encode,decode)
31 where 33 where
32 encode x = vp <> (x - m) 34 encode x = vp #> (x - m)
33 decode x = x <> vp + m 35 decode x = x <# vp + m
34 vp = takeRows n v 36 vp = takeRows n v
35 n = 1 + (length $ fst $ span (< (prec'*sum s)) $ cumSum s) 37 n = 1 + (length $ fst $ span (< (prec'*sum s)) $ cumSum s)
36 cumSum = tail . scanl (+) 0.0 38 cumSum = tail . scanl (+) 0.0
@@ -46,7 +48,7 @@ test :: Stat -> Double -> Vec -> IO ()
46test st prec x = do 48test st prec x = do
47 let (pe,pd) = pca prec st 49 let (pe,pd) = pca prec st
48 let y = pe x 50 let y = pe x
49 print $ dim y 51 print $ size y
50 shdigit (pd y) 52 shdigit (pd y)
51 53
52main = do 54main = do
@@ -63,3 +65,4 @@ main = do
63 let st = stat xs 65 let st = stat xs
64 test st 0.90 x 66 test st 0.90 x
65 test st 0.50 x 67 test st 0.50 x
68