Mvarmodel

MVAR class

Tools for MVAR parameters fitting

class connectivipy.mvarmodel.Mvar[source]

Static class Mvar to multivariete autoregressive model fitting. Possible methods are in fitting_algorithms where key is acronym of algorithm and value is a function from mvar.fitting.

classmethod fit(data, order=None, method='yw')[source]

Mvar model fitting. Args:

data : numpy.array
array with data shaped (k, N), k - channels nr, N-data points)
order = None : int
model order, when default None it estimates order using akaike order criteria.
method = ‘yw’: str
name of mvar fitting algorithm, default Yule-Walker all avaiable methods you can find in fitting_algorithms
Returns:
Av : numpy.array
model coefficients (kXkXorder)
Vf : numpy.array
reflection matrix (kXk)
classmethod order_akaike(data, p_max=None, method='yw')[source]

Akaike criterion of MVAR order estimation.

Args:
data : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
p_max = 5 : int
maximal model order
method = ‘yw’ : str
name of the mvar calculation method
Returns:
best_order : int
minimum of crit array
crit : numpy.array
order criterion values for each value of order p starting from 1

References: .. [1] Blinowska K. J., Zygierewicz J., (2012) Practical

biomedical signal analysis using MATLAB. Boca Raton: Taylor & Francis.
classmethod order_hq(data, p_max=None, method='yw')[source]

Hannan-Quin criterion of MVAR order estimation.

Args:
data : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
p_max = 5 : int
maximal model order
method = ‘yw’ : str
name of the mvar calculation method
Returns:
best_order : int
minimum of crit array
crit : numpy.array
order criterion values for each value of order p starting from 1

References: .. [1] Blinowska K. J., Zygierewicz J., (2012) Practical

biomedical signal analysis using MATLAB. Boca Raton: Taylor & Francis.
classmethod order_schwartz(data, p_max=None, method='yw')[source]

Schwartz criterion of MVAR order estimation.

Args:
data : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
p_max = 5 : int
maximal model order
method = ‘yw’ : str
name of the mvar calculation method
Returns:
best_order : int
minimum of crit array
crit : numpy.array
order criterion values for each value of order p starting from 1

References: .. [1] Blinowska K. J., Zygierewicz J., (2012) Practical

biomedical signal analysis using MATLAB. Boca Raton: Taylor & Francis.
classmethod order_fpe(data, p_max=None, method='yw')[source]

Final Prediction Error criterion of MVAR order estimation. (not recommended) Args:

data : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
p_max = 5 : int
maximal model order
method = ‘yw’ : str
name of the mvar calculation method
Returns:
best_order : int
minimum of crit array
crit : numpy.array
order criterion values for each value of order p starting from 1

References: .. [1] Akaike H, (1970), Statistical predictor identification,

Ann. Inst. Statist. Math., 22 203–217.

Algorithms

connectivipy.mvar.fitting.mvar_gen(Acf, npoints, omit=500)[source]

Generating data point from MVAR coefficiencs matrix Acf. Args:

Acf : numpy.array
array in shape of (p,k,k) where k is number of channels and p is a model order.
npoints : int
number of data points.
Returns:
y : numpy.array
(k, npoints) data points
connectivipy.mvar.fitting.mvar_gen_inst(Acf, npoints, omit=500)[source]

Generating data point from matrix A with MVAR coefficients but it takes into account also zerolag interactions. So here Acf[0] means instantenous interaction not as in mvar_gen one data point lagged. Args:

Acf : numpy.array
array in shape of (p,k,k) where k is number of channels and p is a model order.
npoints : int
number of data points.
Returns:
y : np.array
(k, n) data points
connectivipy.mvar.fitting.stability_mvar(Acf)[source]

Checks stability of MVAR given its parameters matrix Acf.

! Be careful when using ! not tested yet.

Implemented according to: https://sccn.ucsd.edu/wiki/Chapter_3.1._Stationarity_and_Stability

Args:
Acf : numpy.array
array in shape of (p,k,k) where k is number of channels and p is a model order.
Returns:
stable : bool
stability flag - when True it is stable.
connectivipy.mvar.fitting.meanncov(x, y=[], p=0, norm=True)[source]

Wrapper to multichannel case of new covariance ncov. Args:

x : numpy.array
multidimensional data (channels, data points, trials).
y = [] : numpy.array
multidimensional data. If not given the autocovariance of x will be calculated.
p = 0: int
window shift of input data. It can be negative as well.
norm = True: bool
normalization - if True the result is divided by length of x, otherwise it is not.
Returns:
mcov : np.array
covariance matrix
connectivipy.mvar.fitting.ncov(x, y=[], p=0, norm=True)[source]

New covariance. Args:

x : numpy.array
onedimensional data.
y = [] : numpy.array
onedimensional data. If not given the autocovariance of x will be calculated.
p = 0: int
window shift of input data. It can be negative as well.
norm = True: bool
normalization - if True the result is divided by length of x, otherwise it is not.
Returns:
kv : np.array
covariance matrix
connectivipy.mvar.fitting.vieiramorf(y, pmax=1)[source]

Compute multichannel autoregresive model coefficients using Vieira-Morf algorithm. Args:

y : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
pmax: int >0
model order
Returns:
Ar : np.array
matrix with parameters matrix (p, k, k) where p - model order, k - nr of channels
Vr : np.array
prediction error matrix (k,k)

References: .. [1] Marple, Jr S. L., *Digital Spectral Analysis with

Applications*, Prentice-Hall Signal Processing Series, 1987
connectivipy.mvar.fitting.nutallstrand(y, pmax=1)[source]

Compute multichannel autoregresive model coefficients using Nutall-Strand algorithm. Args:

y : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
pmax: int >0
model order
Returns:
Ar : np.array
matrix with parameters matrix (p, k, k) where p - model order, k - nr of channels
Vr : np.array
prediction error matrix (k,k)

References: .. [1] Marple, Jr S. L., *Digital Spectral Analysis with

Applications*, Prentice-Hall Signal Processing Series, 1987
connectivipy.mvar.fitting.yulewalker(y, pmax=1)[source]

Compute multichannel autoregresive model coefficients using Yule-Walker algorithm. Args:

y : numpy.array
multichannel data in shape (k, n) for one trial case and (k, n, tr) for multitrial k - nr of channels, n -data points, tr - nr of trials
pmax: int >0
model order
Returns:
Ar : np.array
matrix with parameters matrix (p, k, k) where p - model order, k - nr of channels
Vr : np.array
prediction error matrix (k,k)

References: .. [1] Marple, Jr S. L., *Digital Spectral Analysis with

Applications*, Prentice-Hall Signal Processing Series, 1987

Additional

Additional tools

connectivipy.mvar.comp.ldl(A)[source]

LDL decomposition (implementation from en.wikipedia.org/wiki/Cholesky_decomposition) Args:

A : numpy.array
matrix kXk
Returns:
L, D, LT : np.array
L is a lower unit triangular matrix, D is a diagonal matrix and LT is a transpose of L.