beginning markov code

constants.py

@@ -1,5 +1,6 @@
 from dataclasses import dataclass
 import numpy as np
+from math import sqrt
 
 
 @dataclass(frozen=True)  # implements constructor among other boilerplate
@@ -10,6 +11,10 @@ class State:
     entry: float  # pi
     exit: float  # eta
 
+    @property
+    def std_dev(self):
+        return sqrt(self.variance)
+
 
 state1 = State(1, 1.44, 0.44, 0.02)
 state2 = State(4, 0.49, 0.56, 0.03)

main.py

@@ -1,3 +0,0 @@
-from constants import *
-
-print(f"state 1: {state1.mean}")

markov.py

@@ -0,0 +1,39 @@
+from dataclasses import dataclass, field
+from typing import List
+import numpy as np
+
+from maths import gaussian
+
+@dataclass
+class Likelihood:
+    forward: float  # forward likelihood
+    backward: float  # backward likelihood
+
+@dataclass
+class TimeStep:
+    states: List[Likelihood] = field(default_factory=list)
+
+@dataclass
+class Transition:
+    pass
+
+class MarkovModel:
+
+    def __init__(self, states: list, observations: list = list(), state_transitions: list = list()):
+        self.observations = observations
+        self.state_transitions = state_transitions
+
+        self.states = states # number of states
+        # self.timesteps = list()
+
+        self.forward = np.zeros((len(states), len(observations)))
+        self.backward = np.zeros((len(states), len(observations)))
+    
+    def populate_forward(self):
+        for t, observation in enumerate(self.observations): # iterate through observations (time)
+            for state_number, state in enumerate(self.states):
+
+                if t == 0: # calcualte initial
+                    self.forward[state_number, t] = self.state_transitions[0, state_number + 1] * gaussian(observation, state.mean, state.std_dev)
+                else:
+                    self.forward[state_number, t] = gaussian(observation, state.mean, state.std_dev)

maths.py

@@ -0,0 +1,15 @@
+from math import sqrt, pi
+
+import numpy as np
+from numpy import exp
+
+root_2_pi = sqrt(2. * pi) # square root is expensive, define as constant here
+
+def gaussian(x: float, mu: float, sd: float):
+    mu_pert = x - mu # mean pertubation
+
+    coefficient = 1. / (sd * root_2_pi)
+
+    return coefficient * exp( - (mu_pert**2)
+                                        /
+                                    (2.*sd**2))