/* This program generates some monte carlo data for use in the forecasting class. The first series is a stationary AR(1) process. The rest of the series are nonstationary. */ options pagesize=60 linesize=74 nodate; goptions device=win gsfname=plot rotate=landscape gsfmode=append target=hpljs3; *border; filename Plot 'c:\sas\work\mcarlo.plt'; /* This segment of the program generates the stationary AR(1) model */ data mc1; x1 = 0; iseed = 12345; do t = 0 to 100; a = rannor(iseed); x2 = 0.5*x1 + a; if t > 0 then output; x1 = x2; end; keep t x2; title ' Monte Carlo AR(1) data with phi(1) = 0.5'; proc print data = mc1; var t x2; proc plot data = mc1; plot x2*t; proc gplot data=mc1; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo AR(1) data with phi(1) = 0.5'; title2 'X=Time Y=AR(1) Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-5 to 5 by 1.0) label=(f=duplex 'AR(1) Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc1; identify var=x2; /* This segment of the program generates a random walk without drift */ data mc2; 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; title 'Monte Carlo Random Walk Without Drift Data'; proc print data = mc2; var t x2; proc plot data = mc2; plot x2*t; proc gplot data=mc2; 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 x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc2; identify var=x2; data mc2; set mc2; x2dif = dif(x2); title 'Differenced Monte Carlo RW Without Drift data'; proc print data = mc2; var t x2dif; proc plot data = mc2; plot x2dif*t; proc gplot data=mc2; 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 'AR(1) Series'); plot x2dif*t / haxis=axis1 vaxis=axis2; proc arima data=mc2; identify var=x2(1); /* This segment of the program generates a Random Walk with Drift */ data mc3; x1 = 0; iseed = 2468; 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; title 'Monte Carlo Random Walk With Drift Data'; proc print data = mc3; var t x2; proc plot data = mc3; plot x2*t; proc gplot data=mc3; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo Random Walk With Drift Data'; title2 'X=Time Y=RW Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-10 to 120 by 10) label=(f=duplex 'RW Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc3; identify var=x2; data mc3; set mc3; x2dif = dif(x2); title 'Differenced Monte Carlo Random Walk With Drift Data'; proc print data = mc3; var t x2dif; proc plot data = mc3; plot x2dif*t; proc gplot data=mc3; 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 'AR(1) Series'); plot x2dif*t / haxis=axis1 vaxis=axis2; proc arima data=mc3; identify var=x2(1); /* This segment of the program shows the effect of incorrectly detrending the Random Walk with Drift data using a time trend. */ title 'Incorrectly Detrended Random Walk with Drift Data'; data mc3; set mc3; t = _n_; proc reg data=mc3; model x2 = t; output out=resid r=resid; 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 from Trend'); plot resid*t / haxis=axis1 vaxis=axis2; proc arima data=resid; identify var=resid; /* This segment of the program generates a Percentage Growth model */ data mc4; x1 = 0; iseed = 333; do t = 0 to 100; a = rannor(iseed); x2 = .03 + 1.0*x1 + a*(.03); if t > 0 then output; x1 = x2; end; keep t x2; data mc4; set mc4; y = exp(x2); title 'Monte Carlo Growth Data'; proc print data = mc4; var t y; proc plot data = mc4; plot y*t; proc gplot data=mc4; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo Growth Data'; title2 'X=Time Y=Growth Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(0 to 25 by 5) label=(f=duplex 'Growth Series'); plot y*t / haxis=axis1 vaxis=axis2; proc arima data=mc4; identify var=y; data mc4; set mc4; ydif = dif(y); title 'Difference in Monte Carlo Growth Data'; proc print data = mc4; var t ydif; proc plot data = mc4; plot ydif*t; proc gplot data=mc4; symbol v=dot c=black i=join h=.8; title1 'Difference of Monte Carlo Growth Data'; title2 'X=Time Y=Diff. of Growth Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-2.0 to 2.0 by .2) label=(f=duplex 'Diff. of Growth Series'); plot ydif*t / haxis=axis1 vaxis=axis2; data mc4; set mc4; x2dif = dif(x2); title 'Percentage Changes in Growth Data'; proc print data = mc4; var t x2dif; proc plot data = mc4; plot x2dif*t; proc gplot data=mc4; symbol v=dot c=black i=join h=.8; title1 'Differenced Growth Series = Percentage Change'; title2 'X=Time Y=Diff. Growth Series = Percentage Change'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-.15 to .20 by .05) label=(f=duplex 'Perc. Change'); plot x2dif*t / haxis=axis1 vaxis=axis2; proc arima data=mc4; identify var=x2(1); /* This segment of the program generates a Deterministic Trend model */ data mc5; x1 = 0; iseed = 7654321; do t = 0 to 100; a = rannor(iseed); x2 = 10.0 + 4.0*t + 10*a; if t > 0 then output; end; keep t x2; title 'Deterministic Trend Data'; proc print data = mc5; var t x2; proc plot data = mc5; plot x2*t; proc gplot data=mc5; symbol v=dot c=black i=join h=.8; title1 'Deterministic Trend Data'; title2 'X=Time Y=Trend Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(0 to 450 by 25) label=(f=duplex 'Trend Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc5; identify var=x2; proc reg data=mc5; model x2 = t; output out=resid r=resid; data mc5; set resid; title 'Residuals from Deterministic Trend'; proc print data = mc5; var t resid; proc plot data = mc5; plot resid*t; proc gplot data=mc5; symbol v=dot c=black i=join h=.8; title1 'Residuals from Deterministic Trend'; title2 'X=Time Y=Residuals from Deterministic Trend'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-80 to 80 by 10) label=(f=duplex 'Residuals from Trend'); plot resid*t / haxis=axis1 vaxis=axis2; proc arima data=mc5; identify var=resid; /* This segment of the program shows the effect of incorrectly differencing deterministic trend. */ data mc5; set mc5; x2dif = dif(x2); proc gplot data=mc5; symbol v=dot c=black i=join h=.8; title1 'Differencing Deterministic Trend Data'; title2 'X=Time Y=Differenced Deterministic Trend Data'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-80 to 80 by 10) label=(f=duplex 'Differenced Deterministic Trend Data'); plot x2dif*t / haxis=axis1 vaxis=axis2; proc arima data=mc5; identify var=x2dif; /* This segment of the program generates the Covariance Switching model */ data mc6; x1 = 0; iseed = 123; do t = 0 to 50; a = rannor(iseed); x2 = 0.5*x1 + a; if t > 0 then output; x1 = x2; end; x1 = 0; iseed = 321; do t = 50 to 100; a = rannor(iseed); x2 = -0.5*x1 + a; if t > 50 then output; x1 = x2; end; keep t x2; title ' Covariance Switching Data'; proc print data = mc6; var t x2; proc plot data = mc6; plot x2*t; proc gplot data=mc6; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo Covariance Switching Data'; title2 'X=Time Y=Cov. Switching Series'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-6 to 6 by 1) label=(f=duplex 'Cov. Switching Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc6; identify var=x2; /* This segment of the program generates the Variance Change model */ data mc7; iseed = 111; do t = 1 to 50; a = rannor(iseed); x2 = a; output; end; iseed = 777; do t = 51 to 100; a = rannor(iseed); x2 = 3*a; output; end; keep t x2; title ' Variance Change Data'; proc print data = mc7; var t x2; proc plot data = mc7; plot x2*t; proc gplot data=mc7; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo Variance Change Data'; title2 'X=Time Y=Variance Change Data'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-10 to 10 by 2) label=(f=duplex 'Var. Change Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc7; identify var=x2; /* This segment of the program generates the Changing Mean model */ data mc8; iseed = 11111; do t = 1 to 50; a = rannor(iseed); x2 = a; output; end; iseed = 77777; do t = 51 to 100; a = rannor(iseed); x2 = 5 + a; output; end; keep t x2; title ' Changing Mean Data'; proc print data = mc8; var t x2; proc plot data = mc8; plot x2*t; proc gplot data=mc8; symbol v=dot c=black i=join h=.8; title1 'Monte Carlo Changing Mean Data'; title2 'X=Time Y=Changing Mean Data'; axis1 order=(0 to 100 by 10) label=(f=duplex 'Obs'); axis2 order=(-10 to 10 by 2) label=(f=duplex 'Changing Mean Series'); plot x2*t / haxis=axis1 vaxis=axis2; proc arima data=mc8; identify var=x2; run;