Значительная разница в ценах при вычислении цены опциона на использование метода Монте-Карло в C++ и Python?

Question

Я пытаюсь рассчитать цену еврозоны с использованием подхода Монте-Карло. Я закодировал алгоритм в C++ и Python. Насколько я знаю, реализация правильная, и по мере увеличения количества N (количество проб) цена должна сходиться к аналогичной ценности в обеих программах.

Моя проблема заключается в том, что по мере роста N, скажем, от 1000 до 10000 испытаний цены сходятся к двум различным значениям. В C++ цена сходится к значению 3,30, а с Python она сходится к 3,70.

Я думаю, что разрыв 0,40 слишком широк, я должен получить более похожие результаты. Почему этот разрыв настолько велик? Что я сделал не так? Кажется, я не могу найти свою ошибку.

Вот код, который я использовал:

Python

import numpy as np 
import matplotlib.pyplot as plt 


def stoc_walk(p,dr,vol,periods): 
    w = np.random.normal(0,1,size=periods) 
    for i in range(periods): 
     p += dr*p + w[i]*vol*p 
    return p 

s0 = 10; 
drift = 0.001502 
volatility = 0.026 
r = 0.02 
days = 255 
N = 10000 
zero_trials = 0 

k=12 
payoffs = [] 

for i in range(N): 
    temp = stoc_walk(s0,drift,volatility,days) 
    if temp > k: 
     payoff = temp-k 
     payoffs.append(payoff*np.exp(-r)) 
    else: 
     payoffs.append(0) 
     zero_trials += 1 

payoffs = np.array(payoffs) 
avg = payoffs.mean() 

print("MONTE CARLO PLAIN VANILLA CALL OPTION PRICING") 
print("Option price: ",avg) 
print("Initial price: ",s0) 
print("Strike price: ",k) 
print("Daily expected drift: ",drift) 
print("Daily expected volatility: ",volatility) 
print("Total trials: ",N) 
print("Zero trials: ",zero_trials) 
print("Percentage of total trials: ",zero_trials/N) 

C++

//Call option Monte Carlo evaluation; 

#include <iostream> 
#include <random> 
#include <math.h> 
#include <chrono> 

using namespace std; 

/* double stoc_walk: returns simulated price after periods 

    p = price at t=t0 
    dr = drift 
    vol = volatility 
    periods (days) 
*/ 
double stoc_walk(double p,double dr,double vol,int periods) 
{ 
    double mean = 0.0; 
    double stdv = 1.0; 

    /* initialize random seed: */ 
    int seed = rand() %1000 + 1; 
    //unsigned seed = std::chrono::system_clock::now().time_since_epoch().count(); 
    std::default_random_engine generator(seed); 
    std::normal_distribution<double> distribution(mean,stdv); 

    for(int i=0; i < periods; i++) 
    { 
     double w = distribution(generator); 
     p += dr*p + w*vol*p; 
    } 
    return p; 
} 

int main() 
{ 
    //Initialize variables 
    double s0 = 10;    //Initial price 
    double drift = 0.001502; //daily drift 
    double volatility = 0.026; //volatility (daily) 
    double r = 0.02;   //Risk free yearly rate 
    int days = 255;    //Days 
    int N = 10000;    //Number of Monte Carlo trials 
    double zero_trials = 0; 

    double k = 12;    //Strike price 
    int temp = 0;    //Temporary variable 
    double payoffs[N];   //Payoff vector 
    double payoff = 0; 

    srand (time(NULL));   //Initialize random number generator 

    //Calculate N payoffs 
    for(int j=0; j < N; j++) 
    { 
     temp = stoc_walk(s0,drift,volatility,days); 
     if(temp > k) 
     { 
      payoff = temp - k; 
      payoffs[j] = payoff * exp(-r); 
     } 
     else 
     { 
      payoffs[j] = 0; 
      zero_trials += 1; 
     } 
    } 

    //Average the results 
    double sum_ = 0; 
    double avg_ = 0; 
    for(int i=0; i<N; i++) 
    { 
     sum_ += payoffs[i]; 
    } 
    avg_ = sum_/N; 

    //Print results 
    cout << "MONTE CARLO PLAIN VANILLA CALL OPTION PRICING" << endl; 
    cout << "Option price: " << avg_ << endl; 
    cout << "Initial price: " << s0 << endl; 
    cout << "Strike price: " << k << endl; 
    cout << "Daily expected drift: " << drift*100 << "%" << endl; 
    cout << "Daily volatility: " << volatility*100 << "%" << endl; 
    cout << "Total trials: " << N << endl; 
    cout << "Zero trials: " << zero_trials << endl; 
    cout << "Percentage of total trials: " << zero_trials/N*100 << "%"; 

    return 0; 
} 

Answer 1

В вашей реализации на C++ есть ошибка.

int seed = rand() %1000 + 1; 

Семена будут находиться в диапазоне [1 ... 1000], что означает ваша гауссовая выборка сильно коррелирует

И вторая ошибка в том, что temp переменных объявлена ​​Int

OK, после того, как простая модификация, см ниже код, avg в C++ версии 3,67

C++: http://codepad.org/D5UZql2P

Python: http://codepad.org/TeAYSwkV

Answer 2

Потому что это Монте-Карло, и потому, что реализация алгоритмов распределения находится вне вашего собственного кода , есть все возможности, что ваша проблема связана с генерацией случайных чисел, а не с вашим кодом, как показано.

Это означает, что вы вряд ли получите код справки об этом здесь.

Что вы должны сделать, так это убедиться, что две реализации работают одинаково для одних и тех же данных. Для этого вы должны сгенерировать набор случайных входных переменных и сохранить их отдельно - затем передать тот же случайный поток в оба алгоритма. Это единственный способ контролировать часть генерации случайных чисел в Монте-Карло и доказать, что реализации действительно логически идентичны.