For example, I just have 2 variables, a, b, in the model and one exogenous shock e_m, which is iid standard normal.

a=0.6*a(-1)+0.4*a(+1)+b;

b=0.88*b(-1)+e_m;

it is easy to compute irf for variable “b”, but how irf is computed for variable “a” since there are lead and lag in the first equations?