clustering 0.1.2 → 0.2.0

raw patch · 15 files changed

+434/−121 lines, 15 filesdep +Rlang-QQdep +paralleldep ~matricesPVP ok

version bump matches the API change (PVP)

Dependencies added: Rlang-QQ, parallel

Dependency ranges changed: matrices

API changes (from Hackage documentation)

Files

• benchmarks/Bench/Hierarchical.hs +44/−0
• benchmarks/Bench/KMeans.hs +31/−0
• benchmarks/Bench/Utils.hs +14/−0
• benchmarks/bench.hs +6/−27
• clustering.cabal +18/−4
• src/AI/Clustering/Hierarchical.hs +10/−11
• src/AI/Clustering/Hierarchical/Internal.hs +12/−0
• src/AI/Clustering/Hierarchical/Types.hs +20/−0
• src/AI/Clustering/KMeans.hs +65/−69
• src/AI/Clustering/KMeans/Internal.hs +101/−0
• src/AI/Clustering/KMeans/Types.hs +33/−0
• tests/Test/Hierarchical.hs +5/−9
• tests/Test/KMeans.hs +57/−0
• tests/Test/Utils.hs +14/−0
• tests/test.hs +4/−1

@@ -0,0 +1,44 @@+module Bench.Hierarchical+    ( benchHierarchical ) where++import Criterion.Main+import qualified Data.Clustering.Hierarchical as C+import qualified Data.Vector as V+import qualified Data.Vector.Unboxed as U+import System.IO.Unsafe (unsafePerformIO)++import AI.Clustering.Hierarchical+import AI.Clustering.Hierarchical.Types++import Bench.Utils++benchHierarchical :: Benchmark+benchHierarchical =+    let dists = computeDists euclidean xs+        fn i j = dists ! (i,j)+        xs = V.fromList $ unsafePerformIO $ randVectors 1000 5+    in bgroup "Hierarchical clustering"+        [ bgroup "AI.Clustering.Hierarchical"+            [ bench "Average Linkage (n = 10)" $+                  whnf (\x -> hclust Average x fn) $! U.enumFromN 0 10+            , bench "Average Linkage (n = 100)" $+                  whnf (\x -> hclust Average x fn) $! U.enumFromN 0 100+            , bench "Average Linkage (n = 500)" $+                  whnf (\x -> hclust Average x fn) $! U.enumFromN 0 500+            ]+        , bgroup "Data.Clustering.Hierarchical"+            [ bench "Average Linkage (n = 10)" $+                  whnf (\x -> C.dendrogram C.UPGMA x fn) $! [0..9]+            , bench "Average Linkage (n = 100)" $+                  whnf (\x -> C.dendrogram C.UPGMA x fn) $! [0..99]+            , bench "Average Linkage (n = 500)" $+                  whnf (\x -> C.dendrogram C.UPGMA x fn) $! [0..499]+            ]++        , bgroup "Distance matrix"+            [ bench "computeDists" $+                  whnf (\x -> computeDists euclidean x) xs+            , bench "computeDists'" $+                  whnf (\x -> computeDists' euclidean x) xs+            ]+        ]

@@ -0,0 +1,31 @@+module Bench.KMeans+    ( benchKMeans ) where++import Criterion.Main+import qualified Data.Matrix.Unboxed as MU+import qualified Data.Vector.Unboxed as U+import System.Random.MWC+import System.IO.Unsafe++import AI.Clustering.KMeans++import Bench.Utils++gen :: GenIO+gen = unsafePerformIO createSystemRandom++dat :: MU.Matrix Double+dat = unsafePerformIO $ fmap MU.fromRows $ randVectors 1000 10++benchKMeans :: Benchmark+benchKMeans = bgroup "KMeans clustering"+    [ bgroup "AI.Clustering.KMeans"+        [ bench "k-means++ (n = 1000, k = 7)" $+            whnfIO $ kmeans' gen KMeansPP 7 dat+        , bench "forgy (n = 1000, k = 7)" $+            whnfIO $ kmeans' gen Forgy 7 dat+        ]+    ]++kmeans' :: GenIO -> Method -> Int -> MU.Matrix Double -> IO (U.Vector Int)+kmeans' g method k = fmap _clusters . kmeans g method k

@@ -0,0 +1,14 @@+module Bench.Utils+    ( randVectors+    ) where++import Control.Monad (replicateM)+import qualified Data.Vector.Unboxed as U+import System.Random.MWC++randVectors :: Int  -- ^ number of samples+            -> Int  -- ^ vector length+            -> IO [U.Vector Double]+randVectors n k = do+    g <- createSystemRandom+    replicateM n $ uniformVector g k

@@ -7,32 +7,11 @@ import AI.Clustering.Hierarchical import AI.Clustering.Hierarchical.Types ((!)) -randSample :: IO [V.Vector Double]-randSample = do-    g <- create-    replicateM 2000 $ uniformVector g 5+import Bench.Hierarchical (benchHierarchical)+import Bench.KMeans (benchKMeans)  main :: IO ()-main = do-    xs <- randSample-    let dists = computeDists euclidean $ V.fromList xs-        fn i j = dists ! (i,j)-    defaultMain-      [ bgroup "AI.Clustering.Hierarchical"-          [ bench "Average Linkage (n = 10)" $-                whnf (\x -> hclust Average x fn) $! V.enumFromN 0 10-          , bench "Average Linkage (n = 100)" $-                whnf (\x -> hclust Average x fn) $! V.enumFromN 0 100-          , bench "Average Linkage (n = 1000)" $-                whnf (\x -> hclust Average x fn) $! V.enumFromN 0 1000-          ]--      , bgroup "Data.Clustering.Hierarchical"-          [ bench "Average Linkage (n = 10)" $-                whnf (\x -> C.dendrogram C.UPGMA x euclidean) $! take 10 xs-          , bench "Average Linkage (n = 100)" $-                whnf (\x -> C.dendrogram C.UPGMA x euclidean) $! take 100 xs-          , bench "Average Linkage (n = 1000)" $-                whnf (\x -> C.dendrogram C.UPGMA x euclidean) $! take 1000 xs-          ]-      ]+main = defaultMain+    [ benchHierarchical+    , benchKMeans+    ]

@@ -2,7 +2,7 @@ -- documentation, see http://haskell.org/cabal/users-guide/  name:                clustering-version:             0.1.2+version:             0.2.0 synopsis:            High performance clustering algorithms description:   Following clutering methods are included in this library:@@ -26,22 +26,26 @@ library   exposed-modules:          AI.Clustering.Hierarchical+    AI.Clustering.Hierarchical.Internal     AI.Clustering.Hierarchical.Types     AI.Clustering.KMeans+    AI.Clustering.KMeans.Internal+    AI.Clustering.KMeans.Types -  other-modules:       -    AI.Clustering.Hierarchical.Internal+--  other-modules:           build-depends:       base >=4.0 && <5.0     , binary     , containers-    , matrices+    , matrices >=0.4.0     , mwc-random+    , parallel     , primitive     , vector    hs-source-dirs:      src+  ghc-options:         -Wall   default-language:    Haskell2010  test-suite test@@ -50,12 +54,15 @@   main-is: test.hs   other-modules:     Test.Hierarchical+    Test.KMeans+    Test.Utils    default-language:    Haskell2010   build-depends:        base     , binary     , mwc-random+    , matrices     , vector     , tasty     , tasty-hunit@@ -63,11 +70,17 @@     , clustering     , hierarchical-clustering     , split+    , Rlang-QQ  benchmark bench   type: exitcode-stdio-1.0   hs-source-dirs: benchmarks+  ghc-options:  -threaded -rtsopts -with-rtsopts=-N2   main-is: bench.hs+  other-modules:+    Bench.Hierarchical+    Bench.KMeans+    Bench.Utils    default-language:    Haskell2010   build-depends: @@ -77,6 +90,7 @@     , vector     , clustering     , hierarchical-clustering+    , matrices  source-repository  head   type: git

@@ -46,8 +46,10 @@     , cutAt     , flatten     , drawDendrogram-    , computeDists++    -- * Distance functions     , euclidean+    , hamming      -- * References     -- $references@@ -56,7 +58,6 @@  import Control.Applicative ((<$>)) import qualified Data.Vector.Generic as G-import qualified Data.Vector.Unboxed as U import Text.Printf (printf)  import AI.Clustering.Hierarchical.Internal@@ -75,7 +76,7 @@ hclust :: G.Vector v a => Linkage -> v a -> DistFn a -> Dendrogram a hclust method xs f = label <$> nnChain dists fn   where-    dists = computeDists f xs+    dists = computeDists' f xs     label i = xs G.! i     fn = case method of         Single -> single@@ -109,17 +110,15 @@     draw (Leaf x) = [x,""]     shift first other = zipWith (++) (first : repeat other) -computeDists :: G.Vector v a => DistFn a -> v a -> DistanceMat-computeDists f vec = DistanceMat n . U.fromList . flip concatMap [0..n-1] $ \i ->-    flip map [i+1..n-1] $ \j -> f (vec `G.unsafeIndex` i) (vec `G.unsafeIndex` j)-  where-    n = G.length vec-{-# INLINE computeDists #-}---- | compute euclidean distance between two points+-- | Compute euclidean distance between two points. euclidean :: G.Vector v Double => DistFn (v Double) euclidean xs ys = sqrt $ G.sum $ G.zipWith (\x y -> (x-y)**2) xs ys {-# INLINE euclidean #-}++-- | Hamming distance.+hamming :: (G.Vector v a, G.Vector v Bool, Eq a) => DistFn (v a)+hamming xs = fromIntegral . G.length . G.filter id . G.zipWith (/=) xs+{-# INLINE hamming #-}  -- $references --

@@ -1,4 +1,16 @@+--------------------------------------------------------------------------------+-- |+-- Module      :  AI.Clustering.Hierarchical.Internal+-- Copyright   :  (c) 2015 Kai Zhang+-- License     :  MIT+--+-- Maintainer  :  kai@kzhang.org+-- Stability   :  experimental+-- Portability :  portable+--+-------------------------------------------------------------------------------- module AI.Clustering.Hierarchical.Internal+{-# WARNING "To be used by developer only" #-}     ( nnChain     , single     , complete

@@ -7,12 +7,16 @@     , DistanceMat(..)     , (!)     , idx+    , computeDists+    , computeDists'     ) where  import Control.Monad (liftM, liftM4)+import Control.Parallel.Strategies (rdeepseq, parMap) import Data.Binary (Binary, put, get, getWord8) import Data.Bits (shiftR) import qualified Data.Vector.Unboxed as U+import qualified Data.Vector.Generic as G import Data.Word (Word8)  type Distance = Double@@ -59,3 +63,19 @@ idx n i j | i <= j = (i * (2 * n - i - 3)) `shiftR` 1 + j - 1           | otherwise = (j * (2 * n - j - 3)) `shiftR` 1 + i - 1 {-# INLINE idx #-}++-- | compute distance matrix+computeDists :: G.Vector v a => DistFn a -> v a -> DistanceMat+computeDists f vec = DistanceMat n . U.fromList . flip concatMap [0..n-1] $ \i ->+    flip map [i+1..n-1] $ \j -> f (vec `G.unsafeIndex` i) (vec `G.unsafeIndex` j)+  where+    n = G.length vec+{-# INLINE computeDists #-}++-- | compute distance matrix in parallel+computeDists' :: G.Vector v a => DistFn a -> v a -> DistanceMat+computeDists' f vec = DistanceMat n . U.fromList . concat . flip (parMap rdeepseq) [0..n-1] $ \i ->+    flip map [i+1..n-1] $ \j -> f (vec `G.unsafeIndex` i) (vec `G.unsafeIndex` j)+  where+    n = G.length vec+{-# INLINE computeDists' #-}

@@ -1,3 +1,4 @@+{-# LANGUAGE FlexibleContexts #-} -------------------------------------------------------------------------------- -- | -- Module      :  AI.Clustering.KMeans@@ -12,54 +13,72 @@ module AI.Clustering.KMeans     ( KMeans(..)     , kmeans+    , kmeansBy     , kmeansWith      -- * Initialization methods-    , Initialization(..)+    , Method(..) +    -- * Useful functions     , decode+    , withinSS++    -- * References+    -- $references     ) where  import Control.Monad (forM_) import Control.Monad.Primitive (PrimMonad, PrimState) import qualified Data.Matrix.Unboxed as MU+import qualified Data.Matrix.Generic as MG import qualified Data.Matrix.Unboxed.Mutable as MM import Data.Ord (comparing) import qualified Data.Vector as V+import qualified Data.Vector.Generic as G import qualified Data.Vector.Mutable as VM import qualified Data.Vector.Unboxed as U import qualified Data.Vector.Unboxed.Mutable as UM-import Data.List (minimumBy, nub)-import System.Random.MWC (uniformR, Gen)+import Data.List (minimumBy, foldl')+import System.Random.MWC (Gen) --- | Results from running kmeans-data KMeans = KMeans-    { _clusters :: U.Vector Int     -- ^ A vector of integers (0 ~ k-1)-                                    -- indicating the cluster to which each-                                    -- point is allocated.-    , _centers :: MU.Matrix Double  -- ^ A matrix of cluster centres.-    } deriving (Show)+import AI.Clustering.KMeans.Types (KMeans(..), Method(..))+import AI.Clustering.KMeans.Internal (sumSquares, forgy, kmeansPP) --- | Lloyd's algorithm, also known as K-means algorithm-kmeans :: PrimMonad m+-- | Perform K-means clustering+kmeans :: (PrimMonad m, MG.Matrix mat U.Vector Double)        => Gen (PrimState m)-       -> Initialization-       -> Int                           -- ^ number of clusters-       -> MU.Matrix Double             -- ^ each row represents a point+       -> Method+       -> Int+       -> mat U.Vector Double        -> m KMeans-kmeans g method k mat = do-    initial <- case method of-        Forgy -> forgy g k mat-        _ -> undefined-    return $ kmeansWith initial mat+kmeans g method k mat = kmeansBy g method k dat (MG.takeRow mat)+  where+    dat = U.enumFromN 0 $ MG.rows mat {-# INLINE kmeans #-} --- | Lloyd's algorithm, also known as K-means algorithm-kmeansWith :: MU.Matrix Double   -- ^ initial set of k centroids-           -> MU.Matrix Double   -- ^ each row represents a point+-- | K-means algorithm+kmeansBy :: (PrimMonad m, G.Vector v a)+         => Gen (PrimState m)+         -> Method+         -> Int                   -- ^ number of clusters+         -> v a                   -- ^ data stores in rows+         -> (a -> U.Vector Double)+         -> m KMeans+kmeansBy g method k dat fn = do+    initial <- case method of+        Forgy -> forgy g k dat fn+        KMeansPP -> kmeansPP g k dat fn+    return $ kmeansWith initial dat fn+{-# INLINE kmeansBy #-}++-- | K-means algorithm+kmeansWith :: G.Vector v a+           => MU.Matrix Double   -- ^ initial set of k centroids+           -> v a                -- ^ each row represents a point+           -> (a -> U.Vector Double)            -> KMeans-kmeansWith initial mat | d /= MU.cols initial || k > n = error "check input"-                       | otherwise = KMeans member centers+kmeansWith initial dat fn | d /= MU.cols initial || k > n = error "check input"+                          | otherwise = KMeans member centers   where     (member, centers) = loop initial U.empty     loop means membership@@ -70,18 +89,20 @@      -- Assignment step     assign means = U.generate n $ \i ->-        let x = MU.takeRow mat i-        in fst $ minimumBy (comparing snd) $ zip [0..k-1] $ map (dist x) $ MU.toRows means+        let x = fn $ G.unsafeIndex dat i+        in fst $ minimumBy (comparing snd) $ zip [0..k-1] $ map (sumSquares x) $ MU.toRows means      -- Update step-    update membership = MM.create $ do+    update membership = MU.create $ do         m <- MM.replicate (k,d) 0.0         count <- UM.replicate k (0 :: Int)         forM_ [0..n-1] $ \i -> do             let x = membership U.! i             UM.unsafeRead count x >>= UM.unsafeWrite count x . (+1)++            let vec = fn $ dat G.! i             forM_ [0..d-1] $ \j ->-                MM.unsafeRead m (x,j) >>= MM.unsafeWrite m (x,j) . (+ mat MU.! (i,j))+                MM.unsafeRead m (x,j) >>= MM.unsafeWrite m (x,j) . (+ (vec U.! j))         -- normalize         forM_ [0..k-1] $ \i -> do             c <- UM.unsafeRead count i@@ -89,48 +110,11 @@                 MM.unsafeRead m (i,j) >>= MM.unsafeWrite m (i,j) . (/fromIntegral c)         return m -    dist xs = U.sum . U.zipWith (\x y -> (x - y)**2) xs--    n = MU.rows mat+    n = G.length dat     k = MU.rows initial-    d = MU.cols mat+    d = MU.cols initial {-# INLINE kmeansWith #-} --- | Different initialization methods-data Initialization = Forgy    -- ^ The Forgy method randomly chooses k unique-                               -- observations from the data set and uses these-                               -- as the initial means-                    | KMeansPP -- ^ K-means++ algorithm, not implemented.--forgy :: PrimMonad m-      => Gen (PrimState m)-      -> Int                 -- number of clusters-      -> MU.Matrix Double    -- data-      -> m (MU.Matrix Double)-forgy g k mat | k > n = error "k is larger than sample size"-              | otherwise = iter-  where-    iter = do-        vec <- sample g k . U.enumFromN 0 $ n-        let xs = map (MU.takeRow mat) . U.toList $ vec-        if length (nub xs) == length xs-           then return . MU.fromRows $ xs-           else iter-    n = MU.rows mat-{-# INLINE forgy #-}---- random select k samples from a population-sample :: PrimMonad m => Gen (PrimState m) -> Int -> U.Vector Int -> m (U.Vector Int)-sample g k xs = do-    v <- U.thaw xs-    forM_ [0..k-1] $ \i -> do-        j <- uniformR (i, lst) g-        UM.unsafeSwap v i j-    U.unsafeFreeze . UM.take k $ v-  where-    lst = U.length xs - 1-{-# INLINE sample #-}- -- | Assign data to clusters based on KMeans result decode :: KMeans -> [a] -> [[a]] decode result xs = V.toList $ V.create $ do@@ -143,4 +127,16 @@     n = U.maximum membership + 1  -- | Compute within-cluster sum of squares---withinSS :: Matrix+withinSS :: KMeans -> MU.Matrix Double -> [Double]+withinSS result mat = zipWith f (decode result [0 .. MU.rows mat-1]) .+                          MU.toRows . _centers $ result+  where+    f c center = foldl' (+) 0 $ map (sumSquares center . MU.takeRow mat) c+++-- $references+--+-- Arthur, D. and Vassilvitskii, S. (2007). k-means++: the advantages of careful +-- seeding. Proceedings of the eighteenth annual ACM-SIAM symposium on Discrete +-- algorithms. Society for Industrial and Applied Mathematics Philadelphia, PA, +-- USA. pp. 1027–1035.

@@ -0,0 +1,101 @@+{-# LANGUAGE BangPatterns #-}+--------------------------------------------------------------------------------+-- |+-- Module      :  AI.Clustering.KMeans.Internal+-- Copyright   :  (c) 2015 Kai Zhang+-- License     :  MIT+--+-- Maintainer  :  kai@kzhang.org+-- Stability   :  experimental+-- Portability :  portable+--+-- <module description starting at first column>+--------------------------------------------------------------------------------+module AI.Clustering.KMeans.Internal+{-# WARNING "To be used by developer only" #-}+    ( forgy+    , kmeansPP+    , sumSquares+    ) where++import Control.Monad (forM_)+import Control.Monad.Primitive (PrimMonad, PrimState)+import Data.List (nub)+import qualified Data.Matrix.Unboxed as MU+import qualified Data.Vector.Generic as G+import qualified Data.Vector.Unboxed as U+import qualified Data.Vector.Unboxed.Mutable as UM+import System.Random.MWC (uniformR, Gen)++forgy :: (PrimMonad m, G.Vector v a)+      => Gen (PrimState m)+      -> Int                 -- number of clusters+      -> v a                 -- data+      -> (a -> U.Vector Double)+      -> m (MU.Matrix Double)+forgy g k dat fn | k > n = error "k is larger than sample size"+                 | otherwise = iter+  where+    iter = do+        vec <- randN g k . U.enumFromN 0 $ n+        let xs = map (\i -> fn $ dat `G.unsafeIndex` i) . U.toList $ vec+        if length (nub xs) == length xs+           then return . MU.fromRows $ xs+           else iter+    n = G.length dat+{-# INLINE forgy #-}++kmeansPP :: (PrimMonad m, G.Vector v a)+         => Gen (PrimState m)+         -> Int+         -> v a+         -> (a -> U.Vector Double)+         -> m (MU.Matrix Double)+kmeansPP g k dat fn+    | k > n = error "k is larger than sample size"+    | otherwise = do+        c1 <- uniformR (0,n-1) g+        loop [c1] 1+  where+    loop centers !k'+        | k' == k = return $ MU.fromRows $ map (\i -> fn $ dat `G.unsafeIndex` i) centers+        | otherwise = do+            c' <- chooseWithProb g $ U.map (shortestDist centers) rowIndices+            loop (c':centers) (k'+1)++    n = G.length dat+    rowIndices = U.enumFromN 0 n+    shortestDist centers x = minimum $ map (\i ->+        sumSquares (fn $ dat `G.unsafeIndex` x) (fn $ dat `G.unsafeIndex` i)) centers+{-# INLINE kmeansPP #-}++chooseWithProb :: PrimMonad m+               => Gen (PrimState m)+               -> U.Vector Double    -- ^ weights, may not be normalized+               -> m Int              -- ^ result/index+chooseWithProb g ws = do+    x <- uniformR (0,sum') g+    return $ loop x 0 0+  where+    loop v !cdf !i | cdf' >= v = i+                   | otherwise = loop v cdf' (i+1)+      where cdf' = cdf + ws `U.unsafeIndex` i++    sum' = U.sum ws+{-# INLINE chooseWithProb #-}++-- | Random select k samples from a population+randN :: PrimMonad m => Gen (PrimState m) -> Int -> U.Vector Int -> m (U.Vector Int)+randN g k xs = do+    v <- U.thaw xs+    forM_ [0..k-1] $ \i -> do+        j <- uniformR (i, lst) g+        UM.unsafeSwap v i j+    U.unsafeFreeze . UM.take k $ v+  where+    lst = U.length xs - 1+{-# INLINE randN #-}++sumSquares :: U.Vector Double -> U.Vector Double -> Double+sumSquares xs = U.sum . U.zipWith (\x y -> (x - y)**2) xs+{-# INLINE sumSquares #-}

@@ -0,0 +1,33 @@+--------------------------------------------------------------------------------+-- |+-- Module      :  AI.Clustering.KMeans.Types+-- Copyright   :  (c) 2015 Kai Zhang+-- License     :  MIT+--+-- Maintainer  :  kai@kzhang.org+-- Stability   :  experimental+-- Portability :  portable+--+-- <module description starting at first column>+--------------------------------------------------------------------------------+module AI.Clustering.KMeans.Types+    ( KMeans(..)+    , Method(..)+    ) where++import qualified Data.Matrix.Unboxed as MU+import qualified Data.Vector.Unboxed as U++-- | Results from running kmeans+data KMeans = KMeans+    { _clusters :: U.Vector Int     -- ^ A vector of integers (0 ~ k-1)+                                    -- indicating the cluster to which each+                                    -- point is allocated.+    , _centers :: MU.Matrix Double  -- ^ A matrix of cluster centers.+    } deriving (Show)++-- | Different initialization methods+data Method = Forgy    -- ^ The Forgy method randomly chooses k unique+                       -- observations from the data set and uses these+                       -- as the initial means.+            | KMeansPP -- ^ K-means++ algorithm.

@@ -13,6 +13,7 @@ import Test.Tasty.QuickCheck  import AI.Clustering.Hierarchical+import Test.Utils  tests :: TestTree tests = testGroup "Hierarchical:"@@ -23,11 +24,6 @@     , testCase "Weighted Linkage" testWeighted     ] -randSample :: IO [V.Vector Double]-randSample = do-    g <- create-    replicateM 500 $ uniformVector g 5- isEqual :: Eq a => Dendrogram a -> C.Dendrogram a -> Bool isEqual (Leaf x) (C.Leaf x') = x == x' isEqual (Branch _ d x y) (C.Branch d' x' y') = abs (d - d') < 1e-8 &&@@ -36,7 +32,7 @@  testSingle :: Assertion testSingle = do-    xs <- randSample+    xs <- randVectors 500 5     let true = C.dendrogram C.SingleLinkage xs euclidean         test = hclust Single (V.fromList xs) euclidean     assertBool (unlines ["Expect: ", show true, "But see: ", show test]) $@@ -44,7 +40,7 @@  testComplete :: Assertion testComplete = do-    xs <- randSample+    xs <- randVectors 500 5     let true = C.dendrogram C.CompleteLinkage xs euclidean         test = hclust Complete (V.fromList xs) euclidean     assertBool (unlines ["Expect: ", show true, "But see: ", show test]) $@@ -52,7 +48,7 @@  testAverage :: Assertion testAverage = do-    xs <- randSample+    xs <- randVectors 500 5     let true = C.dendrogram C.UPGMA xs euclidean         test = hclust Average (V.fromList xs) euclidean     assertBool (unlines ["Expect: ", show true, "But see: ", show test]) $@@ -60,7 +56,7 @@  testWeighted :: Assertion testWeighted = do-    xs <- randSample+    xs <- randVectors 500 5     let true = C.dendrogram C.FakeAverageLinkage xs euclidean         test = hclust Weighted (V.fromList xs) euclidean     assertBool (unlines ["Expect: ", show true, "But see: ", show test]) $

@@ -0,0 +1,57 @@+{-# LANGUAGE QuasiQuotes #-}+{-# LANGUAGE DataKinds #-}+module Test.KMeans+    ( tests+    ) where++import Control.Monad+import qualified Data.Matrix.Unboxed as MU+import qualified Data.Vector.Unboxed as V+import Data.List+import RlangQQ+import System.Random.MWC+import Test.Tasty+import Test.Tasty.HUnit+import Test.Tasty.QuickCheck++import AI.Clustering.KMeans+import AI.Clustering.KMeans.Internal++import Test.Utils++tests :: TestTree+tests = testGroup "KMeans:"+    [ testCase "KMeans" testKMeans+    ]++rKmeans :: Int -> [Double] -> [Double] -> IO [Int]+rKmeans n dat center = do+    o <- [r| x <- matrix(hs_dat, ncol=hs_n,byrow=T);+             y <- matrix(hs_center, ncol=hs_n,byrow=T);+             hs_result <- kmeans(x,y,iter.max=1000000,algorithm="Lloyd")$cluster;+         |]+    let x = Label :: Label "result"+    return $ o .!. x++testKMeans :: Assertion+testKMeans = do+    let n = 2000+        d = 15+        k = 10+    g <- createSystemRandom+    xs <- randVectors n d++    let mat = MU.fromRows xs :: MU.Matrix Double+        dat = V.enumFromN 0 $ MU.rows mat+        fn = MU.takeRow mat++    centers <- kmeansPP g k dat fn++    r <- rKmeans d (MU.toList mat) (MU.toList centers)+    let test = sort $ map sort $ decode result xs+        result = kmeansWith centers dat fn+        true = sort $ map sort $ decode result{_clusters=V.fromList $ map (subtract 1) r} xs+        show' xs = unlines $ map (show . map (unwords . map show . V.toList)) xs++    assertBool ("Expect: " ++ show' true ++ "\nBut saw: " ++ show' test) $+        test == true

@@ -0,0 +1,14 @@+module Test.Utils+    ( randVectors+    ) where++import Control.Monad (replicateM)+import qualified Data.Vector.Unboxed as U+import System.Random.MWC++randVectors :: Int  -- ^ number of samples+            -> Int  -- ^ vector length+            -> IO [U.Vector Double]+randVectors n k = do+    g <- createSystemRandom+    replicateM n $ uniformVector g k

@@ -1,7 +1,10 @@ import Test.Tasty  import qualified Test.Hierarchical as Hierarchical+import qualified Test.KMeans as KMeans  main :: IO () main = defaultMain $ testGroup "Main"-    [ Hierarchical.tests ]+    [ Hierarchical.tests+    , KMeans.tests+    ]