/* Mean Reversion Process */
new; cls;
eta=0.5;
sig=0.1;
xbar=50;
X0=50;
T=1;
n=252;
X=X0|meanrev(eta,sig,xbar,X0,T,n);

library pgraph;
graphset;
xy(seqa(0,T/n,n+1),X);




proc meanrev(eta,sig,xbar,X0,T,n);
   local delt,Z,X,i;
   delt=T/n;
   Z=rndn(n,1);
   X=zeros(n,1);
   X[1]=X0*exp(-eta*delt)+xbar*(1-exp(-eta*delt))+sig*sqrt((1-exp(-2*eta*delt))/(2*eta))*Z[1];
   i=1;
   do while i<=n-1;
      if X[i]<0;
         X[i]=X0/100;   /* for example to avoid negtive X[i] */
      endif;
      X[i+1]=X[i]*exp(-eta*delt)+xbar*(1-exp(-eta*delt))+sig*sqrt((1-exp(-2*eta*delt))/(2*eta))*Z[i+1];
      i=i+1;
   endo;
   retp(X);
endp;