hierarchical-spectral-clustering 0.4.1.3 → 0.5.0.0

raw patch · 7 files changed

+190/−80 lines, 7 filesdep +mwc-randomPVP ok

version bump matches the API change (PVP)

Dependencies added: mwc-random

API changes (from Hackage documentation)

Files

• app/Main.hs +19/−16
• hierarchical-spectral-clustering.cabal +2/−1
• src/Math/Clustering/Hierarchical/Spectral/Dense.hs +54/−30
• src/Math/Clustering/Hierarchical/Spectral/Sparse.hs +55/−32
• src/Math/Clustering/Hierarchical/Spectral/Test.hs +8/−0
• src/Math/Clustering/Hierarchical/Spectral/Types.hs +7/−1
• src/Math/Clustering/Hierarchical/Spectral/Utility.hs +45/−0

@@ -78,6 +78,8 @@                                        <?> "([1] | INT) Number of eigenvectors to use while clustering with kmeans. Takes from the first eigenvector. Recommended to start at 2 and work up from there if needed."                        , minModularity  :: Maybe Double                                        <?> "([0] | DOUBLE) Minimum modularity to be over to continue recursion."+                       , numRuns        :: Maybe Int+                                       <?> "([Nothing] | INT) Number of runs for permutation test at each split. If present, returns p-value for each split in the tree."                        , eigenGroup     :: Maybe String                                        <?> "([SignGroup] | KMeansGroup) Whether to group the eigenvector using the sign or kmeans while clustering. While the default is sign, kmeans may be more accurate (but starting points are arbitrary)."                        , separateComponents :: Bool@@ -117,6 +119,7 @@         minSize'            = fmap MinSize . unHelpful . minSize $ opts         numEigen'           = fmap NumEigen . unHelpful . numEigen $ opts         minModularity'      = fmap Q . unHelpful . minModularity $ opts+        numRuns'            = unHelpful . numRuns $ opts         eigenGroup'         =           maybe SignGroup (readOrErr "Cannot read --eigen-group")             . unHelpful@@ -143,40 +146,40 @@             Dense -> do                 (items, mat) <- readDenseAdjMatrix decodeOpt stdin -                let cluster items = clusteringTreeToGenericClusteringTree+                let cluster items = fmap clusteringTreeToGenericClusteringTree                                   . HD.hierarchicalSpectralClusterAdj                                       eigenGroup'                                       (fmap unNumEigen numEigen')                                       (fmap unMinSize minSize')                                       minModularity'+                                      numRuns'                                       items -                return $-                    if separateComponents'-                        then Multiple-                           . fmap (uncurry cluster)-                           . getComponentMatsItems items-                           $ mat-                        else Single $ cluster items mat+                if separateComponents'+                    then fmap Multiple+                       . mapM (uncurry cluster)+                       . getComponentMatsItems items+                       $ mat+                    else Single <$> cluster items mat             Sparse -> do                 -- (items, mat) <- readEigenSparseAdjMatrix decodeOpt stdin                 (items, mat) <- readSparseAdjMatrix decodeOpt stdin -                let cluster items = clusteringTreeToGenericClusteringTree+                let cluster items = fmap clusteringTreeToGenericClusteringTree                                   . HS.hierarchicalSpectralClusterAdj                                       eigenGroup'                                       (fmap unNumEigen numEigen')                                       (fmap unMinSize minSize')                                       minModularity'+                                      numRuns'                                       items -                return $-                    if separateComponents'-                        then Multiple-                           . fmap (uncurry cluster)-                           . getComponentMatsItems items-                           $ mat-                        else Single $ cluster items mat+                if separateComponents'+                    then fmap Multiple+                       . mapM (uncurry cluster)+                       . getComponentMatsItems items+                       $ mat+                    else Single <$> cluster items mat              (Premade file) ->                 fmap (either error Single . A.eitherDecode) . B.readFile $ file

@@ -1,6 +1,6 @@ cabal-version: >=1.10 name: hierarchical-spectral-clustering-version: 0.4.1.3+version: 0.5.0.0 license: GPL-3 license-file: LICENSE copyright: 2019 Gregory W. Schwartz@@ -43,6 +43,7 @@         managed >=1.0.6,         modularity >=0.2.1.1,         mtl >=2.2.2,+        mwc-random >=0.13.6.0,         safe >=0.3.17,         sparse-linear-algebra >=0.3.1,         spectral-clustering >=0.3.1.2,

@@ -53,37 +53,48 @@                             -> Maybe NumEigen                             -> Maybe Int                             -> Maybe Q+                            -> Maybe Runs                             -> Items a                             -> Either FeatureMatrix B-                            -> ClusteringTree a-hierarchicalSpectralCluster eigenGroup normFlag numEigenMay minSizeMay minModMay initItems initMat =+                            -> IO (ClusteringTree a)+hierarchicalSpectralCluster eigenGroup normFlag numEigenMay minSizeMay minModMay runsMay initItems initMat =     go initItems initB   where     initB = either (getB normFlag) id $ initMat     minMod   = fromMaybe (Q 0) minModMay     minSize  = fromMaybe 1 minSizeMay     numEigen = fromMaybe 1 numEigenMay-    go :: Items a -> B -> ClusteringTree a+    go :: Items a -> B -> IO (ClusteringTree a)     go !items !b =         if (H.rows $ unB b) > 1             && hasMultipleClusters clusters             && ngMod > minMod             && H.rows (unB left) >= minSize             && H.rows (unB right) >= minSize-            then-                Node { rootLabel = vertex-                     , subForest = [ go (subsetVector items leftIdxs) left-                                   , go (subsetVector items rightIdxs) right-                                   ]-                     }--            else-                Node {rootLabel = vertex, subForest = []}+            then do+                goLeft <- go (subsetVector items leftIdxs) left+                goRight <- go (subsetVector items rightIdxs) right+                vertex <- getVertex+                return $ Node { rootLabel = vertex+                              , subForest = [goLeft, goRight]+                              }+            else do+                vertex <- getVertex+                return $ Node {rootLabel = vertex, subForest = []}       where-        vertex      = ClusteringVertex-                        { _clusteringItems = items-                        , _ngMod = ngMod-                        }+        getVertex = do+          pVal <- mapM (\ x -> permutationTest+                                (\l m -> unQ $ getBModularity l m)+                                x+                                clusters+                                b+                       )+                  runsMay+          return $ ClusteringVertex+                 { _clusteringItems = items+                 , _ngMod = ngMod+                 , _pVal = pVal+                 }         clusters :: H.Vector Double         clusters = spectralClustering eigenGroup b         spectralClustering :: EigenGroup -> B -> H.Vector Double@@ -104,43 +115,56 @@                             -> Maybe NumEigen                             -> Maybe Int                             -> Maybe Q+                            -> Maybe Runs                             -> Items a                             -> AdjacencyMatrix-                            -> ClusteringTree a-hierarchicalSpectralClusterAdj !eigenGroup !numEigenMay !minSizeMay !minModMay !items !adjMat =+                            -> IO (ClusteringTree a)+hierarchicalSpectralClusterAdj !eigenGroup !numEigenMay !minSizeMay !minModMay !runsMay !items !adjMat =     if H.rows adjMat > 1         && hasMultipleClusters clusters         && ngMod > minMod         && H.rows left >= minSize         && H.rows right >= minSize-        then-            Node { rootLabel = vertex-                 , subForest =-                    [ hierarchicalSpectralClusterAdj+        then do+            goLeft <- hierarchicalSpectralClusterAdj                         eigenGroup                         numEigenMay                         minSizeMay                         minModMay+                        runsMay                         (subsetVector items leftIdxs)                         left-                    , hierarchicalSpectralClusterAdj+            goRight <- hierarchicalSpectralClusterAdj                         eigenGroup                         numEigenMay                         minSizeMay                         minModMay+                        runsMay                         (subsetVector items rightIdxs)                         right-                    ]-                 }-        else-            Node {rootLabel = vertex, subForest = []}+            vertex <- getVertex+            return $ Node { rootLabel = vertex , subForest = [goLeft, goRight] }+        else do+            vertex <- getVertex+            return $ Node {rootLabel = vertex, subForest = []}   where     minMod      = fromMaybe (Q 0) minModMay     minSize     = fromMaybe 1 minSizeMay     numEigen    = fromMaybe 1 numEigenMay-    vertex      = ClusteringVertex { _clusteringItems = items-                                   , _ngMod = ngMod-                                   }+    getVertex = do+      pVal <- mapM+                (\ x -> permutationTest+                        (\l m -> unQ $ getModularity l m)+                        x+                        clusters+                        adjMat+                )+                runsMay+      return $ ClusteringVertex+             { _clusteringItems = items+             , _ngMod = ngMod+             , _pVal = pVal+             }     clusters    = spectralClustering eigenGroup adjMat     spectralClustering :: EigenGroup -> AdjacencyMatrix -> H.Vector Double     spectralClustering SignGroup   = spectralClusterNorm

@@ -29,6 +29,7 @@ import qualified Data.Sparse.Common as S import qualified Data.Vector as V import qualified Data.Vector.Storable as VS+import qualified Data.Vector.Unboxed as U import qualified Numeric.LinearAlgebra.Sparse as S  -- Local@@ -55,37 +56,48 @@                             -> Maybe NumEigen                             -> Maybe Int                             -> Maybe Q+                            -> Maybe Runs                             -> Items a                             -> Either FeatureMatrix B-                            -> ClusteringTree a-hierarchicalSpectralCluster eigenGroup normFlag numEigenMay minSizeMay minModMay initItems initMat =+                            -> IO (ClusteringTree a)+hierarchicalSpectralCluster eigenGroup normFlag numEigenMay minSizeMay minModMay runsMay initItems initMat =     go initItems initB   where     initB = either (getB normFlag) id $ initMat     minMod   = fromMaybe (Q 0) minModMay     minSize  = fromMaybe 1 minSizeMay     numEigen = fromMaybe 1 numEigenMay-    go :: Items a -> B -> ClusteringTree a+    go :: Items a -> B -> IO (ClusteringTree a)     go !items !b =         if (S.nrows $ unB b) > 1             && hasMultipleClusters clusters             && ngMod > minMod             && S.nrows (unB left) >= minSize             && S.nrows (unB right) >= minSize-            then-                Node { rootLabel = vertex-                     , subForest = [ go (subsetVector items leftIdxs) left-                                   , go (subsetVector items rightIdxs) right-                                   ]-                     }--            else-                Node {rootLabel = vertex, subForest = []}+            then do+                goLeft <- go (subsetVector items leftIdxs) left+                goRight <- go (subsetVector items rightIdxs) right+                vertex <- getVertex+                return $ Node { rootLabel = vertex+                              , subForest = [goLeft, goRight]+                              }+            else do+                vertex <- getVertex+                return $ Node {rootLabel = vertex, subForest = []}       where-        vertex      = ClusteringVertex-                        { _clusteringItems = items-                        , _ngMod = ngMod-                        }+        getVertex = do+          pVal <- mapM (\ x -> permutationTestSparse+                                (\l m -> unQ $ getBModularity l m)+                                x+                                clusters+                                b+                       )+                  runsMay+          return $ ClusteringVertex+                 { _clusteringItems = items+                 , _ngMod = ngMod+                 , _pVal = pVal+                 }         clusters :: S.SpVector Double         clusters = spectralClustering eigenGroup b         spectralClustering :: EigenGroup -> B -> S.SpVector Double@@ -118,36 +130,47 @@                                -> Maybe NumEigen                                -> Maybe Int                                -> Maybe Q+                               -> Maybe Runs                                -> Items a                                -> AdjacencyMatrix-                               -> ClusteringTree a-hierarchicalSpectralClusterAdj eigenGroup numEigenMay minSizeMay minModMay initItems initMat =+                               -> IO (ClusteringTree a)+hierarchicalSpectralClusterAdj eigenGroup numEigenMay minSizeMay minModMay runsMay initItems initMat =     go initItems initMat   where     minMod      = fromMaybe (Q 0) minModMay     minSize     = fromMaybe 1 minSizeMay     numEigen    = fromMaybe 1 numEigenMay-    go :: Items a -> AdjacencyMatrix -> ClusteringTree a+    go :: Items a -> AdjacencyMatrix -> IO (ClusteringTree a)     go !items !mat =         if S.nrows mat > 1             && hasMultipleClusters clusters             && ngMod > minMod             && S.nrows left >= minSize             && S.nrows right >= minSize-            then-                Node { rootLabel = vertex-                     , subForest = [ go (subsetVector items leftIdxs) left-                                   , go (subsetVector items rightIdxs) right-                                   ]-                     }--            else-                Node {rootLabel = vertex, subForest = []}+            then do+                goLeft <- go (subsetVector items leftIdxs) left+                goRight <- go (subsetVector items rightIdxs) right+                vertex <- getVertex+                return $ Node { rootLabel = vertex+                              , subForest = [ goLeft, goRight]+                              }+            else do+                vertex <- getVertex+                return $ Node {rootLabel = vertex, subForest = []}       where-        vertex      = ClusteringVertex-                        { _clusteringItems = items-                        , _ngMod = ngMod-                        }+        getVertex = do+          pVal <- mapM (\ x -> permutationTestSparse+                                (\l m -> unQ $ getModularity l m)+                                x+                                clusters+                                mat+                       )+                  runsMay+          return $ ClusteringVertex+                 { _clusteringItems = items+                 , _ngMod = ngMod+                 , _pVal = pVal+                 }         clusters :: S.SpVector Double         clusters = spectralClustering eigenGroup mat         spectralClustering :: EigenGroup -> AdjacencyMatrix -> S.SpVector Double

@@ -86,6 +86,7 @@                     Nothing                     Nothing                     Nothing+                    Nothing                     exampleItems                     (Left $ exampleMatrix 3 exampleData) @@ -95,6 +96,7 @@                     (Just 2)                     Nothing                     Nothing+                    Nothing                     exampleItems                     (Left $ exampleMatrix 3 exampleData) @@ -103,6 +105,7 @@                         Nothing                         Nothing                         Nothing+                        Nothing                         exampleItems                         adjacencyExample @@ -111,6 +114,7 @@                         (Just 2)                         Nothing                         Nothing+                        Nothing                         exampleItems                         adjacencyExample @@ -136,6 +140,7 @@                           Nothing                           Nothing                           Nothing+                          Nothing                           exampleItems                           denseAdjacencyExample @@ -144,6 +149,7 @@                           (Just 2)                           Nothing                           Nothing+                          Nothing                           exampleItems                           denseAdjacencyExample @@ -153,6 +159,7 @@                         Nothing                         Nothing                         Nothing+                        Nothing                         exampleItems                         (Left denseFeatureExample) @@ -160,6 +167,7 @@                         KMeansGroup                         True                         (Just 2)+                        Nothing                         Nothing                         Nothing                         exampleItems

@@ -22,6 +22,7 @@     , getClusterItemsGenericTree     , Q (..)     , NumEigen (..)+    , Runs (..)     ) where  -- Remote@@ -42,17 +43,20 @@ type ClusteringTree a = Tree (ClusteringVertex a) type GenericClusteringTree a = Tree (GenericClusteringVertex a) type NumEigen = Int+type Runs = Int  data EigenGroup = SignGroup | KMeansGroup deriving (Read, Show, Generic)  data ClusteringVertex a = ClusteringVertex     { _clusteringItems  :: !(Items a)     , _ngMod            :: !Q+    , _pVal           :: !(Maybe Double)     } deriving (Eq, Ord, Read, Show, Generic)  data GenericClusteringVertex a = GenericClusteringVertex     { _item :: !(Maybe (Items a))     , _distance :: !(Maybe Double)+    , _pValue :: !(Maybe Double)     } deriving (Eq, Ord, Read, Show, Generic)  -- | Convert a ClusteringTree to a Dendrogram. Modularity is the distance.@@ -95,6 +99,7 @@       Node { rootLabel = ( GenericClusteringVertex                              { _item = Just $ _clusteringItems n                              , _distance = Nothing+                             , _pValue = Nothing                              }                          )            , subForest = []@@ -103,6 +108,7 @@       Node { rootLabel = ( GenericClusteringVertex                             { _item = Nothing                             , _distance = Just . unQ . _ngMod $ n+                            , _pValue = _pVal n                             }                          )            , subForest = fmap go xs@@ -115,7 +121,7 @@ -- | Gather clusters (leaves) from the tree. getClusterItemsTree :: ClusteringTree a -> [Items a] getClusterItemsTree = fmap _clusteringItems . leaves-                    + -- | Gather clusters (leaves) from the generic tree. getClusterItemsGenericTree :: GenericClusteringTree a -> [Items a] getClusterItemsGenericTree = catMaybes . fmap _item . leaves

@@ -8,12 +8,20 @@  module Math.Clustering.Hierarchical.Spectral.Utility   ( subsetVector+  , permutationTest+  , permutationTestSparse   ) where  -- Remote+import Data.List (genericLength) import Data.Maybe (fromMaybe) import qualified Data.Foldable as F+import qualified Data.Sparse.Common as S import qualified Data.Vector as V+import qualified Data.Vector.Unboxed as U+import qualified Data.Vector.Generic as G+import qualified System.Random.MWC as R+import qualified System.Random.MWC.Distributions as R  -- Local @@ -26,3 +34,40 @@                       $ xs V.!? i                       ) : acc                     ) []+++-- | Permutation test returning p-value.+permutationTest :: (G.Vector v a)+                => (v a -> b -> Double)+                -> Int+                -> v a+                -> b+                -> IO Double+permutationTest f runs labs samples = do+  let obs = f labs samples+      run r = do+        g <- R.initialize . V.singleton $ fromIntegral r+        ls <- R.uniformShuffle labs g+        return $ f ls samples++  n <- genericLength . filter (>= obs) <$> mapM run [1..runs]++  return $ n / fromIntegral runs++-- | Permutation test returning p-value from sparse vector.+permutationTestSparse :: (S.SpVector Double -> b -> Double)+                -> Int+                -> S.SpVector Double+                -> b+                -> IO Double+permutationTestSparse f runs labsInit samples = do+  let obs = f labsInit samples+      labs = U.fromList $ S.toDenseListSV labsInit+      run r = do+        g <- R.initialize . V.singleton $ fromIntegral r+        ls <- S.vr . U.toList <$> R.uniformShuffle labs g+        return $ f ls samples++  n <- genericLength . filter (>= obs) <$> mapM run [1..runs]++  return $ n / fromIntegral runs