[Rd] Problems Simulating (PR#3471)

[smiles at dept.econ.yorku.ca]

Full_Name: Stan Miles
Version: 1.5.0.
OS: Windows 2000 Server
Submission from: (NULL) (130.63.74.220)


Short Description of the Problem:  

I use R to simulate a random variable, with mean 0.14 and variance 0.2 .  I
simulate 20 sets, 30000 realizations/set.
I take the average of each set, and all 20 of the averages are 
much higher than 0.14.  In fact, they are All about 2-3 stdev higher. 
It seems that R is not generating the random variable with the mean I ask for.

The Details:

I simulate Prices [Pt] and Returns that are defined as 
Rt = [P(t)-P(t-1)]/P(t-1).  The returns have Normal distribution
=>  Rt is has N(mu=0.14, sigmaSquared=0.2) distribution.

To do the task above, I use the following files:

c:\ProgramFiles\R\rw1050\Code1.txt

consists of the following text:

====================
gpsim1_function(t,runs,mu,delta,stdev){
jold=0
xx_rep(0,t*runs)
for(i in 1:runs){
ind=(i-1)*t
xx[ind+1]=1
for(j in 2:t) {
xx[ind+j]_xx[ind+j-1]*
exp(rnorm(n=1,mean=mu*delta*delta,sd=stdev*delta))
}
}
xx

}

=======================

File "experiment4commands.txt" has text:

source("code1.txt")
aa=gpsim1(30000,1,0.14,0.052342392,0.4472135955)
sink("experiment4mu14z20Path1.txt")
for(i in 1:length(aa)){print(aa[i])}
sink()

source("code1.txt")
aa=gpsim1(30000,1,0.14,0.052342392,0.4472135955)
sink("experiment4mu14z20Path2.txt")
for(i in 1:length(aa)){print(aa[i])}
sink()

:
:
:

source("code1.txt")
aa=gpsim1(30000,1,0.14,0.052342392,0.4472135955)
sink("experiment4mu14z20Path20.txt")
for(i in 1:length(aa)){print(aa[i])}
sink()
=========================
{Above, 0.052342392 = 1/365, it is my discretization factor.}


I paste the text from the "experiment4commands.txt" file into R.

I now open each of the resulting 20 files in Excel.  The output [with 
variable Pt] is located in column B.  To get back Rt, I paste
the following formula into cells C2:C30000 - 
Ci = [B(i)-B(i-1)]/B(i-1).

My discretization factor was 1/365, thus, when I take the Average 
[and the Variance] of cells C2:C30000, I also have to multiply this 
number by 365, to get back the annual Rt.

I then made a table of the resulting 365*Average(C2:C30000), for each
of the 20 Sets that I have simulated using the methodology described above.


Set     365*Average(C2:C30000) for the set

1	0.284948
2	0.238827892
3	0.176492
4	0.242007
5	0.263882
6	0.211672
7	0.203536
8	0.19211
9	0.225185
10	0.234478
11	0.272792
12	0.205158
13	0.189349
14	0.222317
15	0.220652
16	0.165553
17	0.294211
18	0.333488
19	0.147971
20	0.290742

The average of the realizations of mu [true value = 0.14] over 
20 paths, each with 30000 observations, is about 0.23.

For this case, daily mu = mu over each algorithm step
is 0.14/365 = 0.0003835616438

The daily variance is 0.2/365 = 0.0005479452055

However, according to statistics, as the number of observations increases
variance will fall.  In particular, for 30,000 observations, the 
StDev is calculated as follows:  

SQRT(0.2/365)/SQRT(30000) = 0.0001351474 = {1 STDEV}

The average of the values in the table above, 0.23, is

[(0.23/365) - (0.14/365)]/0.0001351474 = 1.824 StDeviations above the mean.
As you can see Most of the 20 sets had realized values about 2-3 StDeviations
higher than the mean.

Something is wrong.

Thank you for considering this bug report.