Basic Gaussian Process simulation/ white noise - massive problem

Question

I have a certain problem with Mathematica or more precisely, myself handling it. It is correlated with a mathematical problem and its realization/simulation in Mathematica.

I want to start with the mathematical part: I am given a gaussian distributed random variable $\epsilon(t)$ with $E\left[\epsilon(t)\right]=0$ and correlation function $E\left[\epsilon(t)\epsilon(s)\right]=\sigma_{\epsilon}^2\delta(t-s)$, so gaussian white noise. The thing i want to calculate is $E\left[\exp\left(-\frac{1}{2}\int_{0}^{t}\epsilon(t')\mathrm{d}t'\right)\right]$. We see that the integral of a gaussian-distributed random variable is again gaussian distributed, so from the moment generating function of the gaussian distribution we get: \begin{equation*} E\left[\exp\left(k X\right)\right]=G(k)=\exp\left(\mu_X k+\frac{1}{2}k^2\sigma_X^2 \right) \end{equation*}

With $X=\int_{0}^{t}\epsilon(t')\mathrm{d}t'$ it follows that $\mu_X=E\left[X\right]=0$ and $E\left[X^2\right]=\sigma_{\epsilon}^2t$.

So it follows: \begin{equation*} G\left(-\frac{1}{2}\right)=E\left[\exp\left(-\frac{1}{2}\int_{0}^{t}\epsilon(t')\mathrm{d}t'\right)\right]=\exp\left(\frac{1}{8}\sigma^2_\epsilon t\right) \end{equation*}

So far so good, this is the thing I want to 'get as a result' theoretically. Now I wanted to verify it with a little simulation in Mathematica.

Here is what I did in Mathematica (I tried to add comments in my code, please ask if there are questions/ambiguities)

I really hope you guys can help me, I can't work any further without having solved that certain problem...and I'm getting desperate. I also seems that the result depends on the time-steps I choose... But the plots really don't match.

(simulated with )

(theoretically)

Standarddeviation=5;
(*define standard deviation of the gaussian random variable epsilon *)

DeltaT=0.04; 
(* time steps with which i want to evaluate the integral *)

Repetitions=5000; 
(* in order to get a good average, i take 5000 reps*)

Time = 20; 
(*maximal time *)



AverageExactly[t_] := Exp[1/8*Standarddeviation^2*t];
 (*exact average theoretically*)

(*It follows the loop for the 5000 repetitions *)
For[z = 0, z < Repetitions, z++, 
  GaussianDistribution=RandomVariate[NormalDistribution[0,Standarddeviation], 
   Round[Time/DeltaT + 1]];
(*generate Gaussian numbers, there shall be enough that all 
time is covered*)
 SimulatedValue[t_] := Exp[-1/2*Deltat*Sum[GaussianDistribution[[i]],{i,1,t/Deltat}]]; 
 (* Evaluate the integral in form of a sum*)
 ListSimulatedValues = 
  Table[SimulatedValue[t], {t, 0, Time, DeltaT}]; 
(*make a table out of the simulated values*)
 ListAxis = Table[t, {t, 0, Time, Deltat}]; 
(*define the axis-values*)


 ListSimulatedValuesPlot[z] = 
  Transpose@{# & @@@ ListAxis, # & @@@ 
     ListSimulatedValues}; (*merge both tables*)

 ]

SimulatedAveraged = List[1/Repetitions*Sum[ListSimulatedValuesPlot[z],{z, 0,  Repetitions - 1}]]
 (*average over all 5000 repetitions *)
  ListPlot[SimulatedAveraged =     List[1/Repetitions*Sum[ListSimulatedValuesPlot[z], {z, 0, Repetitions - 1}]], 
 PlotRange -> {{0, 20}, {0.9, 2}}] (*Plot simulated values *)
Plot[AverageExactly[t], {t, 0, 20}](*Plot theoretical values values *)

Answer 1

See this for example. The simulation is band-limited, requiring a correction to the variance of the distribution. If you do that and stick to low values for time the agreement is okay.

sd = 5;
dt = 0.001;
reps = 10000;
time = 1;

exact[t_] := Exp[1/8*sd^2*t];

gd = RandomVariate[NormalDistribution[0, sd/Sqrt[dt]], {Round[time/dt + 1], reps}];

sim = Mean[Exp[-(dt/2) Transpose@Accumulate[gd]]];

Show[
 Plot[exact[t], {t, 0, time}, PlotLegends -> {"Exact"}],
 ListLinePlot[sim, DataRange -> {0, time}, PlotStyle -> Orange, 
  PlotLegends -> {"Simulation"}]]