AlssmPolyMeixner

AlssmPolyMeixner(poly_degree: int, segment, **kwargs)

Bases: ModelBase

              flowchart TD
              lmlib.statespace.model.AlssmPolyMeixner[AlssmPolyMeixner]
              lmlib.statespace.model.ModelBase[ModelBase]

                              lmlib.statespace.model.ModelBase --> lmlib.statespace.model.AlssmPolyMeixner
                


              click lmlib.statespace.model.AlssmPolyMeixner href "" "lmlib.statespace.model.AlssmPolyMeixner"
              click lmlib.statespace.model.ModelBase href "" "lmlib.statespace.model.ModelBase"
            

ALSSM whose output basis is the Meixner polynomials, orthogonal under the geometric (exponential) weight \(\gamma^j\) on \(j = 0, 1, 2, \ldots\)

Unlike AlssmPoly (Pascal/monomial basis), which suffers from Gram matrix condition numbers of \(\mathcal{O}(g^{2D})\), and AlssmPolyLegendre, which requires a finite window specification, AlssmPolyMeixner is designed for infinite or semi-infinite exponential windows and keeps the composite Gram matrix near the theoretical minimum.

Mathematical background

The Meixner polynomials \(M_n(j;\,1,\gamma)\) satisfy

\[ \sum_{j=0}^{\infty} \gamma^j\, M_m(j)\, M_n(j) = W_{n}\,\delta_{mn}, \qquad W_n = g\!\left(\frac{g}{g-1}\right)^{\!n} = \frac{1}{(1-\gamma)\,\gamma^n}, \]

where \(\gamma = (g-1)/g < 1\) is the exponential decay of the backward segment and \(g > 1\) is its effective window length. The norms \(W_n\) are exact and closed-form, growing by the constant factor \(g/(g-1)\) per degree. \(\delta_{mn}\) is the kronecker delta.

Accepting a Segment

The recommended constructor is AlssmPolyMeixner(poly_degree, segment) which infers all parameters from the Segment:

• segment.direction → selects \(A\):

• 'bw' (backward): \(A = A_{\rm bw} = I - \tfrac{1}{g-1}\operatorname{triu}(\mathbf{1},1)\) — basis \(C A_{\rm bw}^j x = M_n(j;\gamma)\) at lag \(j \ge 0\).

• 'fw' (forward): \(A = A_{\rm fw} = A_{\rm bw}^{-1}\) — the forward filter's internal \(A^{-1} = A_{\rm bw}\) step recovers the decaying Meixner basis at lags \(j \le 0\).

• segment.a (backward) or segment.b (forward) → shifts \(C\) so that the Gram matrix remains \(W_{\rm ss}\) (diagonal) even when the segment does not start at \(j=0\):

• Backward \([a, \infty)\): \(C \leftarrow [1,\ldots,1]\,A_{\rm bw}^{-a}\) (makes \(C A^j x = M_n(j-a;\gamma)\) for \(j \ge a\)).

• Forward \((-\infty, b]\): \(C \leftarrow [1,\ldots,1]\,A_{\rm bw}^{-b}\) (makes the relative basis start at the boundary \(j=b\)).

The g and direction attributes are always available.

Condition number

\(\kappa(W) = (g/(g-1))^D\), independent of the segment shift.

Notes

The Meixner polynomials used here are the standard (monic normalisation) \(M_n(x;\,1,\gamma) = {}_2F_1(-n,\,-x;\,1;\,1 - 1/\gamma)\).

Parameters:

• poly_degree (int) –

  Polynomial degree \(D \geq 0\). The model order is \(N = D+1\).

• segment (Segment) –

  Target segment. Encodes direction, \(g\), and boundary shift:

  • segment.g → window size \(g\) (giving decay \(\gamma = (g-1)/g\)).
  • segment.direction → 'bw' selects \(A_{\rm bw}\); 'fw' selects \(A_{\rm fw} = A_{\rm bw}^{-1}\).
  • segment.a (BW, if finite) or segment.b (FW, if finite) → shifts \(C\) to restore Gram matrix orthogonality when the segment does not start at the canonical origin.

  Meixner orthogonality requires a semi-infinite support, so the segment must be backward with b=+inf or forward with a=-inf.

• **kwargs –

  Forwarded to ModelBase (e.g. label).

Methods:

• update –

  Recompute \(A\) and \(C\) from all stored parameters.

• eval_output –

  Evaluate the ALSSM output for one or more state vectors.

• dump_tree –

  Return the internal ALSSM tree structure as a string.

• set_state_var_label –

  Register a label for one or more state vector indices.

• get_state_var_labels –

  Return all registered state-variable labels together with their index tuples.

• get_state_var_indices –

  Return the state-vector indices for a state variable identified by its label.

• get_alssm_output_dimension –

  Return the ALSSM output dimension \(Q\) (number of output channels).

Attributes:

• steady_state_basis –
• poly_degree –

  int : Polynomial degree \(D\).

• g –

  float : Effective window size \(g\) (read-only; from the segment).

• direction –

  str : Segment direction 'bw' (backward) or 'fw' (forward).

• shift –

  int : Boundary shift applied to \(C\); 0 for the canonical origin.

• gamma –

  float : Exponential decay \(\gamma = (g-1)/g\).

• label –

  str : Label of the model

• C_init –

  ndarray, shape=([Q,] N) : Initialized Output matrix \(C \in \mathbb{R}^{Q \times N}\)

• force_MC –

  bool : If True, a 1-D output vector C is broadcast to a 2-D array of shape (1, N) (multi-channel form).

• A –

  ndarray, shape=(N, N) : State matrix \(A \in \mathbb{R}^{N \times N}\)

• C –

  ndarray, shape=([Q,] N) : Output matrix \(C \in \mathbb{R}^{Q \times N}\)

• N –

  int : Model order \(N\)

• Q –

  int : Number of output channels \(Q\).

• alssms –

  list : Sub-ALSSMs that compose this model (empty for leaf nodes such as Alssm).

• lambdas –

  ndarray : Per-ALSSM scalar output scaling factors \(\lambda_m\) applied to each sub-model's output matrix \(C_m\).

• is_MC –

  bool : True if the output matrix C is 2-D (multi-channel form), False if 1-D (scalar output).

Source code in lmlib/statespace/model.py

def __init__(self, poly_degree: int, segment, **kwargs):
    r"""
    Parameters
    ----------
    poly_degree : int
        Polynomial degree $D \geq 0$.  The model order is $N = D+1$.
    segment : Segment
        Target segment.  Encodes direction, $g$, and boundary shift:

        * ``segment.g`` → window size $g$ (giving decay
          $\gamma = (g-1)/g$).
        * ``segment.direction`` → ``'bw'`` selects $A_{\rm bw}$;
          ``'fw'`` selects $A_{\rm fw} = A_{\rm bw}^{-1}$.
        * ``segment.a`` (BW, if finite) or ``segment.b`` (FW, if finite)
          → shifts $C$ to restore Gram matrix orthogonality when the
          segment does not start at the canonical origin.

        Meixner orthogonality requires a semi-infinite support, so the
        segment must be backward with ``b=+inf`` or forward with ``a=-inf``.
    **kwargs
        Forwarded to [`ModelBase`][lmlib.statespace.model.ModelBase] (e.g. ``label``).
    """

    super().__init__(**kwargs)
    assert isinstance(poly_degree, int) and poly_degree >= 0, \
        'poly_degree must be a non-negative int'
    self._poly_degree = int(poly_degree)

    assert isinstance(segment, Segment), \
        "'segment' must be a lmlib.statespace.segment.Segment instance."
    # ── Validate: Meixner orthogonality requires a semi-infinite segment ──
    if segment.direction == BACKWARD and np.isfinite(segment.b):
        raise ValueError(
            f"AlssmPolyMeixner requires b=+inf for BACKWARD segments "
            f"(got b={segment.b}). The Meixner polynomials are orthogonal only on "
            f"a semi-infinite support. For finite windows use AlssmPoly or "
            f"AlssmPolyLegendre instead.")
    if segment.direction == FORWARD and np.isfinite(segment.a):
        raise ValueError(
            f"AlssmPolyMeixner requires a=-inf for FORWARD segments "
            f"(got a={segment.a}). The Meixner polynomials are orthogonal only on "
            f"a semi-infinite support. For finite windows use AlssmPoly or "
            f"AlssmPolyLegendre instead.")
    self._g = float(segment.g)
    self._direction = segment.direction
    # ── Determine C shift from the finite boundary ──────────────────────
    #
    # BACKWARD [a, ∞):
    #   C ← [1..1] A_bw^{-a}  →  output at lag j = M_n(j-a; γ).
    #   This keeps W = W_ss regardless of a (the Stein boundary term Q_b
    #   = C^T C = ones at relative lag 0).
    #
    # FORWARD (-∞, b]:
    #   A = A_fw = A_bw^{-1}.
    #   C ← [1..1] A_bw^{-|b|} = [1..1] A_fw^{|b|}
    #   This makes Q_b = outer([1..1],[1..1]) in the Stein equation,
    #   keeping W = W_ss regardless of b (standard b=-1 → shift=1).
    if segment.direction == BACKWARD:
        self._shift = 0 if np.isinf(segment.a) else int(segment.a)
    else:  # FORWARD
        self._shift = 0 if np.isinf(segment.b) else int(abs(segment.b))

    self.update()

Methods

update

update()

Recompute \(A\) and \(C\) from all stored parameters.

Transition matrix \(A\):

• direction 'bw': \(A = A_{\rm bw} = I_N - \tfrac{1}{g-1}\,\operatorname{triu}(\mathbf{1}_N, 1)\)
• direction 'fw': \(A = A_{\rm fw} = A_{\rm bw}^{-1}\)

Output vector \(C\) (before shift): \(C_{\rm base} = [1,\ldots,1]\).

Shift (when shift != 0): \(C \leftarrow C_{\rm base}\,A_{\rm bw}^{-\text{shift}}\).

This ensures that for a backward segment \([a, \infty)\) with shift = a, the output at absolute lag \(j\) is \(M_n(j - a;\,\gamma)\) — orthogonal under the shifted window weight.

Source code in lmlib/statespace/model.py

def update(self):
    r"""
    Recompute $A$ and $C$ from all stored parameters.

    **Transition matrix** $A$:

    * direction ``'bw'``:
      $A = A_{\rm bw} = I_N - \tfrac{1}{g-1}\,\operatorname{triu}(\mathbf{1}_N, 1)$
    * direction ``'fw'``:
      $A = A_{\rm fw} = A_{\rm bw}^{-1}$

    **Output vector** $C$ (before shift):
    $C_{\rm base} = [1,\ldots,1]$.

    **Shift** (when ``shift != 0``):
    $C \leftarrow C_{\rm base}\,A_{\rm bw}^{-\text{shift}}$.

    This ensures that for a backward segment $[a, \infty)$ with
    ``shift = a``, the output at absolute lag $j$ is
    $M_n(j - a;\,\gamma)$ — orthogonal under the shifted window weight.
    """
    N = self._poly_degree + 1
    A_bw = np.eye(N) - (1.0 / (self._g - 1.0)) * np.triu(np.ones((N, N)), 1)

    if self._direction == 'fw':
        self.A = np.linalg.inv(A_bw)
    else:
        self.A = A_bw

    C_base = np.ones(N)
    if self._shift != 0:
        # C @ A_bw^{-shift}: matrix_power handles negative exponents via inverse.
        C_base = C_base @ matrix_power(A_bw, -self._shift)
    self.C = C_base

    self._init_state_var_labels()
    self._broadcast_C_to_multichannel()

eval_output

eval_output(xs, js=None)

Evaluate the ALSSM output for one or more state vectors.

Without evaluation index (js=None):

\[ s(x) = C x \]

With evaluation indices (js provided):

\[ s_j(x) = C A^j x \]

Parameters:

• xs (array_like of shape (..., N)) –

  State vector(s). The last dimension must equal the model order N.

• js (array_like of shape (J,) or None, default: None ) –

  Sequence of integer evaluation indices. If None, evaluates at \(j = 0\) only (i.e. returns \(Cx\)).

Returns:

• s ( ndarray ) –

  If js is None: shape (..., [Q]). If js is provided: shape (J, ..., [Q]). The [Q] dimension is present only when is_MC is True.

Source code in lmlib/statespace/model.py

def eval_output(self, xs, js=None):
    r"""
    Evaluate the ALSSM output for one or more state vectors.

    Without evaluation index (``js=None``):

    $$
    s(x) = C x
    $$

    With evaluation indices (``js`` provided):

    $$
    s_j(x) = C A^j x
    $$

    Parameters
    ----------
    xs : array_like of shape (..., N)
        State vector(s). The last dimension must equal the model order N.
    js : array_like of shape (J,) or None, optional
        Sequence of integer evaluation indices. If None, evaluates at
        $j = 0$ only (i.e. returns $Cx$).

    Returns
    -------
    s : ndarray
        If ``js`` is None: shape ``(..., [Q])``.
        If ``js`` is provided: shape ``(J, ..., [Q])``.
        The ``[Q]`` dimension is present only when [`is_MC`][lmlib.statespace.model.ModelBase.is_MC] is True.
    """
    xs = np.asarray(xs)

    # Ensure last axis is state dimension
    if xs.shape[-1] != self.N:
        raise ValueError(f"Last dimension of xs must be {self.N}")

    # No propagation: s = C x
    if js is None:
        _subscript = 'ln,...n->...l' if self.is_MC else 'n,...n->...'
        return np.einsum(_subscript, self.C, xs)

    # Propagation: s_j = C A^j x
    A_powers = [matrix_power(self.A, int(j)) for j in js]
    _subscript = 'ln,...n->...l' if self.is_MC else 'n,...n->...'
    return np.asarray([np.einsum(_subscript, self.C @ Aj, xs) for Aj in A_powers])

dump_tree

dump_tree() -> str

Return the internal ALSSM tree structure as a string.

Returns:

• out ( str ) –

  Multi-line string representing the nested ALSSM structure.

Example

>>> import lmlib as lm
>>> import numpy as np
>>> alssm_poly = lm.AlssmPoly(4, label="high order polynomial")
>>> A = [[1, 1], [0, 1]]
>>> C = [[1, 0]]
>>> alssm_line = lm.Alssm(A, C, label="line")
>>> stacked_alssm = lm.AlssmStacked((alssm_poly, alssm_line), label='stacked model')
>>> print(stacked_alssm.dump_tree())
└-AlssmStacked, A: (7, 7), C: (2, 7), label: stacked model
  └-AlssmPoly, A: (5, 5), C: (5,), label: high order polynomial
  └-Alssm, A: (2, 2), C: (1, 2), label: line

Source code in lmlib/statespace/model.py

def dump_tree(self) -> str:
    """
    Return the internal ALSSM tree structure as a string.

    Returns
    -------
    out : str
        Multi-line string representing the nested ALSSM structure.

    Example
    --------
    ```python
    >>> import lmlib as lm
    >>> import numpy as np
    >>> alssm_poly = lm.AlssmPoly(4, label="high order polynomial")
    >>> A = [[1, 1], [0, 1]]
    >>> C = [[1, 0]]
    >>> alssm_line = lm.Alssm(A, C, label="line")
    >>> stacked_alssm = lm.AlssmStacked((alssm_poly, alssm_line), label='stacked model')
    >>> print(stacked_alssm.dump_tree())
    └-AlssmStacked, A: (7, 7), C: (2, 7), label: stacked model
      └-AlssmPoly, A: (5, 5), C: (5,), label: high order polynomial
      └-Alssm, A: (2, 2), C: (1, 2), label: line
    ```
    """
    return self._rec_tree(level=0)

set_state_var_label

set_state_var_label(label: str, indices: tuple[int])

Register a label for one or more state vector indices.

Labels allow state components to be referenced by name rather than by numeric index; see get_state_var_indices.

Parameters:

• label (str) –

  Label name to register.

• indices (tuple of int) –

  State vector indices associated with this label.

Example

>>> import lmlib as lm
>>> alssm = lm.AlssmPoly(poly_degree=1, label='slope_with_offset')
>>> alssm.set_state_var_label('slope', (1,))
>>> alssm.get_state_var_indices('slope_with_offset.slope')
(1,)

Source code in lmlib/statespace/model.py

def set_state_var_label(self, label:str, indices:tuple[int]):
    r"""
    Register a label for one or more state vector indices.

    Labels allow state components to be referenced by name rather than by
    numeric index; see [`get_state_var_indices`][lmlib.statespace.model.ModelBase.get_state_var_indices].

    Parameters
    ----------
    label : str
        Label name to register.
    indices : tuple of int
        State vector indices associated with this label.

    Example
    --------
    ```python
    >>> import lmlib as lm
    >>> alssm = lm.AlssmPoly(poly_degree=1, label='slope_with_offset')
    >>> alssm.set_state_var_label('slope', (1,))
    >>> alssm.get_state_var_indices('slope_with_offset.slope')
    (1,)
    ```
    """
    self._state_var_labels[label] = indices

get_state_var_labels

get_state_var_labels()

Return all registered state-variable labels together with their index tuples.

Labels are accumulated recursively from all nested sub-ALSSMs, with each label prefixed by the current model's label. The state indices are adjusted to reflect the position within the combined (block-diagonal) state vector.

Returns:

• out ( list of (str, tuple of int) ) –

  List of (label_string, indices) pairs. label_string is a dot-separated path (e.g. 'stacked.poly.x0') and indices is the corresponding tuple of integer state-vector positions.

Source code in lmlib/statespace/model.py

def get_state_var_labels(self):
    r"""
    Return all registered state-variable labels together with their index tuples.

    Labels are accumulated recursively from all nested sub-ALSSMs, with
    each label prefixed by the current model's [`label`][lmlib.statespace.model.ModelBase.label].  The state
    indices are adjusted to reflect the position within the combined
    (block-diagonal) state vector.

    Returns
    -------
    out : list of (str, tuple of int)
        List of ``(label_string, indices)`` pairs.  ``label_string`` is a
        dot-separated path (e.g. ``'stacked.poly.x0'``) and ``indices`` is
        the corresponding tuple of integer state-vector positions.
    """
    state_list = []
    for var_label, indices in self._state_var_labels.items():
        state_list.append((self.label + '.' + var_label, indices))

    N = 0
    for alssm in self.alssms:
        for var_label, indices in alssm.get_state_var_labels():
            state_list.extend([(self.label + '.' + var_label, tuple(i + N for i in indices))])
        N += alssm.N
    return state_list

get_state_var_indices

get_state_var_indices(label)

Return the state-vector indices for a state variable identified by its label.

Parameters:

• label (str) –

  Fully qualified state label (dot-separated path), as returned by get_state_var_labels.

Returns:

• out ( tuple of int or list of int ) –

  State-vector indices associated with label. Returns an empty list if label is not found.

Source code in lmlib/statespace/model.py

def get_state_var_indices(self, label):
    r"""
    Return the state-vector indices for a state variable identified by its label.

    Parameters
    ----------
    label : str
        Fully qualified state label (dot-separated path), as returned by
        [`get_state_var_labels`][lmlib.statespace.model.ModelBase.get_state_var_labels].

    Returns
    -------
    out : tuple of int or list of int
        State-vector indices associated with ``label``.
        Returns an empty list if ``label`` is not found.
    """

    for l, indices in self.get_state_var_labels():
        if label == l:
            return indices
    return []

get_alssm_output_dimension

get_alssm_output_dimension()

Return the ALSSM output dimension \(Q\) (number of output channels).

Source code in lmlib/statespace/model.py

def get_alssm_output_dimension(self):
    r"""Return the ALSSM output dimension $Q$ (number of output channels)."""
    return self.Q