Simulations of time series with known ARMA structure are often needed as part of other simulations or as learning data sets for developing time series analysis skills. The following program generates a time series by using the IML functions NORMAL, ARMACOV, HANKEL, PRODUCT, RATIO, TOEPLITZ, and ROOT.
*ods trace output;
options ls=80 nostimer nodate ps=60;
proc iml;
reset noname;
start armasim(y,n,phi,theta,seed);
/*-----------------------------------------------------------*/
/* IML Module: armasim */
/* Purpose: Simulate n data points from ARMA process */
/* exact covariance method */
/* Arguments: */
/* */
/* Input: n : series length */
/* phi : AR coefficients */
/* theta: MA coefficients */
/* seed : integer seed for normal deviate generator */
/* Output: y: realization of ARMA process */
/* ----------------------------------------------------------*/
p=ncol(phi)-1;
q=ncol(theta)-1;
y=normal(j(1,n+q,seed));
/* Pure MA or white noise */
if p=0 then y=product(theta,y)[,(q+1):(n+q)];
else do; /* Pure AR or ARMA */
/* Get the autocovariance function */
call armacov(gamma,cov,ma,phi,theta,p);
if gamma[1]<0 then
do;
print 'ARMA parameters not stable.';
print 'Execution terminating.';
stop;
end;
/* Form covariance matrix */
gamma=toeplitz(gamma);
/* Generate covariance between initial y and */
/* initial innovations */
if q>0 then
do;
psi=ratio(phi,theta,q);
psi=hankel(psi[,-((-q):(-1))]);
m=max(1,(q-p+1));
psi=psi[-((-q):(-m)),];
if p>q then psi=j(p-q,q,0)//psi;
gamma=(gamma||psi)//(psi`||i(q));
end;
/* Use Cholesky root to get startup values */
gamma=root(gamma);
startup=y[,1:(p+q)]*gamma;
e=y[,(p+q+1):(n+q)];
/* Generate MA part */
if q>0 then
do;
e=startup[,(p+1):(p+q)]||e;
e=product(theta,e)[,(q+1):(n+q-p)];
end;
y=startup[,1:p];
phi1=phi[,-(-(p+1):(-2))]`;
/* Use difference equation to generate */
/* remaining values */
do ii=1 to n-p;
y=y||(e[,ii]-y[,ii:(ii+p-1)]*phi1);
end;
end;
y=y`;
finish armasim; /* ARMASIM */
run armasim(y,10,{1 -0.8},{1 0.5},1234321);
print ,'Simulated Series:', y;
