Haskell implementations of Hidden Markov Models & related algorithms from AIMA book
forwardbackward and fixed-lag-smoothing
+490
+20
.gitignore
+20
.gitignore
+5
CHANGELOG.md
+5
CHANGELOG.md
+202
LICENSE
+202
LICENSE
···
1
+
2
+
Apache License
3
+
Version 2.0, January 2004
4
+
http://www.apache.org/licenses/
5
+
6
+
TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION
7
+
8
+
1. Definitions.
9
+
10
+
"License" shall mean the terms and conditions for use, reproduction,
11
+
and distribution as defined by Sections 1 through 9 of this document.
12
+
13
+
"Licensor" shall mean the copyright owner or entity authorized by
14
+
the copyright owner that is granting the License.
15
+
16
+
"Legal Entity" shall mean the union of the acting entity and all
17
+
other entities that control, are controlled by, or are under common
18
+
control with that entity. For the purposes of this definition,
19
+
"control" means (i) the power, direct or indirect, to cause the
20
+
direction or management of such entity, whether by contract or
21
+
otherwise, or (ii) ownership of fifty percent (50%) or more of the
22
+
outstanding shares, or (iii) beneficial ownership of such entity.
23
+
24
+
"You" (or "Your") shall mean an individual or Legal Entity
25
+
exercising permissions granted by this License.
26
+
27
+
"Source" form shall mean the preferred form for making modifications,
28
+
including but not limited to software source code, documentation
29
+
source, and configuration files.
30
+
31
+
"Object" form shall mean any form resulting from mechanical
32
+
transformation or translation of a Source form, including but
33
+
not limited to compiled object code, generated documentation,
34
+
and conversions to other media types.
35
+
36
+
"Work" shall mean the work of authorship, whether in Source or
37
+
Object form, made available under the License, as indicated by a
38
+
copyright notice that is included in or attached to the work
39
+
(an example is provided in the Appendix below).
40
+
41
+
"Derivative Works" shall mean any work, whether in Source or Object
42
+
form, that is based on (or derived from) the Work and for which the
43
+
editorial revisions, annotations, elaborations, or other modifications
44
+
represent, as a whole, an original work of authorship. For the purposes
45
+
of this License, Derivative Works shall not include works that remain
46
+
separable from, or merely link (or bind by name) to the interfaces of,
47
+
the Work and Derivative Works thereof.
48
+
49
+
"Contribution" shall mean any work of authorship, including
50
+
the original version of the Work and any modifications or additions
51
+
to that Work or Derivative Works thereof, that is intentionally
52
+
submitted to Licensor for inclusion in the Work by the copyright owner
53
+
or by an individual or Legal Entity authorized to submit on behalf of
54
+
the copyright owner. For the purposes of this definition, "submitted"
55
+
means any form of electronic, verbal, or written communication sent
56
+
to the Licensor or its representatives, including but not limited to
57
+
communication on electronic mailing lists, source code control systems,
58
+
and issue tracking systems that are managed by, or on behalf of, the
59
+
Licensor for the purpose of discussing and improving the Work, but
60
+
excluding communication that is conspicuously marked or otherwise
61
+
designated in writing by the copyright owner as "Not a Contribution."
62
+
63
+
"Contributor" shall mean Licensor and any individual or Legal Entity
64
+
on behalf of whom a Contribution has been received by Licensor and
65
+
subsequently incorporated within the Work.
66
+
67
+
2. Grant of Copyright License. Subject to the terms and conditions of
68
+
this License, each Contributor hereby grants to You a perpetual,
69
+
worldwide, non-exclusive, no-charge, royalty-free, irrevocable
70
+
copyright license to reproduce, prepare Derivative Works of,
71
+
publicly display, publicly perform, sublicense, and distribute the
72
+
Work and such Derivative Works in Source or Object form.
73
+
74
+
3. Grant of Patent License. Subject to the terms and conditions of
75
+
this License, each Contributor hereby grants to You a perpetual,
76
+
worldwide, non-exclusive, no-charge, royalty-free, irrevocable
77
+
(except as stated in this section) patent license to make, have made,
78
+
use, offer to sell, sell, import, and otherwise transfer the Work,
79
+
where such license applies only to those patent claims licensable
80
+
by such Contributor that are necessarily infringed by their
81
+
Contribution(s) alone or by combination of their Contribution(s)
82
+
with the Work to which such Contribution(s) was submitted. If You
83
+
institute patent litigation against any entity (including a
84
+
cross-claim or counterclaim in a lawsuit) alleging that the Work
85
+
or a Contribution incorporated within the Work constitutes direct
86
+
or contributory patent infringement, then any patent licenses
87
+
granted to You under this License for that Work shall terminate
88
+
as of the date such litigation is filed.
89
+
90
+
4. Redistribution. You may reproduce and distribute copies of the
91
+
Work or Derivative Works thereof in any medium, with or without
92
+
modifications, and in Source or Object form, provided that You
93
+
meet the following conditions:
94
+
95
+
(a) You must give any other recipients of the Work or
96
+
Derivative Works a copy of this License; and
97
+
98
+
(b) You must cause any modified files to carry prominent notices
99
+
stating that You changed the files; and
100
+
101
+
(c) You must retain, in the Source form of any Derivative Works
102
+
that You distribute, all copyright, patent, trademark, and
103
+
attribution notices from the Source form of the Work,
104
+
excluding those notices that do not pertain to any part of
105
+
the Derivative Works; and
106
+
107
+
(d) If the Work includes a "NOTICE" text file as part of its
108
+
distribution, then any Derivative Works that You distribute must
109
+
include a readable copy of the attribution notices contained
110
+
within such NOTICE file, excluding those notices that do not
111
+
pertain to any part of the Derivative Works, in at least one
112
+
of the following places: within a NOTICE text file distributed
113
+
as part of the Derivative Works; within the Source form or
114
+
documentation, if provided along with the Derivative Works; or,
115
+
within a display generated by the Derivative Works, if and
116
+
wherever such third-party notices normally appear. The contents
117
+
of the NOTICE file are for informational purposes only and
118
+
do not modify the License. You may add Your own attribution
119
+
notices within Derivative Works that You distribute, alongside
120
+
or as an addendum to the NOTICE text from the Work, provided
121
+
that such additional attribution notices cannot be construed
122
+
as modifying the License.
123
+
124
+
You may add Your own copyright statement to Your modifications and
125
+
may provide additional or different license terms and conditions
126
+
for use, reproduction, or distribution of Your modifications, or
127
+
for any such Derivative Works as a whole, provided Your use,
128
+
reproduction, and distribution of the Work otherwise complies with
129
+
the conditions stated in this License.
130
+
131
+
5. Submission of Contributions. Unless You explicitly state otherwise,
132
+
any Contribution intentionally submitted for inclusion in the Work
133
+
by You to the Licensor shall be under the terms and conditions of
134
+
this License, without any additional terms or conditions.
135
+
Notwithstanding the above, nothing herein shall supersede or modify
136
+
the terms of any separate license agreement you may have executed
137
+
with Licensor regarding such Contributions.
138
+
139
+
6. Trademarks. This License does not grant permission to use the trade
140
+
names, trademarks, service marks, or product names of the Licensor,
141
+
except as required for reasonable and customary use in describing the
142
+
origin of the Work and reproducing the content of the NOTICE file.
143
+
144
+
7. Disclaimer of Warranty. Unless required by applicable law or
145
+
agreed to in writing, Licensor provides the Work (and each
146
+
Contributor provides its Contributions) on an "AS IS" BASIS,
147
+
WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or
148
+
implied, including, without limitation, any warranties or conditions
149
+
of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A
150
+
PARTICULAR PURPOSE. You are solely responsible for determining the
151
+
appropriateness of using or redistributing the Work and assume any
152
+
risks associated with Your exercise of permissions under this License.
153
+
154
+
8. Limitation of Liability. In no event and under no legal theory,
155
+
whether in tort (including negligence), contract, or otherwise,
156
+
unless required by applicable law (such as deliberate and grossly
157
+
negligent acts) or agreed to in writing, shall any Contributor be
158
+
liable to You for damages, including any direct, indirect, special,
159
+
incidental, or consequential damages of any character arising as a
160
+
result of this License or out of the use or inability to use the
161
+
Work (including but not limited to damages for loss of goodwill,
162
+
work stoppage, computer failure or malfunction, or any and all
163
+
other commercial damages or losses), even if such Contributor
164
+
has been advised of the possibility of such damages.
165
+
166
+
9. Accepting Warranty or Additional Liability. While redistributing
167
+
the Work or Derivative Works thereof, You may choose to offer,
168
+
and charge a fee for, acceptance of support, warranty, indemnity,
169
+
or other liability obligations and/or rights consistent with this
170
+
License. However, in accepting such obligations, You may act only
171
+
on Your own behalf and on Your sole responsibility, not on behalf
172
+
of any other Contributor, and only if You agree to indemnify,
173
+
defend, and hold each Contributor harmless for any liability
174
+
incurred by, or claims asserted against, such Contributor by reason
175
+
of your accepting any such warranty or additional liability.
176
+
177
+
END OF TERMS AND CONDITIONS
178
+
179
+
APPENDIX: How to apply the Apache License to your work.
180
+
181
+
To apply the Apache License to your work, attach the following
182
+
boilerplate notice, with the fields enclosed by brackets "[]"
183
+
replaced with your own identifying information. (Don't include
184
+
the brackets!) The text should be enclosed in the appropriate
185
+
comment syntax for the file format. We also recommend that a
186
+
file or class name and description of purpose be included on the
187
+
same "printed page" as the copyright notice for easier
188
+
identification within third-party archives.
189
+
190
+
Copyright [yyyy] [name of copyright owner]
191
+
192
+
Licensed under the Apache License, Version 2.0 (the "License");
193
+
you may not use this file except in compliance with the License.
194
+
You may obtain a copy of the License at
195
+
196
+
http://www.apache.org/licenses/LICENSE-2.0
197
+
198
+
Unless required by applicable law or agreed to in writing, software
199
+
distributed under the License is distributed on an "AS IS" BASIS,
200
+
WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
201
+
See the License for the specific language governing permissions and
202
+
limitations under the License.
+1
README.md
+1
README.md
···
1
+
Haskell implementation of Hidden Markov Models & related algorithms from AIMA book
+1
cabal.project
+1
cabal.project
···
1
+
packages: ./*.cabal
+8
default.nix
+8
default.nix
+33
hmm-haskell.cabal
+33
hmm-haskell.cabal
···
1
+
cabal-version: 2.4
2
+
-- Initial package description 'hmm-haskell.cabal' generated by 'cabal
3
+
-- init'. For further documentation, see
4
+
-- http://haskell.org/cabal/users-guide/
5
+
6
+
name: hmm-haskell
7
+
version: 0.1.0.0
8
+
-- synopsis:
9
+
-- description:
10
+
-- bug-reports:
11
+
license: Apache-2.0
12
+
license-file: LICENSE
13
+
author: Kaushik Chakraborty
14
+
maintainer: kaushik.chakraborty3@cognizant.com
15
+
-- copyright:
16
+
category: Math
17
+
extra-source-files: CHANGELOG.md, README.md
18
+
19
+
library
20
+
exposed-modules: HMM
21
+
-- other-modules:
22
+
-- other-extensions:
23
+
build-depends: base ^>=4.12.0.0 && < 4.13
24
+
, papa >=0.3.1 && < 1
25
+
, linear >= 1.20.8 && < 1.21
26
+
, vector >= 0.12.0 && < 0.13
27
+
, mtl >= 2.2.2 && < 2.3
28
+
, transformers >= 0.5.6 && < 0.6
29
+
, hmatrix >= 0.19.0 && < 0.20
30
+
hs-source-dirs: src
31
+
default-language: Haskell2010
32
+
default-extensions: NoImplicitPrelude
33
+
ghc-options: -Wall -O2
+12
hmm.nix
+12
hmm.nix
···
1
+
{ mkDerivation, base, hmatrix, linear, mtl, papa, stdenv
2
+
, transformers, vector
3
+
}:
4
+
mkDerivation {
5
+
pname = "hmm-haskell";
6
+
version = "0.1.0.0";
7
+
src = ./.;
8
+
libraryHaskellDepends = [
9
+
base hmatrix linear mtl papa transformers vector
10
+
];
11
+
license = stdenv.lib.licenses.asl20;
12
+
}
+7
shell.nix
+7
shell.nix
···
1
+
{ nixpkgs ? import <nixpkgs> {}, compiler ? "default"}:
2
+
let
3
+
inherit (nixpkgs) pkgs;
4
+
drv = import ./default.nix { inherit nixpkgs compiler; };
5
+
drvWithTools = pkgs.haskell.lib.addBuildDepends drv [ pkgs.cabal-install ];
6
+
in
7
+
if pkgs.lib.inNixShell then drvWithTools.env else drvWithTools
+199
src/hmm.hs
+199
src/hmm.hs
···
1
+
{-# LANGUAGE DataKinds #-}
2
+
{-# LANGUAGE KindSignatures #-}
3
+
{-# LANGUAGE ScopedTypeVariables #-}
4
+
{-# LANGUAGE TemplateHaskell #-}
5
+
{-# LANGUAGE TypeOperators #-}
6
+
module HMM where
7
+
8
+
-- import Prelude (undefined)
9
+
-- import System.IO (print)
10
+
import Text.Printf (printf)
11
+
12
+
import Papa
13
+
import GHC.TypeNats
14
+
import Data.Proxy (Proxy(..))
15
+
16
+
import Control.Monad.State (State, runState)
17
+
import Control.Monad.Trans.Maybe (MaybeT, runMaybeT)
18
+
19
+
import qualified Data.Vector as Vec
20
+
21
+
import Linear.V
22
+
import Linear
23
+
24
+
import qualified Numeric.LinearAlgebra.Static as HM
25
+
import qualified Numeric.LinearAlgebra as HMat
26
+
-----------------------------------------------------------------------------------
27
+
28
+
-- | Type synonyms
29
+
type R = Double
30
+
type TransitionModel (s :: Nat) a = V s (V s a)
31
+
type SensorDiagonal (s :: Nat) a = V s (V s a)
32
+
type SensorModel (t :: Nat) (s :: Nat) a = V t (SensorDiagonal s a)
33
+
type Message (s :: Nat) a = V s (V 1 a)
34
+
type Distribution (s :: Nat) a = V s a
35
+
36
+
-- | Utility Functions
37
+
38
+
-- | Multiply each number by a constant such that the sum is 1.0
39
+
40
+
normalise :: (Fractional a, Foldable f, Functor f) => f a -> f a
41
+
normalise xs | null xs = xs
42
+
| otherwise =
43
+
let
44
+
s = sum xs
45
+
in
46
+
(/ s) <$> xs
47
+
48
+
toHM :: (KnownNat s, KnownNat t) => V s (V t R) -> HM.L s t
49
+
toHM = HM.matrix . foldMap (Vec.toList . toVector)
50
+
51
+
fromHM :: (KnownNat s, KnownNat t) => HM.L s t -> V s (V t R)
52
+
fromHM m = V $ V <$> Vec.fromList (fmap Vec.fromList $ HMat.toLists $ HM.extract m)
53
+
54
+
inverse :: forall s. (KnownNat s) => V s (V s R) -> V s (V s R)
55
+
inverse = fromHM . HM.inv . toHM
56
+
57
+
reverseV :: forall s a. (KnownNat s) => V s a -> V s a
58
+
reverseV = V . Vec.fromList . Papa.reverse . toList
59
+
60
+
extract1 :: (KnownNat s) => V s (V 1 a) -> V s a
61
+
extract1 = V . foldMap toVector
62
+
63
+
infix 8 !**!
64
+
(!**!) :: (KnownNat m, KnownNat n) => V m (V n R) -> V m (V n R) -> V m (V n R)
65
+
as !**! bs = fromHM $ toHM as * toHM bs
66
+
67
+
unitColumn :: forall s. KnownNat s => Message s R
68
+
unitColumn = V $ Vec.replicate (fromIntegral $ natVal (Proxy :: Proxy s)) (toV $ V1 1.0)
69
+
70
+
71
+
-- | Hidden Markov Models
72
+
73
+
data HMM (s :: Nat) (t :: Nat) a b = HMM {
74
+
-- prior distribution which is used as initial forward message
75
+
_prior :: Message s a,
76
+
-- transition model with `s` states as `sxs` matrix
77
+
_tModel :: TransitionModel s a,
78
+
-- evidence value vector, each index maps to the corresponding evidence values
79
+
_sDist :: V t b,
80
+
-- sensor model with `t` evidence values having `sxs` diagonal matrix capturing their ditributions for each state `s`
81
+
_sModel :: SensorModel t s a
82
+
} deriving (Show)
83
+
makeLenses ''HMM
84
+
85
+
86
+
-- | Smart HMM constructors
87
+
88
+
-- | Make HMM
89
+
-- `mp :: Maybe (Message s a)` -- Prior distribution on the initial state, `P(X0)`. If nothing then considered `(0.5,0.5)`
90
+
-- `xs :: TransitionModel s a` -- Transition Model as `sxs` matrix, `s` being the number of states
91
+
-- `evs :: V t b` - evidence values vector map
92
+
-- `es :: V t (V s a)` -- a `txs` matrix for each evidence value `t`, a vector capturing conditional probabilities for each state `s`
93
+
mkHMM_ :: (KnownNat s, KnownNat t) => Maybe (Message s R) -> TransitionModel s R -> V t b -> V t (V s R) -> HMM s t R b
94
+
mkHMM_ mp xs evs = let p = fromMaybe (V $ Vec.replicate (dim xs) 0.5) mp
95
+
in
96
+
HMM p xs evs . (scaled <$>)
97
+
98
+
-- | A default HMM where both transition model states and evidence variables have boolean support
99
+
mkHMM :: KnownNat s => TransitionModel s R -> V 2 (V s R) -> HMM s 2 R Bool
100
+
mkHMM ts = mkHMM_ Nothing ts (V $ Vec.fromList [True, False])
101
+
102
+
-- | A default HMM with one state variable having boolean support
103
+
mkHMM1 :: V2 (V2 R) -> V2 (V2 R) -> HMM 2 2 R Bool
104
+
mkHMM1 ts ss = let f = (_V #) . over mapped toV
105
+
in
106
+
uncurry mkHMM $ over both f (ts , ss)
107
+
108
+
sensorDiagonal :: (KnownNat s, KnownNat t, Eq b) => HMM s t a b -> b -> Maybe (SensorDiagonal s a)
109
+
sensorDiagonal hmm e = findIndexOf (sDist.folded) (== e) hmm >>= (\i -> hmm ^? sModel . ix i)
110
+
111
+
-- | Fixed Lag Smoothing
112
+
113
+
-- | Persistent State
114
+
data Persistent (s :: Nat) a b d= Persistent {
115
+
_t :: d,
116
+
_f_msg :: Message s a,
117
+
_b :: V s (V s a),
118
+
_e_td_t :: Vec.Vector b
119
+
} deriving (Show)
120
+
makeLenses ''Persistent
121
+
122
+
-- | State to start with
123
+
persistentInit :: forall s b d. (KnownNat s, Integral d) => Message s R -> Persistent s R b d
124
+
persistentInit p = Persistent { _t = 1, _f_msg = p, _b = identity, _e_td_t = Vec.empty}
125
+
126
+
-- | Filtering message propagated forward
127
+
forward :: (KnownNat s) => HMM s t R b -> SensorDiagonal s R -> Message s R -> Message s R
128
+
forward hmm ot f = normalise (ot !*! Linear.transpose (hmm ^. tModel) !*! f)
129
+
130
+
-- | Smoothing message propagated backward
131
+
backward :: (KnownNat s) => HMM s t R b -> SensorDiagonal s R -> Message s R -> Message s R
132
+
backward hmm ok1 bk2 = (hmm ^. tModel) !*! ok1 !*! bk2
133
+
134
+
135
+
-- | Online Algorithm for smoothing with a fixed time lag of `d` steps
136
+
-- hmm -- HMM model
137
+
-- d -- length of lag
138
+
-- e -- evidence at time t
139
+
fixedLagSmoothing :: forall s u b d. (KnownNat s, KnownNat u, Eq b, Integral d) => HMM s u R b -> d -> b -> MaybeT (State (Persistent s R b d)) (Distribution s R)
140
+
fixedLagSmoothing hmm d e = do
141
+
e_td_t %= flip Vec.snoc e
142
+
143
+
o_t <- uses e_td_t $ (^?! _Just) . sensorDiagonal hmm . Vec.last
144
+
145
+
t' <- use t
146
+
if t' > d then
147
+
do
148
+
e_td_t %= Vec.drop 1
149
+
o_tmd <- uses e_td_t $ (^?! _Just) . sensorDiagonal hmm . Vec.head
150
+
f_msg %= forward hmm o_tmd
151
+
b %= \b' -> inverse o_tmd !*! inverse (hmm ^. tModel) !*! b' !*! (hmm ^. tModel) !*! o_t
152
+
t += 1
153
+
154
+
(f'', b'') <- liftA2 (,) (use f_msg) (use b)
155
+
156
+
return $ extract1 $ normalise (f'' !**! (b'' !*! unitColumn))
157
+
else
158
+
do
159
+
b %= \b' -> b' !*! (hmm ^. tModel) !*! o_t
160
+
t += 1
161
+
mzero
162
+
163
+
164
+
-- | The forward–backward algorithm for smoothing: computing posterior prob- abilities of a sequence of states given a sequence of observations
165
+
-- hmm - HMM model as a way to implement
166
+
-- evs - list of evidences for each time step
167
+
forwardBackward :: forall s t u b . (KnownNat t, KnownNat s, KnownNat u, Eq b) => HMM s u R b -> V t b -> V t (Distribution s R)
168
+
forwardBackward hmm evs = let
169
+
-- forward messages
170
+
fv :: V (t + 1) (Message s R)
171
+
fv = V $ Vec.fromList $ foldl' (\m@(x:_) e -> forward hmm (sensorDiagonal hmm e ^?! _Just) x : m) [hmm ^. prior] evs
172
+
-- reifying the number of evidences
173
+
tNat = dim evs
174
+
-- getting rid of the prior message from the end of the list
175
+
fv_0 :: V t (Message s R)
176
+
fv_0 = V $ Vec.fromList $ fv ^.. taking tNat traversed
177
+
-- backward messages
178
+
bs :: V t (Message s R)
179
+
bs = V $ Vec.fromList $ foldl' (\m@(x:_) e -> backward hmm (sensorDiagonal hmm e ^?! _Just) x : m) [unitColumn] $ reverseV evs
180
+
in
181
+
-- smoothing probabilities in reverse order of the list of evidences
182
+
liftA2 (\f b' -> extract1 $ normalise $ f !**! b') fv_0 (reverseV bs)
183
+
184
+
-- | execution
185
+
umbrellaHMM :: HMM 2 2 R Bool
186
+
umbrellaHMM = mkHMM1 (V2 (V2 0.7 0.3) (V2 0.3 0.7)) (V2 (V2 0.9 0.2) (V2 0.1 0.8))
187
+
188
+
runFLSAlgo :: [Bool] -> [Maybe (Distribution 2 String)]
189
+
runFLSAlgo bs = let hmm = mkHMM1 (V2 (V2 0.7 0.3) (V2 0.3 0.7)) (V2 (V2 0.9 0.2) (V2 0.1 0.8))
190
+
initState = persistentInit (V $ Vec.fromList [V $ Vec.singleton 0.5, V $ Vec.singleton 0.5]) :: Persistent 2 R Bool Integer
191
+
algo = fixedLagSmoothing hmm 1
192
+
a :: ([Maybe (Distribution 2 R)], Persistent 2 R Bool Integer)
193
+
a = foldl' (\(rs, s) x -> runState (runMaybeT $ algo x) s & _1 #%~ ((rs ++) . pure)) ([], initState) bs
194
+
in
195
+
(a ^. _1) & traverse.traverse.traverse %~ \(x :: R) -> printf "%.3f" x :: String
196
+
197
+
198
+
runFBAlgo :: (KnownNat t) =>V t Bool -> [Distribution 2 R]
199
+
runFBAlgo bs = toList $ forwardBackward umbrellaHMM bs