1
$\begingroup$

I am trying to perform Monte Carlo Simulations using quasi random standard normal numbers. I understand that we can use sobol sequences to generate uniform numbers, and then use probability integral transform to convert them to standard normal numbers. My code gives unrealistic values of the simulated asset path:

import sobol_seq
import numpy as np
from scipy.stats import norm

def i4_sobol_generate_std_normal(dim_num, n, skip=1):
    """
    Generates multivariate standard normal quasi-random variables.
    Parameters:
      Input, integer dim_num, the spatial dimension.
      Input, integer n, the number of points to generate.
      Input, integer SKIP, the number of initial points to skip.
      Output, real np array of shape (n, dim_num).
    """

    sobols = sobol_seq.i4_sobol_generate(dim_num, n, skip)

    normals = norm.ppf(sobols)

    return normals

def GBM(Ttm, TradingDaysInAYear, NoOfPaths, UnderlyingPrice, RiskFreeRate, Volatility):
    dt = float(Ttm) / TradingDaysInAYear
    paths = np.zeros((TradingDaysInAYear + 1, NoOfPaths), np.float64)
    paths[0] = UnderlyingPrice
    for t in range(1, TradingDaysInAYear + 1):
        rand = i4_sobol_generate_std_normal(1, NoOfPaths)
        lRand = []
        for i in range(len(rand)):
            a = rand[i][0]
            lRand.append(a)
        rand = np.array(lRand)

        paths[t] = paths[t - 1] * np.exp((RiskFreeRate - 0.5 * Volatility ** 2) * dt + Volatility * np.sqrt(dt) * rand)
    return paths

GBM(1, 252, 8, 100., 0.05, 0.5)

array([[1.00000000e+02, 1.00000000e+02, 1.00000000e+02, ...,
        1.00000000e+02, 1.00000000e+02, 1.00000000e+02],
       [9.99702425e+01, 1.02116774e+02, 9.78688323e+01, ...,
        1.00978615e+02, 9.64128959e+01, 9.72154915e+01],
       [9.99404939e+01, 1.04278354e+02, 9.57830834e+01, ...,
        1.01966807e+02, 9.29544649e+01, 9.45085180e+01],
       ...,
       [9.28295879e+01, 1.88049044e+04, 4.58249200e-01, ...,
        1.14117599e+03, 1.08089096e-02, 8.58754653e-02],
       [9.28019642e+01, 1.92029616e+04, 4.48483141e-01, ...,
        1.15234371e+03, 1.04211828e-02, 8.34842557e-02],
       [9.27743486e+01, 1.96094448e+04, 4.38925214e-01, ...,
        1.16362072e+03, 1.00473641e-02, 8.11596295e-02]])

Values like 8.11596295e-02 should not be generated, hence I think there is something wrong in the code.

References: https://stats.stackexchange.com/questions/27450/best-method-for-transforming-low-discrepancy-sequence-into-normal-distribution, https://stackoverflow.com/questions/9412339/recommendations-for-low-discrepancy-e-g-sobol-quasi-random-sequences-in-pytho, https://github.com/naught101/sobol_seq

Any help is appreciated.

$\endgroup$
  • $\begingroup$ What is a "simple" implementation? What DID you find and why is it not OK? Implementation in pseudo-code or what language? $\endgroup$ – g g Apr 20 at 18:32
  • $\begingroup$ @gg: I edited the post to include the workings $\endgroup$ – SPaul Apr 21 at 5:07
  • $\begingroup$ So maybe "sobol_seq.i4_sobol_generate" is incorrect? What did you do to test correctness? And how can we see from your code that "8.11596295e-02 should not be generated"? $\endgroup$ – g g Apr 22 at 21:29
0
$\begingroup$

Is this simple enough in Python?

from scipy.stats import norm

uniform_rvs = generated_array_as_you_like(method='SOBOL_SEQUENCE')
normal_rvs = norm.ppf(uniform_rvs)
$\endgroup$
  • $\begingroup$ I edited the post to include the workings. Also, I didn't get what you meant by this : generated_array_as_you_like(method='SOBOL_SEQUENCE'). $\endgroup$ – SPaul Apr 21 at 5:09
  • 1
    $\begingroup$ Your workings are exactly what I have specified above. You have replaced generate_array_as_you_like(method='SOBOL_SEQUENCE') with isobol_seq.i4_sobol_generate(*args), i.e. you have generated a random uniform array with the method you like, and then generated normal random variables. If there is a problem it is not with your random number generator. The Normal distribution can contain any number from -inf to +inf $\endgroup$ – Attack68 Apr 21 at 5:48
  • $\begingroup$ If I use standard normal draws from the numpy library rand = np.random.standard_normal(NoOfPaths) then the price matches with the Black Scholes price. Hence I think the problem is with the random number generator. The value 8.11596295e-02 refers to a price in a path, and its very unlikely that the price would come down from 100 (initial price) to 8.11596295e-02. $\endgroup$ – SPaul Apr 21 at 7:19

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.