/*  This segment of the program generates a random walk without drift  */

data mc1;

x1 = 0;
iseed = 8888;

do t = 0 to 100;
  a = rannor(iseed);

x2 = 1.0*x1 + a;

if t > 0 then output;

x1 = x2;

end;

keep t x2;
run;

data mc1;
  set mc1;
  y = x2;
  time = _N_;
  drop t x2;
  run;

proc gplot data=mc1;
  symbol v=dot c=black i=join h=.8;
  title1 'Monte Carlo Random Walk Without Drift Data';
  title2 'X=Time Y=RW Series';
  axis1 order=(0 to 100 by 10)
  label=(f=duplex 'Obs');
  axis2 order=(-10 to 10 by 2.0)
  label=(f=duplex 'RW Series');
  plot y*time / haxis=axis1 vaxis=axis2 overlay vref = 0;
  run;

  title;

proc arima data=mc1;
  identify var=y;
  run;

data mc1;
  set mc1;
  ydif = dif(y);
  run;

proc gplot data=mc1;
  symbol v=dot c=black i=join h=.8;
  title1 'Differenced RW Without Drift';
  title2 'X=Time Y=Diff. RW Series';
  axis1 order=(0 to 100 by 10)
  label=(f=duplex 'Obs');
  axis2 order=(-5 to 5 by 1.0)
  label=(f=duplex 'Dif Series');
  plot ydif*time / haxis=axis1 vaxis=axis2;
  run;

  title;

proc arima data=mc1;
  identify var=y(1);
  run;

/*  This segment of the program generates a Random Walk with Drift  */

data mc2;

x1 = 0;
iseed = 445;

do t = 0 to 100;
  a = rannor(iseed);

x2 = 1.0 + 1.0*x1 + 2*a;

if t > 0 then output;

x1 = x2;

end;

keep t x2;
run;

data mc2;
  set mc2;
  time = _N_;
  y = x2;
  drop x2 t;
  run;

title 'Random Walk with Drift Data';
proc sgplot data=mc2 noautolegend;
   reg x=time y=y/lineattrs=(color=black) markerattrs=(size=5.0);
   run;

  title;

proc arima data=mc2;
  identify var=y;
  run;

data mc2;
  set mc2;
  ydif = dif(y);
  run;

proc gplot data=mc2;
  symbol v=dot c=black i=join h=.8;
  title1 'Differenced RW With Drift';
  title2 'X=Time Y=Diff. RW Series';
  axis1 order=(0 to 100 by 10)
  label=(f=duplex 'Obs');
  axis2 order=(-10 to 10 by 1.0)
  label=(f=duplex 'Dif Series');
  plot ydif*time / haxis=axis1 vaxis=axis2 overlay vref = 0;
  run;

proc arima data=mc2;
  identify var=y(1);
  run;

proc reg data=mc2;
  model y = time;
  output out=resid r=resid;
  run;

proc gplot data=resid;
  symbol v=dot c=black i=join h=.8;
  title1 'Deterministic Trend Residuals for RW with Drift Data';
  title2 'X = Time Y = Residuals from Deterministic Trend';
  axis1 order=(0 to 100 by 10)
  label=(f=duplex 'OBS');
  axis2 order=(-40 to 40 by 10)
  label=(f=duplex 'Residuals');
  plot resid*time / haxis=axis1 vaxis=axis2 overlay vref = 0;;

proc arima data=resid;
   identify var=resid;
   run;