在Python中使用逆變換方法生成隨機變量
目標
在仿真理論中,生成隨機變量是最重要的“構建塊”之一,而這些隨機變量大多是由均勻分布的隨機變量生成的。其中一種可以用來產生隨機變量的方法是逆變換法。在本文中,我將向您展示如何使用Python中的逆變換方法生成隨機變量(包括離散和連續的情況)。
概念
給定隨機變量U,其中U在(0,1)中均勻分布。 假設我們要生成隨機變量X,其中累積分布函數(CDF)為:

逆變換方法的思想是通過如下使用其逆CDF從任何概率分布中生成一個隨機數。

對于離散隨機變量,步驟略有不同。假設我們想生成一個離散隨機變量X的值,它具有一個概率質量函數(PMF)

為了生成X的值,需要生成一個隨機變量U,U在(0,1)中均勻分布,并且定義

通過以上步驟,我們可以按如下方法創建逆變換方法的算法。

連續隨機數代碼實現
首先,我們實現此方法以生成連續隨機變量。 假設我們要模擬一個隨機變量X,該變量遵循均值λ(即X〜EXP(λ))的指數分布。 我們知道指數分布的概率分布函數(PDF)是

CDF如下

然后,我們可以使用以下的方法寫出逆CDF

在Python中,我們可以通過如下編寫這些代碼行來簡單地實現它。
- ### Generate exponential distributed random variables given the mean
- ### and number of random variables
- def exponential_inverse_trans(n=1,mean=1):
- U=uniform.rvs(size=n)
- X=-mean*np.log(1-U)
- actual=expon.rvs(size=n,scale=mean)
- plt.figure(figsize=(12,9))
- plt.hist(X, bins=50, alpha=0.5, label="Generated r.v.")
- plt.hist(actual, bins=50, alpha=0.5, label="Actual r.v.")
- plt.title("Generated vs Actual %i Exponential Random Variables" %n)
- plt.legend()
- plt.show()
- return X
我們可以通過運行以下示例來嘗試上面的代碼。 請注意,由于我們要生成隨機變量,因此結果可能會有所不同。
- cont_example1=exponential_inverse_trans(n=100,mean=4)
- cont_example2=exponential_inverse_trans(n=500,mean=4)
- cont_example3=exponential_inverse_trans(n=1000,mean=4)


看起來很有趣。 如果將其與實際變量進行比較,我們可以看到生成的隨機變量具有非常相似的結果。 可以調整均值(請注意,我為expon.rvs()函數定義的均值是指數分布中的比例參數)和/或 生成的隨機變量的數量,以查看不同的結果。
離散隨機數實現代碼
對于離散隨機變量情況,假設我們要模擬遵循以下分布的離散隨機變量情況X

首先,我們編寫函數以使用這些代碼行為一個樣本生成離散隨機變量。
- ### Generate arbitary discrete distributed random variables given
- ### the probability vector
- def discrete_inverse_trans(prob_vec):
- U=uniform.rvs(size=1)
- if U<=prob_vec[0]:
- return 1
- else:
- for i in range(1,len(prob_vec)+1):
- if sum(prob_vec[0:i])<U and sum(prob_vec[0:i+1])>U:
- return i+1
然后,我們創建一個函數以使用這些代碼行生成許多隨機變量樣本。
- def discrete_samples(prob_vec,n=1):
- sample=[]
- for i in range(0,n):
- sample.append(discrete_inverse_trans(prob_vec))
- return np.array(sample)
最后,我們創建一個函數來模擬結果,并通過這些代碼行將其與實際結果進行比較。
- def discrete_simulate(prob_vec,numbers,n=1):
- sample_disc=discrete_samples(prob_vec,n)
- unique, counts=np.unique(sample_disc,return_counts=True)
- fig=plt.figure()
- ax=fig.add_axes([0,0,1,1])
- prob=counts/n
- ax.bar(numbers,prob)
- ax.set_title("Simulation of Generating %i Discrete Random Variables" %n)
- plt.show()
- data={'X':unique,'Number of samples':counts,'Empirical Probability':prob,'Actual Probability':prob_vec}
- df=pd.DataFrame(data=data)
- return df
我們可以在下面運行一些示例以查看結果。 同樣,請注意,由于我們要生成隨機變量,因此結果可能會有所不同。
- prob_vec=np.array([0.1,0.3,0.5,0.05,0.05])
- numbers=np.array([1,2,3,4,5])
- dis_example1=discrete_simulate(prob_vec, numbers, n=100)
- dis_example2=discrete_simulate(prob_vec, numbers, n=500)
- dis_example3=discrete_simulate(prob_vec, numbers, n=1000)


- In[11]: dis_example1
- Out[11]:
- X Number of samples Empirical Probability Actual Probability
- 0 1 8 0.08 0.10
- 1 2 35 0.35 0.30
- 2 3 50 0.50 0.50
- 3 4 5 0.05 0.05
- 4 5 2 0.02 0.05In[12]: dis_example2
- Out[12]:
- X Number of samples Empirical Probability Actual Probability
- 0 1 53 0.106 0.10
- 1 2 159 0.318 0.30
- 2 3 234 0.468 0.50
- 3 4 30 0.060 0.05
- 4 5 24 0.048 0.05In[13]: dis_example3
- Out[13]:
- X Number of samples Empirical Probability Actual Probability
- 0 1 108 0.108 0.10
- 1 2 290 0.290 0.30
- 2 3 491 0.491 0.50
- 3 4 51 0.051 0.05
- 4 5 60 0.060 0.05
結果很有趣! 我們可以看到,隨著我們增加隨機變量樣本的數量,經驗概率越來越接近實際概率。 嘗試使用不同數量的樣本和/或不同的分布進行實驗,以查看不同的結果。
總結
這種逆變換方法是統計中非常重要的工具,尤其是在仿真理論中,在給定隨機變量均勻分布在(0,1)中的情況下,我們想生成隨機變量。 研究案例本身非常廣泛,您可以使用在生成經驗累積分布函數,預測分析中使用到的這種方法。