proc (1)=forecast(ar,ma,res,armaser,steps);
/*
forecast.src is written by Kristian Jönsson (April 10, 2003), Department of Economics, Lund University.
Contact info:  kristian.jonsson@nek.lu.se

This procedure makes forecasts "#steps" steps ahead, of a specified ARMA process, standing at time t.
The first forecast is made at time max(p,q) which ensures the availability of initial vaules.

Inputs: 
  ar		: A px1 vector with AR coefficients.
  ma	 	: A qx1 vector with MA coefficients.
  res		: A Tx1 vector with estimated residuals for the fitted ARMA-process.
  armaser	: A Tx1 vector with the ARMA-process itself.
  steps		: The number of forecasts to make at each time period (scalar).

Output:
  forecast	: A Txsteps matrix with the forecasts made at each time period.

Comment:
  The data vectors (res and armaser) must be scaled so that they are of equal length.

Note:
This code can be used freely as long as proper reference is given.
No performance guarantee is made.
Bugreports are welcome.
*/

local p,q;
local maxpq;
local forecast;
local y_ser,res_ser;
local i,j;

p=rows(ar);
q=rows(ma);
maxpq=maxc(p|q);

forecast=zeros(rows(armaser),steps);

for i (maxpq,rows(armaser),1); @Starts at a point that ensures initial value availability@

y_ser=armaser[i:i-p+1,1];
res_ser=res[i:i-q+1,1];
forecast[i,1]=ar'*y_ser+ma'*res_ser;

y_ser=forecast[i,1]|y_ser;
y_ser=trimr(y_ser,0,1);
res_ser=0|res_ser;
res_ser=trimr(res_ser,0,1);

if steps>1;  @Forecasts for steps>1@
  for j (2,steps,1);
   forecast[i,j]=ar'*y_ser+ma'*res_ser;
   y_ser=forecast[i,j]|y_ser;
   y_ser=trimr(y_ser,0,1);
   res_ser=0|res_ser;
   res_ser=trimr(res_ser,0,1);
  endfor;
endif;

endfor;


retp(forecast);
endp;