Model sub-package

class TimeSeriesSRC.Model.model.Parameters(epochs: int, goal: float, min_grad: float, mu: float, mu_dec: float, mu_inc: float, mu_max: float, show: int, max_time: float, delta: float)

Bases: object

epochs: int

goal: float

min_grad: float

mu: float

mu_dec: float

mu_inc: float

mu_max: float

show: int

max_time: float

delta: float

class TimeSeriesSRC.Model.model.pmodel(xtype, na=[0, 1], nb=[], nc=[], nd=[], nf=[], delay=[], diff=[0], per=[], upre=[], ypre=[], ypost=[], eFcn='estimlm', indexFcn='pmodmse', initFcn='initrand')

Bases: object

Prediction model for discrete-time ARMA / ARX / ARMAX / BJTF / regression.

Stores the polynomial orders, delay, differencing specification, and parameter arrays for one of the five supported model structures. After construction the object is passed to estimate() to fit the parameters, and to predict() to generate one-step-ahead predictions.

Parameters:

• xtype ({'arma', 'arx', 'armax', 'bjtf', 'regr'}) –

  Model structure.

  • 'arma' — ARMA(nc, nd) noise model, no external input.

  • 'arx' — ARX(na, nb) model with direct input feed-through.

  • 'armax' — ARMAX(na, nb, nc) model.

  • 'bjtf' — Box–Jenkins transfer-function model.

  • 'regr' — Static regression (no dynamics in the output).

• na (int or list of int, optional) – Order of the \(A(q)\) (autoregressive) polynomial. Default [0, 1].

• nb (list of int, optional) – Orders of the \(B(q)\) (input numerator) polynomials; one entry per input channel. Default [].

• nc (int or list of int, optional) – Order(s) of the \(C(q)\) (noise numerator) polynomial(s). Default [].

• nd (list of int, optional) – Order(s) of the \(D(q)\) (noise denominator) polynomial(s). Default [].

• nf (list of int, optional) – Order(s) of the \(F(q)\) polynomial(s) — BJTF only. Default [].

• delay (list of int, optional) – Pure input delay \(k\) for each input channel. Default [].

• diff (list of int, optional) – Differencing orders; diff[0] is the regular order and diff[i] (i > 0) is the seasonal order at period per[i-1]. Default [0].

• per (list of int, optional) – Seasonal periods; must have one fewer element than diff. Default [].

• upre (list of str, optional) – Names of pre-processing functions to apply to each input channel (e.g. ['log']). Default [].

• ypre (list of str, optional) – Names of pre-processing functions to apply to the output before estimation (e.g. ['log']). Default [].

• ypost (list of str, optional) – Names of post-processing functions to apply to the output after simulation (e.g. ['exp']). Default [].

• eFcn (str, optional) – Estimation function identifier. Only 'estimlm' (Levenberg–Marquardt) is currently supported. Default 'estimlm'.

• indexFcn (str, optional) – Performance-index function used during estimation. Default 'pmodmse'.

• initFcn ({'initrand', 'initrandn', 'initzero'}, optional) – Parameter initialisation strategy. Default 'initrand'.

a, b, c, d, f

Estimated polynomial coefficient vectors after calling estimate().

Type:

list of ndarray

estimParams

Levenberg–Marquardt hyper-parameters (epochs, goal, learning-rate schedule, etc.).

Type:

Parameters

Examples

Build a BJTF(1,1,1,1) model with one input and one output:

>>> from TimeSeriesSRC.Model.model import pmodel
>>> pm = pmodel('bjtf', nb=[1], nc=[1], nd=[1], nf=[1], delay=[0])
>>> print(pm)

Build an ARX(2, [1]) model:

>>> pm = pmodel('arx', na=2, nb=[1], delay=[1])

See also

estimate

Fit the model parameters to data.

func_selpmod

Grid search over candidate model structures.

func_pmodsim

Simulate the model output.

set_data(y, u=array([], shape=(1, 0), dtype=float64))

new_model()

newregr()

newbjtf()

newarx()

newarmax()

newarma()

init()

initzero()

initrandn()

initrand()

getmX()

getmXregr()

getmXbjtf()

getmXarx()

getmXarmax()

getmXarma()

setmX(X)

setmXregr(X)

setmXbjtf(X)

setmXarx(X)

setmXarmax(X)

setmXarma(X)

predict(y, u=[])

predregr(y, u=array([], shape=(1, 0), dtype=float64))

PREDICT Compute one-step predictions for regression model.

Parameters:

• y

• u

Returns:

predbjtf(y, u=[])

predarx(y, u=array([], shape=(1, 0), dtype=float64))

predarmax(y, u=[[]])

predarma(y=[])

predictdf(y, u=[[]])

preddfbjtf(y, u)

PREDICT Compute one-step predictions for the Box and Jenkins Transfer Function model. :param y: :param u: :return:

preddfarma(y)

getGHdf()

getGHdfbjtf()

getGHdfarma()

getGH()

getGHbjtf()

getGHarx()

getGHarmax()

getGHarma()

TimeSeriesSRC.Model.estimate.estimate(pmod, y, u=array([], dtype=float64), show_plot=True, show_output=True)

Estimate prediction-model parameters using the Levenberg–Marquardt algorithm.

Pre-processes y (and optionally u) according to the transforms and differencing orders stored in pmod, then calls the LM estimator (func_estimlm()) and returns the fitted model together with the full training record.

Parameters:

• pmod (pmodel) – Prediction model object created by pmodel. Its parameter arrays (a, b, c, d, f) are updated in-place during estimation.

• y (array-like, shape (1, N)) – Output (response) time series. Should be zero-mean (or near zero-mean) after applying the differencing specified in pmod.diff.

• u (array-like, shape (n_inputs, N), optional) – Input time series matrix. Omit for purely ARMA models. Default is an empty array.

• show_plot (bool, optional) – Display the live training-performance (MSE vs epoch) plot. Default True.

• show_output (bool, optional) – Print per-epoch training summary to stdout. Default True.

Returns:

• pmod (pmodel) – The fitted model with updated parameter arrays.

• trec (dict) – Training record. Key 'index' holds the per-epoch MSE vector; 'epoch' holds the number of completed epochs.

• stat (dict) – Final-epoch summary statistics (MSE, gradient norm, Jacobian condition number).

Examples

>>> import numpy as np
>>> from TimeSeriesSRC.Model.model import pmodel
>>> from TimeSeriesSRC.Model.estimate import estimate
>>> y = np.array([-0.19, 0.52, -3.50, 3.01, -3.04, 1.59]).reshape(1, -1)
>>> u = np.array([-0.43, -1.67,  0.13, 0.29, -1.15, 1.19]).reshape(1, -1)
>>> pm = pmodel('bjtf', nb=[1], nc=[1], nd=[1], nf=[1], delay=[0])
>>> pm.estimParams.epochs = 5
>>> pm, trec, stat = estimate(pm, y, u, show_plot=False, show_output=False)

See also

pmodel

Prediction model object with polynomial-order specification.

func_selpmod

Automated grid search over candidate model structures.

func_pmodmse

Compute mean-squared prediction error for a fitted model.

TimeSeriesSRC.Model.selpmod.func_selpmod(filename, y, u=[])

Select the best prediction model by grid-searching AIC and BIC.

Exhaustively searches the order combinations defined in a JSON spec file (or dict), fits each candidate model with estimate(), and returns the best models according to both Akaike (AIC) and Bayesian (BIC) criteria.

Parameters:

• filename (str or dict) – Path to a JSON specification file or a dict with the same structure. The JSON must have a "models" list; each entry is a dict with a "type" key ('arma', 'arx', 'armax', 'bjtf', 'regr') and order keys such as na, nb, nc, nd, nf, diff, delay, per — each holding a list of candidate integer values to try.

• y (array-like) – Observed output time series.

• u (array-like, optional) – Input time series (required for ARX, ARMAX, BJTF, regr models). Default [].

Returns:

result – Keyed by model type (e.g. 'arma', 'arx'). Each value is a dict with:

• 'model' — model type string.

• 'aic' — list of AIC values for every combination tried.

• 'bic' — list of BIC values for every combination tried.

• 'aicstat' — stat dict from estimate() for the AIC winner.

• 'bicstat' — stat dict from estimate() for the BIC winner.

• 'aicmod' — best pmodel selected by AIC.

• 'bicmod' — best pmodel selected by BIC.

Return type:

dict

Examples

>>> import pathlib, json
>>> import TimeSeriesSRC as ts
>>> spec = {
...     "models": [
...         {"type": "arma", "nc": [1, 2], "nd": [1], "diff": [0], "per": []}
...     ]
... }
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> import pandas as pd
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> result = ts.selpmod(spec, y)
>>> best = result['arma']['aicmod']

See also

estimate

Fit a single model.

pmodaic

Akaike Information Criterion.

pmodbic

Bayesian Information Criterion.

TimeSeriesSRC.Model.pmoddisp.func_pmoddisp(pmod, stat)

Print a parameter table and produce an error-bar plot.

Displays each estimated parameter with its ±2σ confidence interval in a formatted table, then plots the same information as a horizontal error-bar chart.

Parameter ordering follows pmod.getmX():

• bjtf — b, c, d, f (one block per input / seasonal period)

• armax — a, b, c

• arx — a, b

• arma — c, d

• regr — b

Parameters:

• pmod (pmodel) – Estimated prediction model as returned by estimate().

• stat (dict) –

  Statistics dictionary returned by estimate():

  • stat['sigma'] — residual variance \(\hat{\sigma}^2\).

  • stat['stdx'] — standard deviation of each parameter (same order as pmod.getmX()).

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> pm = ts.pmodel('arma', nc=[2], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = ts.estimate(pm, y, show_plot=False, show_output=False)
>>> ts.pmoddisp(pm_est, stat)

See also

func_pmodpzplot

Pole-zero map for stability/invertibility check.

estimate

Returns the stat dict consumed here.

TimeSeriesSRC.Model.pmoddisp.func_pmodpzplot(pmod)

Plot poles and zeros of the G and H transfer functions.

Displays an annotated pole-zero map on the complex plane with the unit circle as a stability/invertibility reference.

Polynomial conventions (all in the forward-shift operator \(q\)):

• B polynomial — [b0, b1, ..., b_nb] (zeros of G)

• F polynomial — [1, f1, ..., f_nf], combined with A → poles of G

• A polynomial — [1, a1, ..., a_na]

• C polynomial — [1, c1, ..., c_nc] (zeros of H, regular part)

• D polynomial — [1, d1, ..., d_nd], combined with A → poles of H

For seasonal models (pmod.period non-empty) the H transfer function is decomposed into one panel per seasonal component. For arma models G is omitted. For regr (static) models the function prints a notice and returns immediately.

Parameters:

pmod (pmodel) – Estimated prediction model (result of estimate() or manually constructed).

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> pm = ts.pmodel('arma', nc=[2], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = ts.estimate(pm, y, show_plot=False, show_output=False)
>>> ts.pmodpzplot(pm_est)

See also

func_pmoddisp

Parameter table and error-bar plot.

estimate

Fit model parameters from data.

TimeSeriesSRC.Model.pmodsim.func_pmodsim(pmod, e, u=[])

Simulate the output of a prediction model driven by white noise.

Implements the system:

\[y(t) = G(q)\,u(t) + H(q)\,e(t)\]

where \(G(q)\) and \(H(q)\) are the transfer functions of the fitted model, and \(e(t)\) is a user-supplied white noise sequence.

Parameters:

• pmod (pmodel) – Fitted or manually constructed prediction model.

• e (array-like) – White noise input sequence with the same length as u.

• u (array-like, optional) – External input sequence. Required for ARX, ARMAX, BJTF, and regression models. Default [].

Returns:

y – Simulated output sequence.

Return type:

ndarray

Examples

>>> import numpy as np
>>> import TimeSeriesSRC as ts
>>> pm = ts.pmodel('arma', nc=[1], nd=[1], diff=[0], per=[])
>>> pm.c[0] = np.array([-0.5])
>>> pm.d[0] = np.array([-0.8])
>>> e = np.random.default_rng(0).standard_normal(200)
>>> y = ts.pmodsim(pm, e)

See also

estimate

Fit model parameters from data.

pmodel

Define model structure.

TimeSeriesSRC.Model.pmodaic.func_pmodaic(pmod, y, u=[])

Compute the Akaike Information Criterion (AIC) for a fitted model.

\[\text{AIC} = \ln(\text{MSE}) + \frac{2 k}{N}\]

where \(k\) is the number of free parameters and \(N\) is the number of data points. Lower AIC indicates a better trade-off between fit and model complexity.

Parameters:

• pmod (pmodel) – Fitted prediction model.

• y (array-like) – Desired output sequence.

• u (array-like, optional) – Input sequence. Default [].

Returns:

aic – Akaike Information Criterion.

Return type:

float

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> pm = ts.pmodel('arma', nc=[2], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = ts.estimate(pm, y, show_plot=False, show_output=False)
>>> aic = ts.pmodaic(pm_est, y)

See also

pmodbic

Bayesian Information Criterion (penalises complexity more strongly).

selpmod

Automatic model selection using AIC/BIC grid search.

TimeSeriesSRC.Model.pmodbic.func_pmodbic(pmod, y, u=[])

Compute the Bayesian Information Criterion (BIC) for a fitted model.

\[\text{BIC} = \ln(\text{MSE}) + \frac{k \ln N}{N}\]

where \(k\) is the number of free parameters and \(N\) is the number of data points. BIC penalises model complexity more strongly than AIC for large \(N\). Lower BIC indicates a better model.

Parameters:

• pmod (pmodel) – Fitted prediction model.

• y (array-like) – Desired output sequence.

• u (array-like, optional) – Input sequence. Default [].

Returns:

bic – Bayesian Information Criterion.

Return type:

float

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> pm = ts.pmodel('arma', nc=[1], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = ts.estimate(pm, y, show_plot=False, show_output=False)
>>> bic = ts.pmodbic(pm_est, y)

See also

pmodaic

Akaike Information Criterion.

selpmod

Automatic model selection using AIC/BIC grid search.

TimeSeriesSRC.Model.pmodmse.func_pmodmse(pmod, y, u=[])

Compute the mean squared prediction error (MSE) for a fitted model.

Parameters:

• pmod (pmodel) – Fitted prediction model.

• y (array-like) – Desired output sequence.

• u (array-like, optional) – Input sequence (required for ARX, ARMAX, BJTF models). Default [] (univariate models).

Returns:

• mse (float) – Mean squared one-step-ahead prediction error.

• e (ndarray) – Prediction error sequence y - yhat.

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> y = pd.read_csv(data_dir / 'Series_A_Chemical_Concentration.csv').values.flatten()
>>> pm = ts.pmodel('arma', nc=[2], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = ts.estimate(pm, y, show_plot=False, show_output=False)
>>> mse, e = ts.pmodmse(pm_est, y)

See also

pmodaic

Akaike Information Criterion.

pmodbic

Bayesian Information Criterion.

TimeSeriesSRC.Model.estimlm.plotindex(trec, goal, name, epoch)

Update the live training performance plot.

Parameters:

• trec (dict) – Training record (trec['index'] holds the per-epoch MSE values).

• goal (float) – Target performance value. A horizontal dashed line is drawn if finite and positive.

• name (str) – Algorithm name shown in the figure title.

• epoch (int) – Current epoch index (0-based); only data up to this epoch is plotted.

TimeSeriesSRC.Model.estimlm.reset_plotindex()

Reset the plot figure for a new training session.

TimeSeriesSRC.Model.estimlm.func_estimlm(pmod, y=[], u=[], show_plot=True, show_output=True)

Fit a prediction model with the Levenberg-Marquardt algorithm.

Minimises the sum of squared one-step-ahead prediction errors using an iterative Levenberg-Marquardt update:

\[\Delta X = -(J^T J + \mu I)^{-1} J^T E\]

where \(J\) is the numerical Jacobian of the prediction errors with respect to the model parameters \(X\), \(E\) is the error vector, and \(\mu\) is the adaptive damping factor.

Parameters:

• pmod (pmodel) – Prediction model with initial parameter values and estimation hyper-parameters in pmod.estimParams.

• y (array-like or dict) – Desired output sequence. If a dict, must have keys 'y' (output array) and 'm' (sample-weight vector).

• u (array-like, optional) – Input sequence. Default [] (univariate models).

• show_plot (bool, optional) – Display the live performance-index plot during training. Default True.

• show_output (bool, optional) – Print epoch-by-epoch progress to stdout. Default True.

Returns:

• pmod (pmodel) – Fitted model with updated parameter values.

• trec (dict) – Training record with keys:

  • 'index' — performance index at each epoch.

  • 'mu' — damping factor \(\mu\) at each epoch.

• stat (dict) – Final statistics:

  • 'sigma' — residual variance \(\hat{\sigma}^2\).

  • 'stdx' — standard deviations of parameter estimates.

Notes

Default estimation hyper-parameters (set on pmod.estimParams):

Attribute	Default	Meaning
epochs	10	Maximum iterations
goal	0	Stop when performance ≤ goal
mu	0.001	Initial damping factor
mu_dec	0.1	Multiply mu by this on a successful step
mu_inc	10	Multiply mu by this on a failed step
mu_max	1e10	Stop if mu exceeds this value
min_grad	1e-4	Stop if gradient norm falls below this
max_time	inf	Wall-clock time limit (seconds)
delta	1e-7	Finite-difference step for Jacobian

See also

estimate

Public interface — calls this function internally.

func_jacobian

Numerical Jacobian calculation.

TimeSeriesSRC.Model.jacobian.func_jacobian(pmod, delta, y, u=array([], dtype=float64))

Compute the numerical Jacobian, approximate Hessian, and gradient.

Uses a forward-difference scheme to approximate the Jacobian \(J\) of the one-step-ahead prediction errors with respect to the model parameters \(X\):

\[J_{ij} \approx \frac{\hat{y}(t_i;\, X + \delta e_j) - \hat{y}(t_i;\, X)}{\delta}\]

The approximate Gauss-Newton Hessian and gradient are then:

\[JJ = J M J^T, \qquad je = J M E\]

where \(M\) is a diagonal weight matrix and \(E\) is the error vector.

Parameters:

• pmod (pmodel) – Prediction model with current parameter values.

• delta (float) – Finite-difference step size (e.g. 1e-7).

• y (array-like or dict) – Desired output sequence. If a dict, must have keys 'y' and 'm' (sample weights).

• u (array-like, optional) – Input sequence. Default np.array([]).

Returns:

• je (ndarray, shape (n_params, 1)) – Gradient \(J M E\).

• jj (ndarray, shape (n_params, n_params)) – Approximate Hessian \(J M J^T\).

• normgX (float) – Euclidean norm of the gradient vector.

See also

func_estimlm

Levenberg-Marquardt optimiser that calls this function.