Does simulation time step have an impact on results (Vasicek model for instance)

Question

I am trying to make some computations using Vasicek short rate model. Especially I am trying to compare exact expectation(obtained with the formula) and the expectation from Monte Carlo simulation.

Exact computation

I use:

$E[r_t] = r_0 * \exp(-a*t) + (\theta/a)*(1-\exp(-a*t))$

public override double GetExpectation(double r0, double t)
{
  double expectation = r0 * Math.Exp(-_a * t) + (_theta / _a) * (1 - Math.Exp(-_a * t));
  return expectation;
}

Monte Carlo simulations

I use the following method:

• I compute r from a time t to a time t+dt using:

  public override double ComputeNextValue(double r0, double dt)
  {
      RandomVariableGenerator rvg = RandomVariableGenerator.GetInstance();
  
      double randomGaussian = rvg.GetNextRandomGaussian();
      double r_t_dt = (_theta - _a*r0)*dt + _sigma * Math.Sqrt(dt) * randomGaussian;
  
      return r_t_dt;
  }
  

• then a compute a path from 0 to t with dt as time step using:

  public override double ComputeValue(double r0, double t, double dt)
  {
    double x = r0;
  
    for(double slot = dt; slot <= t; slot += dt)
    {
      x = ComputeNextValue(x, dt);
    }
    return x;
  }
  

• Then I compute the Monte Carlo Expectation using:

  public override double ComputeMonteCarloExpectation(double r0, double t, double dt, int nreps)
  {
     double sum = 0.0;
     double value;
     for (int i = 0; i < nreps; i++)
     {
       value = ComputeValue(r0, t, dt);
       sum += value;
     }
     return sum / nreps;
  }
  

  I use the following parameters:

  double sigma = 0.03;
  double r0 = 0.03;
  double theta = 0.1;
  double a = 0.3;
  
  int nreps = 1000;
  double t = 1;
  

For dt = 0.1:

Exact expectation: 0,108618473059879;  
Monte Carlo Expectation: 0,0101464832161612  

For dt = 1:

Exact expectation: 0,108618473059879;  
Monte Carlo Expectation: 0,092058844704742  

using dt = 1 leads to a result close to exact value while using dt = 0.1 seems to lead to a result having a 0.1 factor difference with exact one.

I think I am doing something wrong but I can't figure it out. Do you have an idea?

Answer 1

It would probably help if you sum up the Euler steps, either in

public override double ComputeNextValue(double r0, double dt)
{
  ...
  return r_0+r_t_dt;
}

or in the integration loop

public override double ComputeValue(double r0, double t, double dt)
{
  double x = r0;

  for(double slot = dt; slot <= t; slot += dt)
  {
    x += ComputeNextValue(x, dt);
  }
  return x;
}

---

The SDE $dr_t=(θ-ar_t)dt+\sigma dW_t$ gets approximated in Euler fashion as $$ \Delta r_t=(θ-ar_t)Δt+\sigma ΔW_t+O(Δt^{3/2}) $$ where $\Delta r_t=r_{t+Δt}-r_t$ so that to compute the next value you need to add the increment to the current value.

One can also directly give an explicit formula for the paths, as $$ r_t=e^{-at}(r_0+σ\widetilde W_{(e^{2at}-1)/(2a)})+\fracθa(1-e^{at}) $$ where $\widetilde W$ is a different Wiener process.