/* GAUSS programm to illustrate realization vs. process */
/* GAUSS is installed on the first rows in the PC lab   */

/***********GAUSS initializations **********/
library pgraph;
graphset;
_pmcolor = {0, 0, 0,0, 0, 0, 0, 0, 15};
_pcolor  = {0 0 0 0 0};
format /rd 5,3; @only three digits after decimal point@
/***********GAUSS initializations **********/

nobs=1000000;              /* Define here the length of the "time series */


x=rndn(nobs,1);       /* Draw standard normal random numbers */
z=(cumsumc(x));       /* Cumulate the numbers drawn          */
times=seqa(1,1,nobs); /* Generate sequence 1,2...nobs        */

/* I plot the time series */
title("Single Realization of a stochastic process: White noise");
xlabel("time");
xy(times,x); 
 
title("Single Realization of a stochastic process: Random walk");   
xy(times,z);

/* I draw 100 independent "time series" (100 replications) */

x=rndn(nobs,30);
z=(cumsumc(x));  /* Again cumulate each of the 100 series */

/* Plots */
ylabel("Many realizations of the stochastic process");
xy(times,x);
xy(times,z);
stop;
/* sample means and variances computed over replications 

  meanc and stdc compute standard deviation and mean over columns of
  a matrix. We want to compute the means and standard deviation over
  the rows of the matrix, hence the matrix is transposed */

meancx=meanc(x');
stdx=stdc(x');
meancz=meanc(z');
stdz=stdc(z'); 

ylabel("ensemble mean");
title("Ensemble means: White noise process");
xy(times,meancx);
title("Ensemble means: Random walk");
xy(times,meancz);
ylabel("ensemble std");
title("Ensemble standard deviations: White noise process");
xy(times,stdx);
title("Ensemble standard deviations: Random walk");
xy(times,stdz); 
stop;

/* Alternative process: */

title("An alternative process");

x=rndn(nobs,1)/20+0.05;
z=exp(cumsumc(x));
times=seqa(1,1,nobs);
xy(times,x);
xy(times,z);

/* More replications */

x=rndn(nobs,100)/20+0.05;
xy(times,x);

z=exp(cumsumc(x));
xy(times,z);