spectral-clustering 0.2.1.4 → 0.2.2.0

raw patch · 3 files changed

+164/−14 lines, 3 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• spectral-clustering.cabal +2/−2
• src/Math/Clustering/Spectral/Dense.hs +162/−11
• src/Math/Clustering/Spectral/Sparse.hs +0/−1

@@ -1,9 +1,9 @@ cabal-version: >=1.10 name: spectral-clustering-version: 0.2.1.4+version: 0.2.2.0 license: GPL-3 license-file: LICENSE-copyright: 2018 Gregory W. Schwartz+copyright: 2019 Gregory W. Schwartz maintainer: gsch@mail.med.upenn.edu author: Gregory W. Schwartz homepage: http://github.com/GregorySchwartz/spectral-clustering#readme

@@ -4,6 +4,8 @@ Collects the functions pertaining to spectral clustering. -} +{-# LANGUAGE BangPatterns #-}+ module Math.Clustering.Spectral.Dense     ( spectralClusterKNorm     , spectralClusterNorm@@ -11,6 +13,15 @@     , getDegreeMatrix     , AdjacencyMatrix (..)     , LabelVector (..)+    , B (..)+    , B1 (..)+    , B2 (..)+    , spectral+    , spectralCluster+    , spectralClusterK+    , getB+    , b1ToB2+    , getSimilarityFromB2     ) where  -- Remote@@ -21,15 +32,165 @@ import Safe (headMay) import qualified AI.Clustering.KMeans as K import qualified Data.Vector as V+import qualified Data.Vector.Storable as VS import qualified Data.Vector.Unboxed as U import qualified Numeric.LinearAlgebra as H+import qualified Numeric.LinearAlgebra.Devel as H import qualified Statistics.Quantile as S+import qualified Numeric.LinearAlgebra.SVD.SVDLIBC as SVD  -- Local +-- | Output vector containing cluster assignment (0 or 1). type LabelVector     = H.Vector Double++-- | Adjacency matrix input. type AdjacencyMatrix = H.Matrix Double +-- | B1 observation by feature matrix.+newtype B1 = B1 { unB1 :: H.Matrix Double } deriving (Show)+-- | B2 term frequency-inverse document frequency matrix of B1.+newtype B2 = B2 { unB2 :: H.Matrix Double } deriving (Show)+-- | Diagonal matrix from \(diag(B(B^{T}1))\).+newtype D  = D { unD :: H.Matrix Double } deriving (Show)+-- | Matrix from \(D^{-1/2}B}\).+newtype C  = C { unC :: H.Matrix Double } deriving (Show)+-- | Normed rows of B2. For a complete explanation, see Shu et al., "Efficient+-- Spectral Neighborhood Blocking for Entity Resolution", 2011.+newtype B  = B { unB :: H.Matrix Double } deriving (Show)++-- | Map hmatrix with indices.+cimap :: (Int -> Int -> Double -> Double) -> H.Matrix Double -> H.Matrix Double+cimap f mat = H.assoc (H.size mat) 0+            . concatMap (\ (!i, xs)+                        -> fmap (\ (!j, !x)+                                -> ( (i, j)+                                  , f i j x+                                  )+                                )+                                xs+                        )+            . zip [0..]+            . fmap (zip [0..])+            . H.toLists+            $ mat++-- | Normalize the input matrix by column. Here, columns are features.+b1ToB2 :: B1 -> B2+b1ToB2 (B1 b1) =+    B2+      . cimap (\ !i !j !x -> (log (fromIntegral n / (fromMaybe 0 $ dVec VS.!? j))) * x)+      $ b1+  where+    dVec :: H.Vector Double+    dVec = H.fromList+         . fmap (H.sumElements . H.step)+         . H.toColumns+         $ b1+    n = H.rows b1+    m = H.cols b1++-- | Euclidean norm each row.+b2ToB :: B2 -> B+b2ToB (B2 b2) =+    B . cimap (\ !i !j !x -> x / (fromMaybe 0 $ eVec VS.!? i)) $ b2+  where+    eVec :: H.Vector Double+    eVec = H.fromList . fmap H.norm_2 . H.toRows $ b2+    n = H.rows b2+    m = H.cols b2++-- | Get the diagonal transformed B matrix.+bToD :: B -> D+bToD (B b) = D+           . H.diag+           . H.flatten+           $ b+        H.<> ((H.tr b) H.<> ((n H.>< 1) [1,1..]))+  where+    n = H.rows b++-- | Get the matrix C as input for SVD.+bdToC :: B -> D -> C+bdToC (B b) (D d) = C $ (H.diag . H.cmap (\x -> x ** (- 1 / 2)) . H.takeDiag $ d) H.<> b++-- | Obtain the second left singular vector (or N earlier) and E on of a sparse+-- matrix.+secondLeft :: Int -> Int -> H.Matrix Double -> [H.Vector Double]+secondLeft n e m =+  fmap (VS.drop (n - 1))+    . H.toColumns+    . (\(!x, _, _) -> x)+    . SVD.sparseSvd (e + (n - 1))+    . H.mkCSR+    . filter (\((_, _), x) -> x /= 0)+    . concatMap (\(!i, xs) -> fmap (\(!j, !x) -> ((i, j), x)) xs)+    . zip [0..]+    . fmap (zip [0..])+    . H.toLists+    $ m++-- | Get the normalized matrix B from an input matrix where the features are+-- columns and rows are observations. Optionally, do not normalize.+getB :: Bool -> H.Matrix Double -> B+getB True = b2ToB . b1ToB2 . B1+getB False = b2ToB . B2++-- | Returns the second left singular vector (or from N) and E on of a sparse+-- spectral process. Assumes the columns are features and rows are observations.+-- B is the normalized matrix (from getB). See Shu et al., "Efficient Spectral+-- Neighborhood Blocking for Entity Resolution", 2011.+spectral :: Int -> Int -> B -> [H.Vector Double]+spectral n e b+    | e < 1     = error "Less than 1 eigenvector chosen for clustering."+    | n < 1 = error "N < 1, cannot go before first eigenvector."+    | otherwise = secondLeft n e . unC . bdToC b . bToD $ b++-- | Returns a vector of cluster labels for two groups by finding the second+-- left singular vector of a special normalized matrix. Assumes the columns are+-- features and rows are observations. B is the normalized matrix (from getB).+-- See Shu et al., "Efficient Spectral Neighborhood Blocking for Entity+-- Resolution", 2011.+spectralCluster :: B -> LabelVector+spectralCluster (B b)+  | H.rows b < 1  = H.fromList []+  | H.rows b == 1 = H.fromList [0]+  | otherwise     = H.cmap (bool 0 1 . (>= 0))+                  . mconcat+                  . spectral 2 1+                  $ B b++-- | Returns a vector of cluster labels for two groups by finding the second+-- left singular vector and on of a special normalized matrix and running kmeans.+-- Assumes the columns are features and rows are observations. B is the+-- normalized matrix (from getB). See Shu et al., "Efficient Spectral+-- Neighborhood Blocking for Entity Resolution", 2011.+spectralClusterK :: Int -> Int -> B -> LabelVector+spectralClusterK e k (B b)+  | H.rows b < 1  = H.fromList []+  | H.rows b == 1 = H.fromList [0]+  | otherwise     = kmeansVec k . spectral 1 e $ B b++-- | Executes kmeans to cluster a vector.+kmeansVec :: Int -> [H.Vector Double] -> LabelVector+kmeansVec k = V.convert+            . U.map fromIntegral+            . K.membership+            . (\x -> K.kmeansBy k x id K.defaultKMeansOpts)+            . V.fromList+            . fmap V.convert+            . H.toRows+            . H.fromColumns+            . fmap H.normalize -- Normalize within eigenvectors (columns).+            . H.toColumns+            . H.fromRows++-- | Get the cosine similarity between two rows using B2.+getSimilarityFromB2 :: B2 -> Int -> Int -> Double+getSimilarityFromB2 (B2 b2) i j =+    H.dot (H.flatten $ b2 H.? [i]) (H.flatten $ b2 H.? [j])+        / (H.norm_2 (H.flatten $ b2 H.? [i]) * H.norm_2 (H.flatten $ b2 H.? [j]))+ -- | Returns the clustering of eigenvectors with the second smallest eigenvalues -- and on of the symmetric normalized Laplacian L. Computes real symmetric part -- of L, so ensure the input is real and symmetric. Diagonal should be 0s for@@ -39,17 +200,7 @@ spectralClusterKNorm e k mat   | H.rows mat < 1  = H.fromList []   | H.rows mat == 1 = H.fromList [0]-  | otherwise       = V.convert-                    . U.map fromIntegral-                    . K.membership-                    . (\x -> K.kmeansBy k x id K.defaultKMeansOpts)-                    . V.fromList-                    . fmap V.convert-                    . H.toRows-                    . H.fromColumns-                    . fmap H.normalize -- Normalize within eigenvectors (columns).-                    . H.toColumns-                    . H.fromRows+  | otherwise       = kmeansVec k                     . spectralNorm 1 e                     $ mat 

@@ -36,7 +36,6 @@ import qualified Numeric.LinearAlgebra as H import qualified Numeric.LinearAlgebra.Devel as H import qualified Numeric.LinearAlgebra.SVD.SVDLIBC as SVD-import Debug.Trace  -- Local