packages feed

learning-hmm 0.3.2.0 → 0.3.2.1

raw patch · 6 files changed

+114/−23 lines, 6 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

CHANGES.md view
@@ -1,6 +1,9 @@ Revision history for Haskell package learning-hmm === +## Version 0.3.2.1+- Bug fix release+ ## Version 0.3.2.0 - Add function `euclideanDistance` which measures the Euclidean distance between   two models
learning-hmm.cabal view
@@ -1,5 +1,5 @@ name:                learning-hmm-version:             0.3.2.0+version:             0.3.2.1 stability:           experimental  synopsis:            Yet another library for hidden Markov models
src/Learning/HMM.hs view
@@ -172,8 +172,8 @@ --   'outputs' of the model. checkDataIn :: Eq o => String -> HMM s o -> [o] -> () checkDataIn fun HMM {..} xs-  | all (`elem` outputs) xs = ()-  | otherwise               = errorIn fun "illegal data"+  | any (`notElem` outputs) xs = errorIn fun "illegal data"+  | otherwise                  = ()  -- | Convert internal 'HMM' to 'HMM'. fromInternal :: (Eq s, Eq o) => [s] -> [o] -> I.HMM -> HMM s o
src/Learning/HMM/Internal.hs view
@@ -28,7 +28,7 @@ import qualified Data.Vector.Mutable              as MV  ( unsafeNew, unsafeRead, unsafeWrite ) import qualified Data.Vector.Unboxed              as U   ( Vector, fromList, length, map, sum, unsafeFreeze, unsafeIndex, unsafeTail, zip ) import qualified Data.Vector.Unboxed.Mutable      as MU  ( unsafeNew, unsafeRead, unsafeWrite )-import qualified Numeric.LinearAlgebra.Data       as H   ( (!), Matrix, Vector, diag, fromColumns, fromList, fromLists, fromRows, konst, maxElement, maxIndex, toColumns, tr )+import qualified Numeric.LinearAlgebra.Data       as H   ( (!), Matrix, Vector, diag, fromColumns, fromList, fromLists, fromRows, ident, konst, maxElement, maxIndex, toColumns, tr ) import qualified Numeric.LinearAlgebra.HMatrix    as H   ( (<>), (#>), sumElements ) import           Prelude                          hiding ( init ) @@ -60,10 +60,12 @@   phi <- H.fromLists <$> replicateM k (stdSimplex (l-1))   return HMM { nStates          = k              , nOutputs         = l-             , initialStateDist = pi0-             , transitionDist   = w-             , emissionDistT    = H.tr phi+             , initialStateDist = q_ H.#> pi0+             , transitionDist   = w H.<> q_+             , emissionDistT    = q_ H.<> H.tr phi              }+  where+    q_ = q k -- Error matrix  withEmission :: HMM -> U.Vector Int -> HMM withEmission (model @ HMM {..}) xs = model'@@ -167,17 +169,34 @@     -- posterior distribution, i.e., gamma and xi values.     (gammas, xis) = posterior model n xs alphas betas cs +    -- Error matrix+    q_ = q nStates+     -- Using the gamma and xi values, we obtain the optimal initial state     -- probability vector, transition probability matrix, and emission     -- probability matrix.-    pi0  = V.unsafeIndex gammas 0+    pi0  = let g0  = V.unsafeIndex gammas 0+               g0_ = g0 / H.konst (H.sumElements g0) nStates+           in q_ H.#> g0_     w    = let ds = V.foldl1' (+) xis         -- denominators                ns = ds H.#> H.konst 1 nStates -- numerators-           in H.diag (H.konst 1 nStates / ns) H.<> ds+               w_ = H.diag (H.konst 1 nStates / ns) H.<> ds+           in w_ H.<> q_+           {- in H.fromRows $ zipWith3 (\n_ t t0 -> if n_ > eps then t else t0)+            -                          (H.toList ns)+            -                          (H.toRows $ w_ H.<> q_)+            -                          (H.toRows transitionDist)+            -}     phi' = let gs' o = V.map snd $ V.filter ((== o) . fst) $ V.zip (G.convert xs) gammas-               ds    = V.foldl' (+) 0 . gs'  -- denominators+               ds    = V.foldl' (+) (H.konst 0 nStates) . gs'  -- denominators                ns    = V.foldl1' (+) gammas -- numerators-           in H.fromRows $ map (\o -> ds o / ns) [0..(nOutputs - 1)]+               phi_  = H.fromRows $ map (\o -> ds o / ns) [0..(nOutputs - 1)]+           in q_ H.<> phi_+           {- in H.fromColumns $ zipWith3 (\n_ e e0 -> if n_ > eps then e else e0)+            -                             (H.toList ns)+            -                             (H.toColumns $ q_ H.<> phi_)+            -                             (H.toColumns emissionDistT)+            -}      -- We finally obtain the new model and the likelihood for the old model.     model' = model { initialStateDist = pi0@@ -235,3 +254,26 @@                alphas betas (G.convert cs)     xis    = V.zipWith3 (\a b x -> H.diag a H.<> transitionDist H.<> H.diag (b * (emissionDistT H.! x)))                alphas (V.unsafeTail betas) (G.convert $ U.unsafeTail xs)++-- | Global error threshold.+{-# INLINE eps #-}+eps :: Double+eps = 1e-4++-- | Error matrix @q k@ is required to guarantee that the elements of initial+--   states vector and emission/transition matrix are all larger than zero.+--   @k@ is assumed to be the number of states. @q k@ is given by+--       [    1 - eps, (1/k') eps,        ..., (1/k') eps ]+--       [ (1/k') eps,    1 - eps,        ..., (1/k') eps ]+--       [                    ...                         ]+--       [ (1/k') eps,        ..., (1/k') eps,    1 - eps ],+--   where the diagonal elements are @1 - eps@ and the remains are @(1/k')+--   eps@. Here @eps@ is a small error value (given by @1e-4@) and+--   @k' = k - 1@.+q :: Int -> H.Matrix Double+{-# INLINE q #-}+q k = H.konst (1 - eps) (k, k) * e + H.konst (eps / k') (k, k) * (one - e)+  where+    e   = H.ident k+    one = H.konst 1 (k, k)+    k'  = fromIntegral (k - 1)
src/Learning/IOHMM.hs view
@@ -198,8 +198,10 @@ --   in the 'inputs' ('outputs') of the model. checkDataIn :: (Eq i, Eq o) => String -> IOHMM i s o -> [i] -> [o] -> () checkDataIn fun IOHMM {..} xs ys-  | all (`elem` inputs) xs && all (`elem` outputs) ys = ()-  | otherwise                                         = errorIn fun "illegal data"+  | any (`notElem` inputs) xs  = errorIn fun "illegal input data"+  | any (`notElem` outputs) ys = errorIn fun "illegal output data"+  | any (`notElem` xs) inputs  = errorIn fun "insufficient input data"+  | otherwise                  = ()  -- | Convert internal 'IOHMM' to 'IOHMM'. fromInternal :: (Eq i, Eq s, Eq o) => [i] -> [s] -> [o] -> I.IOHMM -> IOHMM i s o
src/Learning/IOHMM/Internal.hs view
@@ -28,7 +28,7 @@ import qualified Data.Vector.Mutable              as MV  ( unsafeNew, unsafeRead, unsafeWrite ) import qualified Data.Vector.Unboxed              as U   ( Vector, fromList, length, map, sum, unsafeFreeze, unsafeIndex, unsafeTail, unzip, zip ) import qualified Data.Vector.Unboxed.Mutable      as MU  ( unsafeNew, unsafeRead, unsafeWrite )-import qualified Numeric.LinearAlgebra.Data       as H   ( (!), Matrix, Vector, diag, fromColumns, fromList, fromLists, fromRows, konst, maxElement, maxIndex, toColumns, tr )+import qualified Numeric.LinearAlgebra.Data       as H   ( (!), Matrix, Vector, diag, fromColumns, fromList, fromLists, fromRows, ident, konst, maxElement, maxIndex, toColumns, tr ) import qualified Numeric.LinearAlgebra.HMatrix    as H   ( (<>), (#>), sumElements ) import           Prelude                          hiding ( init ) @@ -63,10 +63,12 @@   return IOHMM { nInputs          = m                , nStates          = k                , nOutputs         = l-               , initialStateDist = pi0-               , transitionDist   = w-               , emissionDistT    = H.tr phi+               , initialStateDist = q_ H.#> pi0+               , transitionDist   = V.map (H.<> q_) w+               , emissionDistT    = q_ H.<> H.tr phi                }+  where+    q_ = q k -- Error matrix  withEmission :: IOHMM -> U.Vector (Int, Int) -> IOHMM withEmission (model @ IOHMM {..}) xys = model'@@ -175,18 +177,37 @@     -- posterior distribution, i.e., gamma and xi values.     (gammas, xis) = posterior model n xys alphas betas cs +    -- Error matrix+    q_ = q nStates+     -- Using the gamma and xi values, we obtain the optimal initial state     -- probability vector, transition probability matrix, and emission     -- probability matrix.-    pi0  = V.unsafeIndex gammas 0+    -- Each simplex in pi0, w, and phi' remains old if their numerators are+    -- zero.+    pi0  = let g0  = V.unsafeIndex gammas 0+               g0_ = g0 / H.konst (H.sumElements g0) nStates+           in q_ H.#> g0_     w    = let xis' i = V.map snd $ V.filter ((== i) . fst) $ V.zip (G.convert $ U.unsafeTail xs) xis-               ds     = V.foldl1' (+) . xis'  -- denominators+               ds     = V.foldl1' (+) . xis'        -- denominators                ns i   = ds i H.#> H.konst 1 nStates -- numerators-           in V.map (\i -> H.diag (H.konst 1 nStates / ns i) H.<> ds i) (V.generate nInputs id)+               w_ i   = H.diag (H.konst 1 nStates / ns i) H.<> ds i+           in flip V.map (V.generate nInputs id) $ \i -> w_ i H.<> q_+                {- H.fromRows $ zipWith3 (\n_ t t0 -> if n_ > eps then t else t0)+                 -                       (H.toList $ ns i)+                 -                       (H.toRows $ w_ i H.<> q_)+                 -                       (H.toRows $ V.unsafeIndex transitionDist i)+                 -}     phi' = let gs' o = V.map snd $ V.filter ((== o) . fst) $ V.zip (G.convert ys) gammas-               ds    = V.foldl' (+) 0 . gs'  -- denominators-               ns    = V.foldl1' (+) gammas -- numerators-           in H.fromRows $ map (\o -> ds o / ns) [0..(nOutputs - 1)]+               ds    = V.foldl' (+) (H.konst 0 nStates) . gs' -- denominators+               ns    = V.foldl1' (+) gammas                   -- numerators+               phi_  = H.fromRows $ map (\o -> ds o / ns) [0..(nOutputs - 1)]+           in q_ H.<> phi_+           {- in H.fromColumns $ zipWith3 (\n_ e e0 -> if n_ > eps then e else e0)+            -                             (H.toList ns)+            -                             (H.toColumns $ q_ H.<> phi_)+            -                             (H.toColumns emissionDistT)+            -}      -- We finally obtain the new model and the likelihood for the old model.     model' = model { initialStateDist = pi0@@ -247,3 +268,26 @@     xis    = V.zipWith3 (\a b (x, y) -> H.diag a H.<> w x H.<> H.diag (b * (emissionDistT H.! y)))                alphas (V.unsafeTail betas) (G.convert $ U.unsafeTail xys)     w = V.unsafeIndex transitionDist++-- | Global error threshold.+{-# INLINE eps #-}+eps :: Double+eps = 1e-4++-- | Error matrix @q k@ is required to guarantee that the elements of initial+--   states vector and emission/transition matrix are all larger than zero.+--   @k@ is assumed to be the number of states. @q k@ is given by+--       [    1 - eps, (1/k') eps,        ..., (1/k') eps ]+--       [ (1/k') eps,    1 - eps,        ..., (1/k') eps ]+--       [                    ...                         ]+--       [ (1/k') eps,        ..., (1/k') eps,    1 - eps ],+--   where the diagonal elements are @1 - eps@ and the remains are @(1/k')+--   eps@. Here @eps@ is a small error value (given by @1e-4@) and+--   @k' = k - 1@.+q :: Int -> H.Matrix Double+{-# INLINE q #-}+q k = H.konst (1 - eps) (k, k) * e + H.konst (eps / k') (k, k) * (one - e)+  where+    e   = H.ident k+    one = H.konst 1 (k, k)+    k'  = fromIntegral (k - 1)