%% Work with Lag Operator Polynomials
%% Create Multivariate Lag Operator Polynomial
% Create the 2-D, degree 3 lag operator polynomial
% \[A(L) = \left[ {\begin{array}{*{20}{c}}
% {1}&{0}\\
% { 0}&{1}
% \end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
% {0.5}&{0.25}\\
% { - 0.1}&{0.4}
% \end{array}} \right]L + \left[ {\begin{array}{*{20}{c}}
% {0.05}&{0.025}\\
% { - 0.01}&{0.04}
% \end{array}} \right]{L^3}.\]

m = 2;
A0 = eye(m);
A1 = [0.5 0.25; -0.1 0.4];
A3 = 0.1*A1;
Coeffs = {A0 A1 A3};
lags = [0 1 3];

A = LagOp(Coeffs,'Lags',lags)

%
%% Determine Polynomial Stability
% A lag operator polynomial is stable if the magnitudes of all eigenvalues of its characteristic polynomial are less than 1.
% Determine whether the polynomial is stable, and return the eigenvalues of its characteristic polynomial.

[tf,evals] = isStable(A)

%
%% Invert Polynomial
% Compute the inverse of A\left(L\right) by right dividing the 2-by-2 identity matrix by A\left(L\right).

Ainv = mrdivide(eye(A.Dimension),A)

% Ainv is a LagOp object representing the inverse of A\left(L\right), a degree 9 lag operator polynomial. In theory, the inverse of a lag operator polynomial is of infinite degree, but the mrdivide coefficient-magnitude tolerances truncate the polynomial.
% Multiply A\left(L\right) and its inverse.

checkinv = Ainv*A

% Because the inverse computation returns the truncated theoretical inverse, the product checkinv contains lags representing the remainder.  You can decrease the mrdivide coefficient-magnitude tolerances to obtain an inverse polynomial that is more precise.
%
%% Filter Time Series
% Produce the 2-D time series y_t=A(L)e_t=A_0e_t+A_1e_{t-1}+A_2e_{t-2}+A_3e_{t-3} by filtering the 2-D series of 100 random standard Gaussian deviates e_t \;through the polynomial.

T = 100;
e = randn(T,m);
y = filter(A,e);
plot((A.Degree + 1):T,y)
title('Filtered Series')

% y is a 97-by-2 matrix representing y_t \ldotp y has p fewer observations than e because filter requires the first p observations of e to initialize the dynamic series when producing y(1:4,:).
% Copyright 2020 The MathWorks, Inc.