wisesummer

Thank you for your codes and efforts!
I think "r1 = gsl_ran_gaussian (r, mu)" and "r2 = r1+ gsl_ran_gaussian (r, mu)" should be replaced with "r1 = sqrt(1/1000)*gsl_ran_gaussian (r, mu)" and "r2 = r1+ sqrt(1/1000)*gsl_ran_gaussian (r, mu)", respectively, since you divided each process into 1000 segments.
And if the there is no points falling in interval (-0.001,0.001) during the whole process, you just assign 0 to "time" right?
So here is the problem, what should we assign to "time" when there is no points falling in interval (-0.001,0.001), 0 or 1?
My opinion, intuitively,  is 0 or 1 with prob. 0.5. provided "s" is in [0,1].

[ 本帖最后由 wisesummer 于 2007-12-24 01:48 AM 编辑 ]

luqi

And if the there is no points falling in interval (-0.001,0.001) during the whole process, you just assign 0 to "time" right?

这个问题不存在我代码中有一个if（time！＝0） 的判断，如果整个过程都没有落在归入0的区间，那么这个过程是不计入统计的。

你只要看看我的图就知道了。我模拟了5000次路径，0那一点对应的个数也就是十几，如果每次都记0，那么0就是5000了。（我记录的是真实的数据，没有经过normalization。）

luqi

I think "r1 = gsl_ran_gaussian (r, mu)" and "r2 = r1+ gsl_ran_gaussian (r, mu)" should be replaced with "r1 = sqrt(1/1000)*gsl_ran_gaussian (r, mu)" and "r2 = r1+ sqrt(1/1000)*gsl_ran_gaussian (r, mu)", respectively, since you divided each process into 1000 segments.

我不是先产生整个过程，再把整个过程分成1000份，我的整个过程是一步一步产生的
先产生人r1
下一步r2 ＝ r1 ＋ 随机数
判断r2是不是接近于0
是就把时间记录下来

然后把r2传递给r1
如此循环，如果记录下每次r1的直，就形成一条路径

wisesummer

I think the most badly defect of Luqi's simulation is setting a "zero zone" to select the zero points.

Although only discrete sample paths could be generated, we still can find the "exact" zero point. The trick is we should simulate each whole process prior to checking its zero points, rather than checking as simulating.
Once the whole process is generated, all the zero points are midpoints between any two most adjacent Brownian Motions with negative product of values, i.e.
t is "zero point" if w(t+dt)*w(t-dt)<0.
And if the whole process fails to have any zero points, just assign 1 to t, by definition L=sup(0<s<1, w(s)=0).
Following is my simulation results:

the red light is nonparametric density estimation, since this is a bounded support density, nonparametric estimation without adjustment would be bothered by the boundary effect, this is why the red line is systemetically higher than blue line within this interval. But the shapes of these two lines match well. At this point, we could not say the theoretical density is wrong.
Another evidence is following Histogram

from this draft, the shape of theoretical density is more apparent.

wisesummer

Luqi, please don't mind I mention following two minor defects of your simulation, I mean nothing more than discussion.

1. "这个问题不存在我代码中有一个if（time！＝0） 的判断，如果整个过程都没有落在归入0的区间，那么这个过程是不计入统计的。"
You should not ignore any simulated path, i.e. never only pick what you want in simulation. In that way, you could only reach a truncated density rather than the true simulation density.

2. "我不是先产生整个过程，再把整个过程分成1000份，我的整个过程是一步一步产生的"
You don't get me at this point. "gsl_ran_gaussian (const gsl_rng * r, double sigma)" returns a Gaussian random variate, with mean zero and standard deviation sigma, so in your setup each step you generate a standard normal, finally you get a process of W(1000) rather than W(1). But I bet you know dW~ N(0,sqrt(dt)).

luqi

受教了，你说的对，应该检查临近的两点，一点大于0，一点小与0，那么中直就是0。已经修改了程序，这下对了。有讨论才有进步，谢谢！！！！

[ 本帖最后由 luqi 于 2007-12-24 03:32 PM 编辑 ]