Take the 2-minute tour ×
Stack Overflow is a question and answer site for professional and enthusiast programmers. It's 100% free, no registration required.

I generate a random sample from the two dimensional distribution, But I got this run time error. Could someone tell me how to fix this error? This is the error message:


ValueError Traceback (most recent call last) in () 25 # Do the accept/reject comparison 26 if x < P(x,y): ---> 27 samples[accepted] = x,y 28 accepted += 1 29

ValueError: setting an array element with a sequence.

P = lambda x, y: np.exp(-x**2-y**2-x**2*y**2)

# domain limits
xmin = -2 # the lower limit of our domain
xmax = 2 # the upper limit of our domain

# range limit (supremum) for y
ymin = -2
ymax = 2

N = 10000 # the total of samples we wish to generate
accepted = 0 # the number of accepted samples
samples = np.zeros(N)
count = 0 # the total count of proposals

# generation loop
while (accepted < N):

    # pick a uniform number on [xmin, xmax) (e.g. 0...10)
    x = np.random.uniform(xmin, xmax)

    # pick a uniform number on [0, ymax)
    y = np.random.uniform(ymin,ymax)

    # Do the accept/reject comparison
    if x < P(x,y):
        samples[accepted] = x,y
        accepted += 1

    count +=1

print count, accepted

# get the histogram info
hinfo = np.histogram(samples,30)

# plot the histogram
plt.hist(samples,bins=30, label=u'Samples', alpha=0.5);

# plot our (normalized) function
xvals=np.linspace(xmin, xmax, 1000)
plt.plot(xvals, hinfo[0][0]*P(xvals), 'r', label=u'P(x)')

# turn on the legend
plt.legend()
share|improve this question
    
Please include the complete error message in your question. –  BrenBarn Feb 14 at 2:37

1 Answer 1

The line here:

samples[accepted] = x,y

is the problem. samples is of type np.ndarray with a specific dtype. Probably np.float. numpy is expecting a single number but getting a tuple(x,y). Probably you meant to create a 2D array:

samples = np.zeros((N,2))

Then your code should work, though I haven't looked deeper into what it all is actually trying to do.

An Aside

I see you are looping over N while generating numbers. This is typical of someone just picking up numpy/scipy from another programming language. I would like to urge you to "think in arrays". This is far more elegant and much much faster:

x = np.random.uniform(xmin,xmax,(N,2))
y = np.random.uniform(ymin,ymax,(N,2))
accepted = x < P(x,y)
samples[accepted] = x[accepted]

The above code replaces your entire while loop.

share|improve this answer
    
Thank you, my code works. I m new to python :) –  Dido Du Feb 14 at 3:14

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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