Trouble using aesara scan to implement error in AR(1)

Questions v5
,
1

Hello all,

I was trying to model the error in an AR(1) process but was having trouble with the mechanics of implementing the scan function.

I am using aesara version 2.6.6. I have also tried it on google colabs and I get the same error in colabs.

I generated some fake data and tried to use aesara scan to implement the error

Here is my sample code:

#Code to look at how to model AR(1) in PyMC

# AR(1) Process:

# y_t = mu + phi * (y_t-1) + eps_t

# Where eps_t = Normal(0, 1)

# Import libraries
import numpy as np
import aesara
import aesara.tensor as at
import pymc as pm
import matplotlib.pyplot as plt
import arviz as az

# Set style
az.style.use('arviz-darkgrid')

# Set Random Seed
RANDOM_SEED = 1234

# Generate the data

# Inputs
Npts = 100
mu_sim = 10.0 # sim stands for simulation
phi_sim = 0.55 # Autoregressive coefficient

# Initialize the data
y_sim = np.zeros(Npts) 
y_sim[0] = mu_sim

# Generate the fake noise
np.random.seed(RANDOM_SEED)
eps_sim = np.random.normal(0, 1, size = Npts - 1) # Create size of N - 1
eps_sim = np.concatenate([[0], eps_sim]) # Then concatenate [0] with the above step

for i in range(1, Npts):
    y_sim[i] = mu_sim + phi_sim * y_sim[i - 1] + eps_sim[i]

# Calculate the error using aesara scan
def ar_1_error(y_tm1, y_t, mu, phi):
    '''
    ----------------------------------------
    Calculate the error from a AR(1) process
    error_t = y_t - mu - phi * y_tm1
    ---------------------------------------
    Inputs:
      y_tm1: lagged 1 value of time series
      y_t: current time series value
      mu: mean
      phi: AR(1) coefficient
    Returns:
      error_t = y_t - mu - phi * y_tm1
    '''

    error_t = y_t - mu - phi * y_tm1

    return error_t

err, _ = aesara.scan(
    fn = ar_1_error,
    sequences = dict(input = y_sim, taps = [-1, 0]),
    non_sequences = [mu_sim, phi_sim]
)

I get the following error:
--------------------------------------------------------------------------- TypeError Traceback (most recent call last) c:\Users\Prabh\OneDrive\Desktop\error_AR.ipynb Cell 1 in <cell line: 120>() 116 error_t = y_t - mu - phi * y_tm1 118 return error_t → 120 err, _ = aesara.scan( 121 fn = ar_1_error, 122 sequences = dict(input = y_sim, taps = [-1, 0]), 123 non_sequences = [mu_sim, phi_sim] 124 ) File c:\Users\Prabh\Anaconda3\envs\intuitive_bayes\lib\site-packages\aesara\scan\basic.py:681, in scan**(fn, sequences, outputs_info, non_sequences, n_steps, truncate_gradient, go_backwards, mode, name, profile, allow_gc, strict, return_list)** 678 else: 679 actual_n_steps = at.as_tensor(n_steps) → 681 scan_seqs = [seq[:actual_n_steps] for seq in scan_seqs] 682 # Conventions : 683 # mit_mot = multiple input taps, multiple output taps ( only provided 684 # by the gradient function ) (…) 688 689 # MIT_MOT – not provided by the user only by the grad function 690 n_mit_mot = 0 File c:\Users\Prabh\Anaconda3\envs\intuitive_bayes\lib\site-packages\aesara\scan\basic.py:681, in (.0) 678 else:

…

688
689 # MIT_MOT – not provided by the user only by the grad function
690 n_mit_mot = 0

TypeError: slice indices must be integers or None or have an index method

Thanks!

2

I implemented the AR generative process (not just the errors) here in case it helps. The focus of that notebook is a bit more dense though: AR PyMC Aeppl.ipynb · GitHub

1 Like
3

Thanks! I will take a look.

It seems that it may be easier to model an AR(1) process as
y[t] ~ Normal(y_tm1*rho, sigma) rather than calculate the error as I did and model that at Normal(0, sigma).

I think this give me a good guide to get started.

Thanks again!

4

Any way to generalize to a ARMA (1, 1)?

This was my attempt, it gives the same error (coding error not error term in time series) as before.
Is there an example notebook for ARMA(1, 1) using aesara scan?

I will still keep the reply Ricardo gave as the solution, this is only for my curiosity.

# Calculate the error using aesara scan
def arma_1_error(y_tm1, y_t, error_tm1, mu, phi, gamma):
    '''
    ----------------------------------------
    Calculate the error from a ARMA(1, 1) process
    error_t = y_t - mu - phi * y_tm1 - gamma * error_tm1
    ---------------------------------------
    Inputs:
      y_tm1: lagged 1 value of time series
      y_t: current time series value
      error_tm1: lagged prior error
      mu: mean
      phi: AR(1) coefficient
      gamma: MA(1) coefficient
    Returns:
      error_t = y_t - mu - phi * y_tm1 - gamma * error_tm1
    '''

    error_t = y_t - mu - phi * y_tm1 - gamma * error_tm1

    return error_t

gamma_sim = 0.000001
init_value = y_sim[0] - mu_sim - phi_sim * mu_sim
err, _ = aesara.scan(
    fn = arma_1_error,
    sequences = dict(input = y_sim, taps = [-1, 0]),
    outputs_info = dict(initial = init_value, taps = [-1]),
    non_sequences = [mu_sim, phi_sim, gamma_sim]
)
5

Maybe @junpenglao has seen such an example

6

I do but only in TensorFlow: 6. Time Series — Bayesian Modeling and Computation in Python
(@Prabh’s implementation looks pretty good)