/* CEV Process (Constant Elasticity of Variance) */
new; cls;
beta=1.5;
mu=0.1;
sig=0.1;
S0=50;
T=1;
n=252;
S=S0|cev(beta,mu,sig,S0,T,n);

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




proc cev(beta,mu,sig,S0,T,n);
   local X0,delt,Z,X,i,S;
   X0=S0^(2-beta);
   delt=T/n;
   Z=rndn(n,1);
   X=zeros(n,1);
   X[1]=X0+(sig^2/2*(2-beta)*(1-beta)+mu*(2-beta)*X0)*delt+sig*(2-beta)*sqrt(X0)*sqrt(delt)*Z[1];
   i=1;
   do while i<=n-1;
      X[i+1]=X[i]+(sig^2/2*(2-beta)*(1-beta)+mu*(2-beta)*X[i])*delt+sig*(2-beta)*sqrt(X[i])*sqrt(delt)*Z[i+1];
      i=i+1;
   endo;
   S=X^(1/(2-beta));
   retp(S);
endp;