numpy.random.normal distribusi berbeda: memilih nilai dari distribusi

Saya memiliki distribusi energi hukum pangkat dan saya ingin memilih n energi acak berdasarkan distribusinya. Saya mencoba melakukan ini secara manual menggunakan nomor acak tetapi itu terlalu tidak efisien untuk apa yang ingin saya lakukan. Saya bertanya-tanya apakah ada metode di numpy (atau lainnya) yang berfungsi seperti numpy.random.normal, kecuali alih-alih menggunakan distribusi normal, distribusinya dapat ditentukan. Jadi menurut saya contohnya mungkin terlihat seperti (mirip dengan numpy.random.normal):

import numpy as np

# Energies from within which I want values drawn
eMin = 50.
eMax = 2500.

# Amount of energies to be drawn
n = 10000

photons = []

for i in range(n):

    # Method that I just made up which would work like random.normal,
    # i.e. return an energy on the distribution based on its probability,
    # but take a distribution other than a normal distribution
    photons.append(np.random.distro(eMin, eMax, lambda e: e**(-1.)))

print(photons)

Pencetakan photons akan memberi saya daftar sepanjang 10.000 yang diisi oleh energi dalam distribusi ini. Jika saya membuat histogram ini, nilai bin akan jauh lebih besar pada energi yang lebih rendah.

Saya tidak yakin apakah metode seperti itu ada tetapi sepertinya memang demikian adanya. Saya harap jelas apa yang ingin saya lakukan.

Sunting:

Saya telah melihat numpy.random.power tetapi eksponen saya -1 jadi menurut saya ini tidak akan berhasil.


python numpy random distribution
person davly    schedule 07.07.2014    source sumber
comment
pdf apa sebenarnya yang kamu inginkan? distribusi daya adalah kasus khusus beta, dapatkah Anda menggunakannya docs.scipy.org/doc/numpy/reference/generated/ ?   -  person wim    schedule 07.07.2014
comment
@wim , Saya yakin saya menginginkan fungsi sepotong-sepotong yang f(x)=0 di luar rentang energi saya dan f(x)=x**a (di mana a dapat berupa nilai dari -5 hingga 5) di dalamnya. Saya tidak bisa melihat cara kerja beta di sini.   -  person davly    schedule 08.07.2014
comment
@davly memperbarui jawaban saya dengan cuplikan kode jika ini bermanfaat   -  person John Greenall    schedule 08.07.2014

Jawaban (3)


arrow_upward
1
arrow_downward

Mengambil sampel dari PDF sembarangan sebenarnya cukup sulit. Ada buku besar dan padat tentang cara mengambil sampel secara efisien dan akurat dari kelompok distribusi standar.

Sepertinya Anda mungkin bisa menggunakan metode inversi khusus untuk contoh yang Anda berikan.

person Robert Kern    schedule 07.07.2014
comment
Bagaimana cara menerapkan metode inversi khusus? Saya tidak melihatnya di sini - person davly; 08.07.2014
comment
Turunkan fungsi kebalikan dari CDF. Gunakan random_sample() untuk mendapatkan nilai yang terdistribusi secara merata antara 0 dan 1. Lewatkan nilai tersebut melalui CDF terbalik untuk mendapatkan nilai yang mengikuti distribusi yang Anda inginkan. Dalam kasus Anda, CDF terbalik adalah lambda u: eMin*(eMax/eMin)**u. - person Robert Kern; 09.07.2014

arrow_upward
1
arrow_downward

Jika Anda ingin mengambil sampel dari distribusi arbitrer, Anda memerlukan kebalikan dari fungsi kepadatan kumulatif (bukan pdf).

Anda kemudian mengambil sampel probabilitas secara seragam dari rentang [0,1] dan memasukkannya ke dalam kebalikan dari cdf untuk mendapatkan nilai yang sesuai.

Seringkali tidak mungkin memperoleh cdf dari pdf secara analitis. Namun, jika Anda ingin memperkirakan distribusinya, Anda dapat melakukannya dengan menghitung f(x) secara berkala pada domainnya, kemudian melakukan cumsum pada vektor ini untuk mendapatkan perkiraan cdf dan dari sini perkiraan kebalikannya.

Cuplikan kode kasar:

import matplotlib.pyplot as plt
import numpy as np
import scipy.interpolate

def f(x):
   """
   substitute this function with your arbitrary distribution
   must be positive over domain
   """
   return 1/float(x)


#you should vary inputVals to cover the domain of f (for better accurracy you can
#be clever about spacing of values as well). Here i space them logarithmically
#up to 1 then at regular intervals but you could definitely do better
inputVals = np.hstack([1.**np.arange(-1000000,0,100),range(1,10000)])

#everything else should just work
funcVals = np.array([f(x) for x in inputVals])
cdf = np.zeros(len(funcVals))
diff = np.diff(funcVals)
for i in xrange(1,len(funcVals)):
   cdf[i] = cdf[i-1]+funcVals[i-1]*diff[i-1]
cdf /= cdf[-1]

#you could also improve the approximation by choosing appropriate interpolator
inverseCdf = scipy.interpolate.interp1d(cdf,inputVals)

#grab 10k samples from distribution
samples = [inverseCdf(x) for x in np.random.uniform(0,1,size = 100000)]

plt.hist(samples,bins=500)
plt.show()
person John Greenall    schedule 07.07.2014

arrow_upward
-1
arrow_downward

Mengapa Anda tidak menggunakan eval dan memasukkan distribusinya ke dalam sebuah string?

>>> cmd = "numpy.random.normal(500)"
>>> eval(cmd)

Anda dapat memanipulasi string sesuai keinginan untuk mengatur distribusinya.

person BigBrownBear00    schedule 07.07.2014
comment
maaf, saya salah memahami pertanyaan Anda. - person BigBrownBear00; 07.07.2014