for impulse response functions (this is only non-orthogonal right now).

Here is some quick code to get VAR IRF MC standard errors:

import numpy as np import scikits.statsmodels.api as sm from scikits.statsmodels.tsa.api import VAR import matplotlib.pyplot as plt mdata = sm.datasets.macrodata.load().data mdata = mdata[['realgdp','realcons','realinv']] #mdata = mdata[['realgdp','realcons','realinv', 'pop', 'realdpi', 'tbilrate']] names = mdata.dtype.names data = mdata.view((float,3)) data = np.diff(np.log(data), axis=0) mymodel = VAR(data,names=names) res = mymodel.fit(maxlags=3,ic=None) #first generate asymptotic standard errors (to compare) res.irf(periods=30).plot(orth=False, stderr_type='asym') #then generate monte carlo standard errors res.irf(periods=30).plot(orth=False, stderr_type='mc', seed=29, repl=1000) plt.show()

Produces the following plots of a shock to (delta) realgdp.

Asymptotic:

Monte Carlo:

I added functions to the VARResults, IRAnalysis, and also added to

some of the pre-existing functions in these classes and also to

plotting.py. Because Monte Carlo standard errors are in general not

symmetric, I had to alter the plot_with_error function from the

tsa.vector_ar.plotting.py file and number of other functions.

Functions added:

VARResults.stderr_MC_irf(self, orth=False, repl=1000, T=25,

signif=0.05, seed=None)

This function generates a tuple that holds the lower and upper error

bands generated by Monte Carlo simulations.

IRAnalysis.cov_mc(self, orth=False, repl=1000, signif=0.05, seed=None)

This just calls the stderr_MC_irf function from the original model

using the number of periods specified when the irf() class is called

from the VAR class.

Modified functions:

BaseIRAnalysis.plot

Added specification of error type stderr_type='asym' or 'mc' (see

example). Also repl (replications) and seed options added if 'mc'

errors are specified.

in tsa.vector_ar.plotting.py

plot_with_error()

irf_grid_plot()

These functions now take the error type specified in the plot()

function above and treats the errors accordingly. Because of how all

of the VAR work (esp. plotting) assumed asymptotic standard errors,

all of these irf plot functions assumed that errors were passed in as

a single matrix with each standard error depending on the MA lag

length and shock variable. Now, depending on which error is specified,

the function will take a tuple of arrays as the standard error if

(stderr_type='mc') rather than a single array.

A serious issue right now is speed. While the asymptotic standard

errors take about a half second on my home laptop to run, the Monte

Carlo standard errors with 3 variables and 1000 replications takes

about 13 seconds to run. Each replication discards the first 100

observations. The most taxing aspect of generating the errors is

re-simulating the data assuming normally distributed standard errors

1000 times (using the util.varsim() function).

Bart