Hidden Markov source
A hidden Markov source, whose symbols are emitted at random from the states of a hidden Markov chain.
- Name:
hidden_markov_source - Field:
mathematics/stochastic_processes - Formalism: discrete map
- References: Lindgren 2024
Factory
@register_model(
field='mathematics/stochastic_processes',
tags=('hidden markov model', 'information theory', 'markov chain', 'discrete map'),
references=('Lindgren 2024',),
defaults={'until': 300, 'dt': 1},
)
def hidden_markov_source(transition: Matrix = CHAIN, emission: Matrix = EMISSION) -> System:
"""Return a hidden Markov source, whose symbols are emitted at random from the states of a hidden Markov chain.
The hidden state follows the transition matrix `P`, and in hidden state `i` the source emits symbol `j` with
probability `E_ij`. The frequency of symbol `j` is `Σ_i π_i E_ij` for the stationary distribution `π` of the hidden
chain, and the symbols alone are not a Markov chain of any finite order in general.
Args:
transition:
The transition matrix `P` of the hidden chain.
emission:
The emission matrix `E`, one row per hidden state and one column per symbol.
Returns:
system:
The model, with the hidden states `source.hidden` and the symbols `source.symbol`.
Examples:
>>> from skcomplex.models.mathematics.stochastic_processes import hidden_markov_source
>>> hidden_markov_source().interactions[0].observed
'symbol'
"""
chain = TransitionInteraction(
'chain',
target='source',
state='hidden',
matrix=lambda s: transition,
emission=lambda s: emission,
observed='symbol',
)
return System(elements=[HiddenMarkovSourceElement('source')], interactions=[chain])
Components
class HiddenMarkovSourceElement(Element):
"""Source of symbols emitted from the states of a hidden Markov chain, one symbol per step.
Args:
label:
The label of the source.
visual:
How the element is drawn, or `None` for `default_visual_`.
hidden:
The initial hidden state.
symbol:
The initial symbol.
Examples:
>>> from skcomplex.models.mathematics.stochastic_processes import HiddenMarkovSourceElement
>>> HiddenMarkovSourceElement('source', hidden=1).hidden
1
"""
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
hidden: Annotated[float, StateSpec(domain=NN)] = 0,
symbol: Annotated[float, StateSpec(domain=NN)] = 0,
) -> None:
super().__init__(label, visual=visual)
self.hidden = hidden
self.symbol = symbol
class TransitionInteraction(Interaction):
"""Move of a discrete state of the target from one value to another, or between all its values by a matrix.
A member whose state holds `source` moves to `destination` at a `rate`, as a jump in continuous time, with a
`probability` in every step of discrete time, or in the step in which a `condition` holds. With `over='neighbors'`
the rate and the probability are those of one neighbour, written with `namespace.target.other` for the neighbour:
the rates of the neighbours add up, and their probabilities combine as `1 - prod(1 - p)`, the probability that at
least one of them moves the member. With a matrix instead, the state takes the values `0` to `n - 1` and steps by
the rows of a transition matrix in discrete time, the row of its current value. An `emission` matrix then draws the
`observed` state of the target from the row that the next value of the state selects, as a hidden Markov source.
With a residence time `after`, a member moves by a probability or a condition only once it has held the source
value for `after` steps, at the earliest at the end of the step that completes them.
In the `'symbolic'` form, every function reads `namespace.target`, the states and parameters of the target, and the
time `namespace.t`. In the `'agent'` form, the target is an agent group and every function takes the agent and
the model view, `f(agent, model)`, and returns a number, a truth value, or the rows of a matrix, and a rate `r`
moves an agent with probability `1 - exp(-r)` in a step. With `over='cell'`, every member acts on the other members
of its cell: it moves each one that holds the source value with the probability, the rate, or the condition
`f(agent, other, model)` of that member `agent` and the acting member `other`, each by its own draw.
The `schedule` of the agent form orders the members. Under `'synchronous'` every member decides from the states at
the start of the step, under `'sequential'` the members act one after the other in their order, and under
`'random-sequential'` in a new random order every step, each reading the values the members before it set. Under
`'member'` every member acts in its own activation, together with the rules of its element and the other member
interactions, in the order the system lists them, so one rule of a member runs its movement, transfers,
transitions, and life cycle.
Args:
label:
The label of the interaction, unique among the parts of its system.
target:
The path of the element or group whose state moves.
state:
The name of the discrete state of the target.
source:
The value the state moves from, a collection of such values, or a function of the value of the state that
returns whether the state moves from it, such as `lambda value: value >= 0`, a condition in the symbolic
form. A state that moves from several values moves to the destination from each.
destination:
The value the state moves to.
rate:
The function of the rate of the move, or `None`.
probability:
The function of the probability of the move in a step, or `None`.
condition:
The function of the condition under which the state moves in a step, or `None`.
matrix:
The function of the transition matrix, one row and one column per value of the state, or `None`.
emission:
With a matrix, the function of the emission matrix, one row per value of the state and one column per value
of the observed state, or `None`.
observed:
The name of the state of the target that the emission matrix draws, or `None`.
over:
`'member'` for a rate or a probability of the member, `'neighbors'` for that of one neighbour, or `'cell'`
for the agent form whose members act on the members of their cell.
neighborhood:
The neighbourhood of `over='neighbors'`, or `None` for the neighbours of the network or the lattice.
after:
The number of steps a member holds the source value before it can move, or `None`.
schedule:
The order of the members of the agent form, `'synchronous'`, `'sequential'`, `'random-sequential'`, or
`'member'`.
form:
`'symbolic'` for laws of the states of the target, or `'agent'` for a rule over an agent group.
phase:
The phase of the step the agent form runs in, or `None` for the default phase.
Raises:
ValidationError:
The transition has not exactly one kind, the source and the destination are equal, an emission has no
matrix or no observed state, `over` is unknown or goes with a form or a kind it cannot take, the
residence time is negative or goes with a rate, a matrix, or members acting on their cell, the form or the
schedule is unknown, the symbolic form has a schedule, or a function is not callable.
Examples:
>>> from skcomplex.core.interactions import TransitionInteraction
>>> recovery = TransitionInteraction(
... 'recovery', target='hosts', state='state', source=1, destination=0, rate=lambda s: 0.5
... )
>>> recovery.check().over
'member'
>>> agents = TransitionInteraction(
... 'recovery', target='hosts', state='state', source=1, destination=2, probability=lambda a, m: 0.1,
... after=3, form='agent'
... )
>>> type(agents.check().prepare().run(SimpleNamespace()).laws_[0]).__name__
'RuleLaw'
>>> separation = TransitionInteraction(
... 'separation', target='workers', state='employer', source=lambda value: value >= 0, destination=-1,
... probability=lambda a, m: 0.05, form='agent'
... )
>>> separation.check().destination
-1
"""
def __init__(
self,
label: str,
*,
target: Annotated[str, RoleSpec()],
state: str = 'state',
source: SourceValues = 0.0,
destination: float = 1.0,
rate: Callable[..., SupportsFloat] | None = None,
probability: Callable[..., SupportsFloat] | None = None,
condition: Callable[..., object] | None = None,
matrix: Callable[..., Sequence[Sequence[object]]] | None = None,
emission: Callable[..., Sequence[Sequence[object]]] | None = None,
observed: str | None = None,
over: str = 'member',
neighborhood: BaseNeighborhood | None = None,
after: float | None = None,
schedule: str = 'synchronous',
form: str = 'symbolic',
phase: str | None = None,
) -> None:
super().__init__(label)
self.target = target
self.state = state
self.source = source
self.destination = destination
self.rate = rate
self.probability = probability
self.condition = condition
self.matrix = matrix
self.emission = emission
self.observed = observed
self.over = over
self.neighborhood = neighborhood
self.after = after
self.schedule = schedule
self.form = form
self.phase = phase
def _check_kind(self) -> str:
"""Raise an error for a transition without exactly one callable kind, or with an emission it cannot take."""
kinds = [kind for kind in KINDS if getattr(self, kind) is not None]
if len(kinds) != 1 or not callable(getattr(self, kinds[0])):
error_msg = (
f'The transition {self.label!r} has the forms {kinds}. Give exactly one of {", ".join(KINDS)}, as a '
'function of the namespace.'
)
raise ValidationError(error_msg)
unpaired = (self.emission is None) != (self.observed is None)
if unpaired or (self.emission is not None and self.matrix is None):
error_msg = (
f'The transition {self.label!r} has the emission {self.emission!r} of the state {self.observed!r}. '
'Give both with a matrix, or neither.'
)
raise ValidationError(error_msg)
return kinds[0]
def _check_form(self) -> None:
"""Raise an error for a transition without exactly one callable kind, or with arguments its kind ignores."""
kind = self._check_kind()
if self.over not in OVERS:
error_msg = f'The transition {self.label!r} runs over {self.over!r}. Use one of {", ".join(OVERS)}.'
raise ValidationError(error_msg)
if self.form not in FORMS:
error_msg = f'The form of transition {self.label!r} is {self.form!r}. Use one of {", ".join(FORMS)}.'
raise ValidationError(error_msg)
if self.schedule not in SCHEDULES or (self.form == 'symbolic' and self.schedule != 'synchronous'):
error_msg = (
f'The schedule of transition {self.label!r} is {self.schedule!r}. Use one of {", ".join(SCHEDULES)}, '
'and the agent form for another than `synchronous`.'
)
raise ValidationError(error_msg)
if self.over == 'cell' and (self.form != 'agent' or kind == 'matrix' or self.after is not None):
error_msg = (
f'The members of transition {self.label!r} act on their cell, which only the agent form does, by a '
'rate, a probability, or a condition and without a residence time.'
)
raise ValidationError(error_msg)
if self.over == 'neighbors' and (kind not in {'rate', 'probability'} or self.form == 'agent'):
error_msg = (
f'The transition {self.label!r} runs over neighbours, which only a rate or a probability of the '
'symbolic form does. In the agent form, read the neighbours of the agent in the function.'
)
raise ValidationError(error_msg)
single = not (callable(self.source) or isinstance(self.source, Collection))
if self.matrix is None and single and self.source == self.destination:
error_msg = f'The transition {self.label!r} moves from {self.source} to itself. Give distinct values.'
raise ValidationError(error_msg)
if self.after is not None and (kind not in {'probability', 'condition'} or not self.after >= 0):
error_msg = (
f'The transition {self.label!r} has the residence time {self.after!r}. Give a number of steps at least '
'0 with a probability or a condition.'
)
raise ValidationError(error_msg)
def _check_transition(self) -> None:
"""Raise an error for arguments the transition cannot take."""
self._check_form()
def _check(self) -> None:
"""Raise an error for a transition without exactly one kind, or with arguments its kind ignores.
Raises:
ValidationError:
The transition has not exactly one kind, the source and the destination are equal, an emission has
no matrix or no observed state, `over` is unknown or goes with a condition, a matrix, or the agent
form, the residence time is negative or goes with a rate or a matrix, the form is unknown, or a
function is not callable.
"""
super()._check()
self._check_transition()
def _get_moving(self, value: Expression, source: SourceValues) -> Expression:
"""Return 1 for a member whose state holds a source value, and 0 otherwise."""
if callable(source):
return where(SR(source(value)), 1, 0)
if isinstance(source, Collection):
return sum((where(value == item, 1, 0) for item in set(source)), SR(0))
return where(value == source, 1, 0)
def _get_jump(self, value: Expression, source: SourceValues, destination: float, rate: Expression) -> JumpLaw:
"""Return the jump of the state from a source value to the destination at a rate."""
if callable(source) or isinstance(source, Collection):
return JumpLaw(rate * self._get_moving(value, source), {value: destination - value})
return JumpLaw(rate * self._get_moving(value, source), {value: destination - source})
def _get_step(self, value: Expression, source: SourceValues, destination: float, moves: Expression) -> UpdateLaw:
"""Return the update of the state that moves a member from a source value where a condition holds."""
if callable(source) or isinstance(source, Collection):
return UpdateLaw(value, where(self._get_moving(value, source) > 0, where(moves, destination, value), value))
return UpdateLaw(value, where(value == source, where(moves, destination, source), value), source=source)
def _holds_source(self, value: float) -> bool:
"""Return whether the value of the state of an agent is a source value."""
if callable(self.source):
return bool(self.source(value))
if isinstance(self.source, Collection):
return value in set(self.source)
return value == self.source
def _combine_rates(self, rate: Expression) -> Expression:
"""Return the rate of the member, or the sum of the rates of its neighbours."""
return neighbors(rate, over=self.neighborhood) if self.over == 'neighbors' else rate
def _combine_probabilities(self, probability: Expression) -> Expression:
"""Return the probability of the member, or the probability that at least one neighbour moves it."""
if self.over == 'member':
return probability
return 1 - exp(neighbors(log(1 - probability), over=self.neighborhood))
def _get_draw(self) -> Expression:
"""Return the uniform draw of the transition in a step."""
return uniform(re.sub(r'\W', '_', f'transition_{self.label}'))
def _get_chain(self, namespace: SimpleNamespace, value: Expression) -> list[BasePart]:
"""Return the transition law of the state by its matrix, with the emission law of the observed state."""
laws: list[BasePart] = []
if self.matrix is not None:
laws.append(TransitionLaw(value, self.matrix(namespace), label=self.label))
if self.emission is not None and self.observed is not None:
observed = getattr(self._get_member(namespace, 'target'), self.observed)
laws.append(EmissionLaw(observed, value, self.emission(namespace), label=self.label))
return laws
def _get_held(self, agent: AgentViewProtocol, model: ModelViewProtocol) -> int:
"""Return the number of whole steps an agent has held its value, counted in the storage of the model."""
value = float(getattr(agent, self.state))
record = model.storage.get(f'{self.label}.{agent.id}')
held = record[1] + 1 if isinstance(record, tuple) and record[0] == value else 0
model.storage[f'{self.label}.{agent.id}'] = (value, held)
return held
def _draw_row(self, row: Sequence[object], model: ModelViewProtocol) -> float:
"""Return a value drawn by the probabilities of a row of a matrix."""
weights = np.asarray(row, dtype=float)
return float(model.rng.choice(weights.size, p=weights / weights.sum()))
def _evaluate[T](
self,
function: Callable[..., T],
agent: AgentViewProtocol,
other: AgentViewProtocol | None,
model: ModelViewProtocol,
) -> T:
"""Return a function of an agent, or of an agent and the member of its cell that acts on it."""
return function(agent, model) if other is None else function(agent, other, model)
def _moves_agent(
self,
agent: AgentViewProtocol,
model: ModelViewProtocol,
other: AgentViewProtocol | None = None,
) -> bool:
"""Return whether an agent that holds the source value moves in this step, moved by another one or alone."""
if self.condition is not None:
return bool(self._evaluate(self.condition, agent, other, model))
if self.rate is not None:
probability = 1.0 - math.exp(-float(self._evaluate(self.rate, agent, other, model)))
elif self.probability is not None:
probability = float(self._evaluate(self.probability, agent, other, model))
else:
probability = 0.0
return bool(model.rng.random() < probability)
def _get_next(self, agent: AgentViewProtocol, model: ModelViewProtocol) -> float | None:
"""Return the next value of the state of an agent, or `None` when it keeps its value."""
value = float(getattr(agent, self.state))
if self.matrix is not None:
return self._draw_row(self.matrix(agent, model)[int(value)], model)
held = self._get_held(agent, model) if self.after is not None else 0
if not self._holds_source(value) or (self.after is not None and held + 1 < self.after):
return None
return float(self.destination) if self._moves_agent(agent, model) else None
def _get_moves(self, agent: AgentViewProtocol, model: ModelViewProtocol) -> list[tuple[AgentViewProtocol, float]]:
"""Return the members an agent moves in this step, itself or the members of its cell, with their values."""
if self.over != 'cell':
value = self._get_next(agent, model)
return [] if value is None else [(agent, value)]
return [
(other, float(self.destination))
for other in agent.agents_here(agent.group)
if other is not agent
and self._holds_source(float(getattr(other, self.state)))
and self._moves_agent(other, model, agent)
]
def _apply(self, moves: list[tuple[AgentViewProtocol, float]], model: ModelViewProtocol) -> None:
"""Set the moved states, and draw the observed state of an emission."""
for agent, value in moves:
setattr(agent, self.state, value)
if self.emission is not None and self.observed is not None:
setattr(agent, self.observed, self._draw_row(self.emission(agent, model)[int(value)], model))
def _move_member(self, agent: AgentViewProtocol, model: ModelViewProtocol) -> None:
"""Move the state of one member, or of the members of its cell, from the states the members before it set."""
if agent.alive:
self._apply(self._get_moves(agent, model), model)
def move(self, model: ModelViewProtocol) -> None:
"""Move the state of every agent of the target that moves in this step, by the schedule.
Args:
model:
The model view of the interaction.
"""
agents = list(model.agents(self._get_role('target')))
if self.schedule == 'synchronous':
moves = [move for agent in agents for move in self._get_moves(agent, model)]
self._apply(moves, model)
return
order = model.rng.permutation(len(agents)) if self.schedule == 'random-sequential' else range(len(agents))
for index in order:
self._move_member(agents[int(index)], model)
def _get_laws(self, namespace: SimpleNamespace) -> list[BasePart]:
"""Return the jump, the update, the transition or emission law of the state, or the rule of the agent form.
Args:
namespace:
The namespace of the symbols the laws can use, with `namespace.target`.
Returns:
laws:
The law of the move, or the rule that moves the agents.
"""
if self.form == 'agent' and self.schedule == 'member':
return [RuleLaw(self._move_member, phase=self.phase, members=self._get_role('target'))]
if self.form == 'agent':
return [RuleLaw(self.move, phase=self.phase)]
value = getattr(self._get_member(namespace, 'target'), self.state)
if self.matrix is not None:
return self._get_chain(namespace, value)
if self.rate is not None:
rate = self._combine_rates(SR(self.rate(namespace)))
return [self._get_jump(value, self.source, self.destination, rate)]
if self.probability is not None:
probability = self._combine_probabilities(SR(self.probability(namespace)))
moves = self._get_draw() < probability
else:
moves = SR(self.condition(namespace)) if self.condition is not None else SR(0)
if self.after is not None:
moves = both(residence(value) + 1 >= self.after, moves)
return [self._get_step(value, self.source, self.destination, moves)]
Equations
source.hidden: hidden_source(t + 1) == where(categorical__transition_hidden_source_chain*(where(hidden_source(t) == 0, 0.9, 0.4) + where(hidden_source(t) == 0, 0.1, 0.6)) >= where(hidden_source(t) == 0, 0.9, 0.4), 1, 0)
source.symbol: symbol_source(t + 1) == where(categorical__emission_symbol_source_chain*(where(where(categorical__transition_hidden_source_chain*(where(hidden_source(t) == 0, 0.9, 0.4) + where(hidden_source(t) == 0, 0.1, 0.6)) >= where(hidden_source(t) == 0, 0.9, 0.4), 1, 0) == 0, 0.7, 0.2) + where(where(categorical__transition_hidden_source_chain*(where(hidden_source(t) == 0, 0.9, 0.4) + where(hidden_source(t) == 0, 0.1, 0.6)) >= where(hidden_source(t) == 0, 0.9, 0.4), 1, 0) == 0, 0.3, 0.8)) >= where(where(categorical__transition_hidden_source_chain*(where(hidden_source(t) == 0, 0.9, 0.4) + where(hidden_source(t) == 0, 0.1, 0.6)) >= where(hidden_source(t) == 0, 0.9, 0.4), 1, 0) == 0, 0.7, 0.2), 1, 0)