python: simulation de mouvement brownien géométrique

Question

Une simulation de base de GBM ne semble pas fonctionner. Qu'est-ce que je fais mal? Le code suivant renvoie toujours des valeurs inférieures à 1E-20, au lieu de quelque chose distribués de façon aléatoire autour de 1.0:

import math 
import random 

p = 1 
dt = 1 
mu = 0 
sigma = 1 
for k in range(100): 
    p *= math.exp((mu - sigma * sigma/2) * dt + 
     sigma * random.normalvariate(0, dt * dt)) 
print(p) 

Je suis en cours d'exécution:

ActivePython 3.1.2.3 (ActiveState Software Inc.) basée sur Python 3.1.2 (r312: 79147, 22 mars 2010, 12:30:45) [MSC v.1500 64 bits (AMD64)] sur win32

Mon système d'exploitation est Windows 7 Professionnel sur le processeur i7-930 (64- bit).

Je serais heureux d'exécuter d'autres tests sur ma machine pour isoler le problème.

Answer 1

J'ai trouvé la réponse. Il n'y a pas de problème avec le code. C'est juste que la distribution lognormale résultante a un paramètre d'échelle énorme = 1 * sqrt (100) = 10. Avec une échelle de 10, l'asymétrie est folle.

Ainsi, même si la moyenne de la distribution est 1.0, il me faudra des milliards d'itérations (sinon milliards de milliards) pour voir un seul nombre supérieur à 1.0.

Answer 2

Semble bien:

import math 
import random 

def compute(): 
    p = 1 
    dt = 1 
    mu = 0 
    sigma = 1 
    for k in range(100): 
     p *= math.exp((mu - sigma * sigma/2) * dt + 
         sigma * random.normalvariate(0, dt * dt)) 
    return p 

print [compute() for x in range(20)] 

donne:

[118.85952235362008, 7.3312246638859059e-14, 29.509674994408684, 1.8720575806397408, 1.5882398997219834e-05, 2967.524471541024, 0.0031130343677571093, 19942.669293314699, 0.00011878818261757498, 5382.80288111769, 0.22867624175360118, 0.028535167622775418, 12.6324011631628, 20.604456159054738, 0.0034504567371058613, 6.5804828930878056e-06, 6398.0493448486704, 0.0014978345496292245, 721546.38343724841, 1285.546939393781] 

Ce utilise python 2.6.1

Answer 3

L'exécution du code en Python 2.6.6 obtient des réponses raisonnables, alors qu'il en cours d'exécution Python 3.1.2 donne les petits nombres que vous décrivez. Je pense que c'est parce que dans Python 2.x l'opérateur de division donne un résultat entier en divisant deux entiers, alors qu'en Python 3.x il donne un résultat à virgule flottante dans la même situation. Par conséquent, la division par deux dans le calcul donne des résultats différents selon la version.

Pour la rendre cohérente, la force la sortie de la division à un entier:

p *= math.exp(int(mu - sigma * sigma/2)) * dt + 
    sigma * random.normalvariate(0, dt * dt)) 

Cela rend la sortie de la même dans les deux 2.x et 3.x. Un ensemble d'échantillons de sortie sur ma machine est:

0,0243898032782, 6126,78299771, ,00450839758621, 1,17316856812, ,00479489258202, 4.88995369021e-06, 0,033957530608, 29,9492464423, 3,16953460691

Cela semble être dans l'ordre de grandeur de ce que vous cherchez pour.

Answer 4

Utilisation division flottante (//) plutôt que la division entière (/) en Python 3.1 fonctionne:

import math 
import random 

p = 1 
dt = 1 
mu = 0 
sigma = 1 
for k in range(100): 
    p *= math.exp((mu - sigma * sigma // 2) * dt + 
     sigma * random.normalvariate(0, dt * dt)) 
print(p) 

Sur une course exemple, j'ai obtenu les numéros suivants:

0,0989269233704
2536660,91466
2146,09989782
0,502233504924
0,43052439984
14,1156450335