def varsim(coefs, intercept, sig_u, steps=100, initvalues=None, seed=None):
 rs = np.random.RandomState(seed=seed)
 rmvnorm = rs.multivariate_normal
 (p, k, k) = coefs.shape
 ugen = rmvnorm(np.zeros(len(sig_u)), sig_u, steps)
 result = np.zeros((steps, k))
 result[p:] = (intercept + ugen[p:])
 for t in range(p, steps):
  ygen = result[t]
  for j in range(p):
   ygen += np.dot(coefs[j], result[((t - j) - 1)])
 return result

def varsim ( coefs , intercept , sig_u , steps = 100 , initvalues = None , seed = None ) : rs = np . random . RandomState ( seed = seed ) rmvnorm = rs . multivariate_normal ( p , k , k ) = coefs . shape ugen = rmvnorm ( np . zeros ( len ( sig_u ) ) , sig_u , steps ) result = np . zeros ( ( steps , k ) ) result [ p : ] = ( intercept + ugen [ p : ] ) for t in range ( p , steps ) : ygen = result [ t ] for j in range ( p ) : ygen += np . dot ( coefs [ j ] , result [ ( ( t - j ) - 1 ) ] ) return result
simulate simple var(p) process with known coefficients .