basefunctions sub-package

TimeSeriesSRC.basefunctions.uniAnal.func_uniAnal(y, na=20, nump=10, nrg=5, ncg=0, diff=[0], per=[], perdsp=1)[source]

Compute and plot the ACF, PACF, and GPAC for a univariate time series.

Differences the series as requested, computes the autocorrelation (ACF), partial autocorrelation (PACF), and generalized partial autocorrelation (GPAC) functions, and produces stem/GPAC plots for order identification.

Parameters:

  • y (array-like) – 1-D time series.

  • na (int, optional) – Maximum lag for ACF; lags range from -na to +na. Default 20.

  • nump (int, optional) – Number of PACF terms to compute (lags 1 .. nump). Default 10.

  • nrg (int, optional) – Number of GPAC rows (numerator orders). Default 5.

  • ncg (int, optional) – Number of GPAC columns (denominator orders); 0 sets it equal to nrg. Default 0.

  • diff (list of int, optional) – Differencing orders to apply before analysis. Must satisfy len(diff) == len(per) + 1. Default [0].

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

  • perdsp (int, optional) – Display period — ACF and PACF are sampled at every perdsp-th lag. Default 1.

Returns:

  • yacf (ndarray, shape (1, 2*na+1)) – Biased ACF from lag -na to +na.

  • ypacf (ndarray, shape (1, nump)) – Partial ACF from lag 1 to nump.

  • ygpac (ndarray, shape (nrg, ncg)) – GPAC table (row = numerator order, column = denominator order).

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.uniAnal import func_uniAnal
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> yacf, ypacf, ygpac = func_uniAnal(y, na=20, nump=10, nrg=5)

See also

multiAnal

Analysis for input-output (transfer function) data.

uniChi

Box-Pierce chi-square whiteness test on residuals.

TimeSeriesSRC.basefunctions.multiAnal.func_multiAnal(u, y, nng=5, ndg=5, nnh=5, ndh=5, lg=12, lh=12)[source]

Multivariate analysis: impulse response, residual ACF, and GPAC tables.

Pre-whitens the input u with a BIC-selected ARMA model, estimates the impulse response between u and y via the Wiener-Hopf equations, and computes GPAC tables for both the G (transfer function) and H (noise) components. Produces plots of the impulse response, residual ACF, and the two GPAC arrays.

Parameters:

  • u (array-like) – 1-D input (exogenous) sequence.

  • y (array-like) – 1-D output sequence (same length as u).

  • nng (int, optional) – Maximum numerator order for the G GPAC table. Default 5.

  • ndg (int, optional) – Maximum denominator order for the G GPAC table. Default 5.

  • nnh (int, optional) – Maximum numerator order for the H GPAC table. Default 5.

  • ndh (int, optional) – Maximum denominator order for the H GPAC table. Default 5.

  • lg (int, optional) – Number of impulse-response lags to compute. Default 12.

  • lh (int, optional) – Accepted for API compatibility; not used in the current implementation. The residual-ACF lag count is determined by nnh + ndh + 1. Default 12.

Returns:

  • g (ndarray, shape (lg+1,)) – Estimated impulse response from u to y.

  • rv (ndarray, shape (1, 2*(nnh+ndh+1)+1)) – Residual autocorrelation function.

  • g_gpac (ndarray, shape (nng, ndg)) – GPAC table for the G transfer function.

  • h_gpac (ndarray, shape (nnh, ndh)) – GPAC table for the H noise model.

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.multiAnal import func_multiAnal
>>> rng = np.random.default_rng(0)
>>> u = rng.standard_normal(300)
>>> e = rng.standard_normal(300) * 0.2
>>> y = lfilter([1], [1, 0.5], u) + lfilter([1], [1, -0.8], e)
>>> g, rv, g_gpac, h_gpac = func_multiAnal(u, y, lg=10)

See also

uniAnal

Univariate ACF/PACF/GPAC analysis.

multiChi

Chi-square test for residual/input independence.

impest

Lower-level impulse-response estimator.

TimeSeriesSRC.basefunctions.uniChi.func_uniChi(pmod, y, k=20, alpha=0.05)[source]

Box-Pierce portmanteau chi-square test on one-step-ahead residuals.

Tests the null hypothesis that the residuals of a fitted model are white noise. Large values of the Q statistic (small p-values) indicate that the residuals are autocorrelated and the model does not adequately fit the data.

Parameters:

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

  • y (array-like) – 1-D desired output sequence used to compute residuals.

  • k (int, optional) – Number of residual-ACF lags used in the Q statistic. Default 20.

  • alpha (float, optional) – Significance level for the test (probability of Type I error). Default 0.05.

Returns:

  • passed (int) – 1 if the test is passed (residuals appear white), 0 otherwise.

  • q (float) – Box-Pierce Q statistic.

  • n (int) – Degrees of freedom (k minus the number of free parameters).

  • pval (float) – p-value of the Q statistic under the chi-square distribution.

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.Model.model import pmodel
>>> from TimeSeriesSRC.Model.estimate import estimate
>>> from TimeSeriesSRC.basefunctions.uniChi import func_uniChi
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1, 0.5], [1, -0.8], e)
>>> pm = pmodel('arma', nc=[1], nd=[1], diff=[0], per=[])
>>> pm_est, trec, stat = estimate(pm, y, show_plot=False, show_output=False)
>>> passed, q, n, pval = func_uniChi(pm_est, y)

See also

multiChi

Chi-square test for transfer function models (residual + cross-correlation).

uniAnal

ACF/PACF/GPAC analysis for order identification.

TimeSeriesSRC.basefunctions.multiChi.func_multiChi(pmod, y, u, k1=20, k2=20, alpha1=0.05, alpha2=0.05)[source]

Multivariate chi-square tests for transfer function model validation.

Performs two hypothesis tests on a fitted transfer function model:

  1. Residual whiteness — Box-Pierce Q statistic on the one-step-ahead residuals e. Tests whether the noise model H(q) is adequate.

  2. Residual/input independence — cross-correlation statistic S between the pre-whitened input u and the residuals e. Tests whether the transfer function G(q) is adequate.

Parameters:

  • pmod (pmodel) – Estimated prediction model (ARX, ARMAX, or BJTF).

  • y (array-like) – 1-D desired output sequence.

  • u (array-like) – 1-D or 2-D input sequence, shape (N,) or (1, N).

  • k1 (int, optional) – Residual-ACF lags for the Q statistic. Default 20.

  • k2 (int, optional) – Cross-correlation lags for the S statistic. Default 20.

  • alpha1 (float, optional) – Significance level for the residual whiteness test. Default 0.05.

  • alpha2 (float, optional) – Significance level for the cross-correlation test. Default 0.05.

Returns:

  • pass_arr (list of int) – [pass_Q, pass_S] — 1 if the test passes, 0 otherwise.

  • q (float) – Box-Pierce Q statistic (residual autocorrelation).

  • pvalq (float) – p-value for Q.

  • s (float) – Cross-correlation chi-square statistic.

  • pvals (float) – p-value for S.

  • nq (int) – Degrees of freedom for Q.

  • ns (int) – Degrees of freedom for S.

Examples

>>> import pathlib, pandas as pd
>>> import TimeSeriesSRC as ts
>>> data_dir = pathlib.Path(ts.__file__).parent / 'TestData'
>>> df  = pd.read_csv(data_dir / 'Series_J_Gas_Furnace.csv')
>>> u   = df.iloc[:, 0].values
>>> y   = df.iloc[:, 1].values
>>> pm  = ts.pmodel('arx', na=[2], nb=[2], delay=[3])
>>> pm_est, trec, stat = ts.estimate(pm, y, u=u,
...                                  show_plot=False, show_output=False)
>>> pass_arr, q, pvalq, s, pvals, nq, ns = ts.multiChi(pm_est, y, u)

See also

uniChi

Univariate chi-square test (ARMA models without input).

multiAnal

Impulse response and GPAC analysis for input-output data.

TimeSeriesSRC.basefunctions.xcorr.func_xcorr(a, b, maxlag=20, flag='none')[source]

Compute the cross-correlation (or autocorrelation) between two sequences.

Produces a symmetric array of length 2*maxlag + 1 covering lags -maxlag to +maxlag.

Parameters:

  • a (array-like) – First 1-D sequence (length N).

  • b (array-like) – Second 1-D sequence (same length as a).

  • maxlag (int, optional) – Maximum lag to compute; must be less than N. Default 20.

  • flag ({'biased', 'unbiased', 'coeff', 'none'}, optional) –

    Normalization method:

    • 'biased' — divide by N (biased estimate).

    • 'unbiased' — divide by N - |k| (unbiased estimate).

    • 'coeff' — normalize so that the zero-lag value is 1.0.

    • 'none' — no scaling (raw inner product). Default.

Returns:

c – Cross-correlation values at lags -maxlag, ..., 0, ..., +maxlag. The zero-lag value is at index maxlag.

Return type:

ndarray, shape (1, 2*maxlag+1)

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.xcorr import func_xcorr
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> acf = func_xcorr(y, y, maxlag=20, flag='biased')
>>> acf.shape
(1, 41)

See also

parcor

Partial autocorrelation via Levinson-Durbin.

gpac

Generalized partial autocorrelation table.

TimeSeriesSRC.basefunctions.parcor.func_parcor(acf, nump)[source]

Compute the partial autocorrelation function via the Levinson-Durbin algorithm.

Parameters:

  • acf (array-like) – Full (two-sided) autocorrelation sequence with the zero-lag value at the center — i.e., the output of func_xcorr(). Shape (1, 2*L+1) or (2*L+1,).

  • nump (int) – Number of PACF terms to compute (orders 1 through nump).

Returns:

  • pacf (ndarray, shape (1, nump)) – Partial autocorrelation values at orders 1 through nump.

  • phi (ndarray, shape (nump+1, 1)) – Final AR parameter vector from the Levinson recursion.

  • sigma (float) – Residual variance of the AR(nump) model.

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.xcorr import func_xcorr
>>> from TimeSeriesSRC.basefunctions.parcor import func_parcor
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> acf = func_xcorr(y, y, 20, 'biased')
>>> pacf, phi, sigma = func_parcor(acf, 10)

See also

xcorr

Compute the autocorrelation sequence.

gpac

Generalized partial autocorrelation table.

TimeSeriesSRC.basefunctions.partoacf.func_partoacf(phi, theta, lagmax, var_a)[source]

Compute the theoretical autocovariance function of an ARMA process.

Uses the Yule-Walker method: solves a linear system for the first p lags, then extends the autocovariance sequence by the AR recursion.

The ARMA convention is:

\[y(t) + \phi_1 y(t-1) + \cdots + \phi_p y(t-p) = a(t) + \theta_1 a(t-1) + \cdots + \theta_q a(t-q)\]
Parameters:

  • phi (array-like) – AR polynomial including the leading 1: [1, phi1, phi2, ...].

  • theta (array-like) – MA polynomial including the leading 1: [1, theta1, theta2, ...].

  • lagmax (int) – Number of lags to compute; output covers lags 0 through lagmax - 1.

  • var_a (float) – Variance \(\sigma_a^2\) of the white-noise input a(t).

Returns:

  • acf (ndarray, shape (lagmax,)) – Theoretical autocovariance function. acf[0] is the variance of y.

  • imp (ndarray, shape (p,)) – First p coefficients of the impulse response of \(C(B) / D(B)\).

Examples

>>> from TimeSeriesSRC.basefunctions.partoacf import func_partoacf
>>> acf, imp = func_partoacf([1, 0.8], [1], 10, 1.0)
>>> round(acf[0], 3)   # variance = 1 / (1 - 0.64) ≈ 2.778
2.778

See also

func_partoacf_pmod

Wrapper that reads polynomials from a pmodel.

uniAnal

Computes the sample ACF from data.

TimeSeriesSRC.basefunctions.partoacf.func_partoacf_pmod(pmod, var_a, lagmax)[source]

Compute the theoretical autocovariance from a fitted pmodel.

Extracts the AR (D) and MA (C) polynomials from pmod, composes seasonal factors when present, and calls func_partoacf().

Parameters:

  • pmod (pmodel) – Fitted ARMA prediction model. Seasonal components are handled automatically using pmod.period.

  • var_a (float) – Variance \(\sigma_a^2\) of the white-noise driving process.

  • lagmax (int) – Number of autocovariance lags to compute (lags 0 .. lagmax-1).

Returns:

  • acf (ndarray, shape (lagmax,)) – Theoretical autocovariance function. acf[0] is the variance of y.

  • imp (ndarray, shape (p,)) – First p impulse-response coefficients of \(C(B)/D(B)\).

  • g_ir (ndarray) – Impulse response of the G transfer function (empty for ARMA models).

Examples

>>> import numpy as np
>>> from TimeSeriesSRC.Model.model import pmodel
>>> from TimeSeriesSRC.basefunctions.partoacf import func_partoacf_pmod
>>> pm = pmodel('arma', nc=[1], nd=[1], diff=[0], per=[])
>>> pm.c[0] = np.array([-0.5])
>>> pm.d[0] = np.array([-0.8])
>>> acf, imp, g_ir = func_partoacf_pmod(pm, 1.0, 10)

See also

func_partoacf

Lower-level function that takes explicit polynomials.

TimeSeriesSRC.basefunctions.gpac.func_gpac(acf, nrows, ncols)[source]

Compute the Generalized Partial Autocorrelation (GPAC) table.

Each cell (j, m) of the GPAC is the ratio of two determinants built from the ACF, following the Woodward & Gray (1981) construction. The table simultaneously identifies the MA order (row index + 1) and AR order (column index + 1) of an ARMA process: cells in an MA(q) column pattern or AR(p) row pattern help select model orders.

Parameters:

  • acf (array-like, shape (1, 2*L+1) or (2*L+1,) or (2*L,)) – Autocorrelation sequence. If two-sided (zero lag at center), the positive-lag half is extracted automatically. If one-sided (zero lag first), it is used directly.

  • nrows (int) – Number of GPAC rows (numerator orders 1 .. nrows).

  • ncols (int) – Number of GPAC columns (denominator orders 1 .. ncols).

Returns:

gpac_array – GPAC table. Row j corresponds to MA order j+1; column m corresponds to AR order m+1.

Return type:

ndarray, shape (nrows, ncols)

Notes

If nrows + ncols exceeds the half-length of the ACF, ncols is silently truncated and a UserWarning is issued.

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.xcorr import func_xcorr
>>> from TimeSeriesSRC.basefunctions.gpac import func_gpac
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> acf = func_xcorr(y, y, 25, 'biased')
>>> G = func_gpac(acf, nrows=5, ncols=5)
>>> G.shape
(5, 5)

See also

plotgpac

Visualize the GPAC table.

uniAnal

Computes and plots ACF, PACF, and GPAC together.

TimeSeriesSRC.basefunctions.plotgpac.func_plotgpac(gpac, gtitle='Gpac Array', ax=None)[source]

Visualize a GPAC table as a grid of colored squares.

Each cell’s area encodes the magnitude of the GPAC value and its color encodes the sign: green for positive, red for negative.

Parameters:

  • gpac (ndarray, shape (nrows, ncols)) – GPAC table as returned by func_gpac().

  • gtitle (str, optional) – Plot title. Default 'Gpac Array'.

  • ax (matplotlib.axes.Axes or None, optional) – Axes on which to draw. If None a new figure is created and displayed. Default None.

Returns:

ax – The axes containing the plot.

Return type:

matplotlib.axes.Axes

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.xcorr import func_xcorr
>>> from TimeSeriesSRC.basefunctions.gpac import func_gpac
>>> from TimeSeriesSRC.basefunctions.plotgpac import func_plotgpac
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> acf = func_xcorr(y, y, 25, 'biased')
>>> G = func_gpac(acf, 5, 5)
>>> ax = func_plotgpac(G, 'My GPAC')

See also

gpac

Compute the GPAC table.

uniAnal

Calls this function automatically.

TimeSeriesSRC.basefunctions.plotacf.func_plotacf(acf, maxlag=None, gtitle='Autocorrelation Function', ax=None)[source]

Plot the autocorrelation function as a stem plot.

Parameters:

  • acf (array-like) – Autocorrelation sequence as returned by func_xcorr(), shape (1, 2*maxlag+1) or (2*maxlag+1,).

  • maxlag (int or None, optional) – Maximum lag displayed. If None, inferred from the length of acf as (len(acf) - 1) // 2. Default None.

  • gtitle (str, optional) – Plot title. Default 'Autocorrelation Function'.

  • ax (matplotlib.axes.Axes or None, optional) – Axes on which to draw. A new figure is created when None. Default None.

Returns:

ax – The axes containing the stem plot.

Return type:

matplotlib.axes.Axes

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.xcorr import func_xcorr
>>> from TimeSeriesSRC.basefunctions.plotacf import func_plotacf
>>> e = np.random.default_rng(0).standard_normal(500)
>>> y = lfilter([1], [1, -0.8], e)
>>> acf = func_xcorr(y, y, 20, 'biased')
>>> ax = func_plotacf(acf, maxlag=20)

See also

xcorr

Compute the autocorrelation sequence.

uniAnal

Higher-level function that plots ACF, PACF, and GPAC together.

TimeSeriesSRC.basefunctions.impest.func_impest(u, y, k)[source]

Estimate the impulse response between two sequences (Wiener-Hopf method).

Solves the Wiener-Hopf equations \(R_{uu}\, g = R_{uy}\) for the finite impulse response g of length k+1.

Parameters:

  • u (array-like, shape (1, N)) – Input (exogenous) sequence.

  • y (array-like, shape (1, N)) – Output sequence (same length as u).

  • k (int) – Number of lags of the impulse response to compute; the result covers lags 0 through k.

Returns:

g – Estimated impulse response coefficients g[0], g[1], …, g[k].

Return type:

ndarray, shape (k+1,)

Examples

>>> import numpy as np
>>> from scipy.signal import lfilter
>>> from TimeSeriesSRC.basefunctions.impest import func_impest
>>> rng = np.random.default_rng(0)
>>> u = rng.standard_normal((1, 300))
>>> e = rng.standard_normal((1, 300)) * 0.1
>>> y = lfilter([1], [1, 0.5], u[0]).reshape(1, -1) + e
>>> g = func_impest(u, y, k=10)
>>> g.shape
(11,)

See also

multiAnal

Higher-level multivariate analysis that calls this function.

TimeSeriesSRC.basefunctions.sdiff.func_sdiff(y, d, p=1)[source]

Apply seasonal (or regular) differencing to a time series.

Computes the d-th order difference at period p:

\[\nabla_p^d\, y(t) = \nabla_p^{d-1}\, y(t) - \nabla_p^{d-1}\, y(t-p)\]

For p = 1 this reduces to ordinary differencing (numpy.diff).

Parameters:

  • y (array-like) – 1-D or row-vector input series.

  • d (int) – Number of differences to apply. d = 0 returns y unchanged.

  • p (int, optional) – Seasonal period. Default 1 (ordinary differencing).

Returns:

yd – Differenced series. Each application of the operator reduces the length by p.

Return type:

ndarray, shape (1, N - d*p)

Raises:

Exception – If d * p is larger than or equal to the length of y.

Examples

>>> import numpy as np
>>> from TimeSeriesSRC.basefunctions.sdiff import func_sdiff
>>> y = np.arange(1.0, 21.0)
>>> yd = func_sdiff(y, d=1, p=4)   # seasonal difference at lag 4
>>> yd.shape
(1, 16)

See also

uniAnal

Passes diff / per arguments directly to this function.

TimeSeriesSRC.basefunctions.makerow.func_makerow(y)[source]

Ensure a numpy array has shape (1, N) (row vector).

If the array is 1-D it is reshaped to (1, N). If it is already 2-D with more rows than columns it is transposed. Otherwise it is returned unchanged.

Parameters:

y (ndarray) – Input array of shape (N,), (N, 1), or (1, N).

Returns:

yr – Row-vector form of y.

Return type:

ndarray, shape (1, N)

Examples

>>> import numpy as np
>>> from TimeSeriesSRC.basefunctions.makerow import func_makerow
>>> func_makerow(np.array([1, 2, 3, 4])).shape
(1, 4)
TimeSeriesSRC.basefunctions.chisqrdf.func_chisqrdf(q, n)[source]

Compute the chi-square cumulative distribution function (CDF).

Returns the probability that a chi-square random variable with n degrees of freedom falls in the interval [0, q].

Parameters:

  • q (float) – Chi-square statistic (upper integration limit), must be > 0.

  • n (int) – Degrees of freedom.

Returns:

pr – \(P(X \leq q)\) where \(X \sim \chi^2(n)\).

Return type:

float

Examples

>>> from TimeSeriesSRC.basefunctions.chisqrdf import func_chisqrdf
>>> func_chisqrdf(18.3, 10)   # should be ≈ 0.95
0.9499...

See also

uniChi

Uses this function to convert Q statistic to p-value.