## Autoregressive processes : learning with an interactive tool

I recently got interested in second order autoregressive processes ("AR(2)" processes), after spending quite some time working only with first order ones... It was time to get an extra degree of freedom/complexity ! However, I didn't have much intuition about the behavior of these models and neither about how to "tune the knobs" that set this behavior so I spent some time on this subject.

For the sake of clarity, let's write down a definition/notation for an AR(2) process $X_k$: $X_k = a_1 X_{k-1} + a_2 X_{k-2} + Z_k$

with the real numbers a1 and a2 being the two parameters that tune the "behavior" of the process. $Z_k$ may be call the "input" in the context of linear filters, while "innovation" (i.e. a white noise input) would better fit in the context of random processes. A nice property of AR(2) processes, as opposed to AR(1) is their ability to capture slowly oscillating behaviors, like ocean waves. Also, within the AR(2) model, AR(1) is included as the special case a2 = 0.

## Bits of a non-interactive study

In order to understand the interplay of the two parameters a1 and a2, I spent extensive time doing paper-based analytical/theoretical studies. I also created a IPython Notebook which I should publish when there is proper Notebook -> HTML conversion system (edit: now available on nbviewer). Most of the AR(2) behavior is governed by the characteristic polynomial $P(z) = z^2 - a_1 z - a_2$. This polynomial helps in showing that there are two regimes : the one with two real roots (when $a_2 > - a_1^2/4$) and the one with two complex-conjugate roots (when $a_2 < - a_1^2/4$) where I can get my long-awaited slow oscillations. A third regime is the special case of a double root $a_2 = - a_1^2/4$ which I don't think is of a special interest.

Beyond the study of the root, I also computed the variance of the process and its spectral density. I found both not so easy to compute !

## Creating an interactive tool

After the static paper-based study and the still semi-static computer-based study, where I was able to derive and plot the impulse response and some realization of random trajectories, I wanted to get a real-time interactive experience of the interplay of the two parameters a1 and a2. To this end I wrote a tiny Python script, based on the event handling machinery of Matplotlib. More specifically, I started from the "looking glass" example script which I altered progressively as my understanding of event handling in Matplotlib was getting better. Here is a frozen screen capture of the resulting plots:

1. the root locus diagram, where one can see the transition between the two real roots and the two complex-conjugate roots as the $a_2 = - a_1^2/4$ boundary is crossed.
2. the Impulse Response (output of the filter for a discrete impulse $Z_k = \delta_k$ ), implemented simply using scipy.signal.lfilter.