setting seed for Monte Carlo

thethinker
107 ●3 ●6 ●13

Does sage seriously not have the ability to change the seed for a Monte Carlo execution? I assumed there is at some way to pass the time, or some other seed generator, to the algorithm.

f(x)=ln(x)
for i in range(10):
    I=monte_carlo_integral(f,[0],[1],100000, algorithm='plain')
    print(I[0])

-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281
-1.0002232192476281

EDIT: Based on a helpful comment: randstate is apparently not impacting the Monte Carlo routines.

I=monte_carlo_integral(f,[0],[1],100000, algorithm='plain')
print(I[0])
from sage.misc.randstate import randstate
r = randstate(54321)
print(r.seed())
I=monte_carlo_integral(f,[0],[1],100000, algorithm='plain')
print(I[0])
r = randstate()
print(r.seed())
I=monte_carlo_integral(f,[0],[1],100000, algorithm='plain')
print(I[0])

-1.0002232192476281
54321
-1.0002232192476281
305159320933427511743223216889072415483
-1.0002232192476281

Any other thoughts?

Comments

The random component is initialized inside the code of the function, the relevant lines are:

    try:
        # Initialize the random number generator
        gsl_rng_env_setup()
        type_rng = gsl_rng_default
        _rng = gsl_rng_alloc(type_rng)

So in order to have access to the randomizer, we need to go deep into the cython tunnel.

Why is it important to have a random component?

dan_fulea ( 0 years ago )

It's one of the basic methods to determine the statistical errors around a MC calculation. Run N times (with N different seeds) to find an average and standard deviation. You can check for Guassianity also, for example. I can't argue for how complicated implementation would be, but assuming Gaussian statistics limits the relevance of the implementation. Mathematica allows seed changes right at the cli, if that matters. Shall I file a feature request?

thethinker ( 0 years ago )

add a comment

Docstring: WARNING: the enclosing module is marked 'needs sage.groups sage.libs.gap sage.libs.ntl sage.libs.pari', so doctests may not pass. The "randstate" class. This class keeps track of random number states and seeds. Type "sage.misc.randstate?" for much more information on random numbers in Sage. Init docstring: WARNING: the enclosing module is marked 'needs sage.groups sage.libs.gap sage.libs.ntl sage.libs.pari', so doctests may not pass. Initialize a new "randstate" object with the given seed (which must be coercible to a Python long). If no seed is given, then a seed is automatically selected using "os.urandom()" if it is available, or the current time otherwise. EXAMPLES: sage: from sage.misc.randstate import randstate sage: r = randstate(54321); r

Comments

Good suggestion, but apparently sage is not accessing the randstate routines when doing monte carlo - see my edit.

thethinker ( 0 years ago )

add a comment

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

setting seed for Monte Carlo

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

setting seed for Monte Carlo savecancel

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

setting seed for Monte Carlo