python 計(jì)算概率密度、累計(jì)分布、逆函數(shù)的例子
計(jì)算概率分布的相關(guān)參數(shù)時(shí)，一般使用 scipy 包，常用的函數(shù)包括以下幾個(gè)：
pdf：連續(xù)隨機(jī)分布的概率密度函數(shù)
pmf：離散隨機(jī)分布的概率密度函數(shù)
cdf：累計(jì)分布函數(shù)
百分位函數(shù)（累計(jì)分布函數(shù)的逆函數(shù)）
生存函數(shù)的逆函數(shù)（1 - cdf 的逆函數(shù)）
函數(shù)里面不僅能跟一個(gè)數(shù)據(jù)，還能跟一個(gè)數(shù)組。下面用正態(tài)分布舉例說明：
>>> import scipy.stats as st >>> st.norm.cdf(0) # 標(biāo)準(zhǔn)正態(tài)分布在 0 處的累計(jì)分布概率值 0.5 >>> st.norm.cdf([-1, 0, 1])# 標(biāo)準(zhǔn)正態(tài)分布分別在 -1， 0， 1 處的累計(jì)分布概率值 array([0.15865525, 0.5, 0.84134475]) >>> st.norm.pdf(0) # 標(biāo)準(zhǔn)正態(tài)分布在 0 處的概率密度值 0.3989422804014327 >>> st.norm.ppf(0.975)# 標(biāo)準(zhǔn)正態(tài)分布在 0.975 處的逆函數(shù)值 1.959963984540054 >>> st.norm.lsf(0.975)# 標(biāo)準(zhǔn)正態(tài)分布在 0.025 處的生存函數(shù)的逆函數(shù)值 1.959963984540054
對(duì)于非標(biāo)準(zhǔn)正態(tài)分布，通過更改參數(shù) loc 與 scale 來改變均值與標(biāo)準(zhǔn)差：
>>> st.norm.cdf(0, loc=2, scale=1) # 均值為 2，標(biāo)準(zhǔn)差為 1 的正態(tài)分布在 0 處的累計(jì)分布概率值 0.022750131948179195
對(duì)于其他隨機(jī)分布，可能更改的參數(shù)不一樣，具體需要查官方文檔。下面我們舉一些常用分布的例子：
>>> st.binom.pmf(4, n=100, p=0.05) # 參數(shù)值 n=100, p=0.05 的二項(xiàng)分布在 4 處的概率密度值 0.17814264156968956 >>> st.geom.pmf(4, p=0.05) # 參數(shù)值 p=0.05 的幾何分布在 4 處的概率密度值 0.04286875 >>> st.poisson.pmf(2, mu=3) # 參數(shù)值 mu=3 的泊松分布在 2 處的概率密度值 0.22404180765538775 >>> st.chi2.ppf(0.95, df=10) # 自由度為 10 的卡方分布在 0.95 處的逆函數(shù)值 18.307038053275146 >>> st.t.ppf(0.975, df=10) # 自由度為 10 的 t 分布在 0.975 處的逆函數(shù)值 2.2281388519649385 >>> st.f.ppf(0.95, dfn=2, dfd=12) # 自由度為 2, 12 的 F 分布在 0.95 處的逆函數(shù)值 3.8852938346523933
補(bǔ)充拓展：給定概率密度,生成隨機(jī)數(shù) python實(shí)現(xiàn)
實(shí)現(xiàn)的方法可以不止一種：
rejection sampling
invert the cdf
Metropolis Algorithm (MCMC)
本篇介紹根據(jù)累積概率分布函數(shù)的逆函數(shù)(2：invert the CDF)生成的方法。
自己的理解不一定正確，有錯(cuò)誤望指正。
目標(biāo)：
已知 y=pdf(x），現(xiàn)想由給定的pdf, 生成對(duì)應(yīng)分布的x
PDF是概率分布函數(shù)，對(duì)其積分或者求和可以得到CDF（累積概率分布函數(shù)），PDF積分或求和的結(jié)果始終為1
步驟（具體解釋后面會(huì)說）：
1、根據(jù)pdf得到cdf
2、由cdf得到inverse of the cdf
3、對(duì)于給定的均勻分布[0,1),帶入inverse cdf，得到的結(jié)果即是我們需要的x
求cdf逆函數(shù)的具體方法：
對(duì)于上面的第二步，可以分成兩類：
1、當(dāng)CDF的逆函數(shù)好求時(shí)，直接根據(jù)公式求取，
2、反之當(dāng)CDF的逆函數(shù)不好求時(shí)，用數(shù)值模擬方法
自己的理解：為什么需要根據(jù)cdf的逆去獲得x？
原因一：
因?yàn)閏df是單調(diào)函數(shù)因此一定存在逆函數(shù)（cdf是s型函數(shù)，而pdf則不一定，例如正態(tài)分布，不單調(diào)，對(duì)于給定的y，可能存在兩個(gè)對(duì)應(yīng)的x，就不可逆）
原因二：
這僅是我自己的直觀理解，根據(jù)下圖所示（左上為pdf，右上為cdf）
由步驟3可知，我們首先生成[0，1)的均勻隨機(jī)數(shù)，此隨機(jī)數(shù)作為cdf的y，去映射到cdf的x（若用cdf的逆函數(shù)表示則是由x映射到y(tǒng)），可以參考上圖的右上，既然cdf的y是均勻隨機(jī)的，那么對(duì)于cdf中同樣范圍的x，斜率大的部分將會(huì)有更大的機(jī)會(huì)被映射，因?yàn)閷?duì)應(yīng)的y范圍更大（而y是隨即均勻分布的），那么，cdf的斜率也就等同于pdf的值，這正好符合若x的pdf較大，那么有更大的概率出現(xiàn)（即重復(fù)很多次后，該x會(huì)出現(xiàn)的次數(shù)最多）
代碼實(shí)現(xiàn)——方法一，公式法
import numpy as np
import math
import random
import matplotlib.pyplot as plt
import collections
count_dict = dict()
bin_count = 20
def inverseCDF():
"""
return the x value in PDF
"""
uniform_random = random.random()
return inverse_cdf(uniform_random)
def pdf(x):
return 2 * x
# cdf = x^2, 其逆函數(shù)很好求，因此直接用公式法
def inverse_cdf(x):
return math.sqrt(x)
def draw_pdf(D):
global bin_count
D = collections.OrderedDict(sorted(D.items()))
plt.bar(range(len(D)), list(D.values()), align='center')
# 因?yàn)橛成鋌in的時(shí)候采用的floor操作，因此加上0.5
value_list = [(key + 0.5) / bin_count for key in D.keys()]
plt.xticks(range(len(D)), value_list)
plt.xlabel('x', fontsize=5)
plt.ylabel('counts', fontsize=5)
plt.title('counting bits')
plt.show()
for i in range(90000):
x = inverseCDF()
# 用bin去映射，否則不好操作
bin = math.floor(x * bin_count) # type(bin): int
count_dict[bin] = count_dict.get(bin, 0) + 1
draw_pdf(count_dict)
結(jié)果：
代碼實(shí)現(xiàn)——方法二，數(shù)值法
數(shù)值模擬cdf的關(guān)鍵是創(chuàng)建lookup table，
table的size越大則結(jié)果越真實(shí)（即區(qū)間劃分的個(gè)數(shù)）
import numpy as np
import math
import random
import matplotlib.pyplot as plt
import collections
lookup_table_size = 40
CDFlookup_table = np.zeros((lookup_table_size))
count_dict = dict()
bin_count = 20
def inverse_cdf_numerically(y):
global lookup_table_size
global CDFlookup_table
value = 0.0
for i in range(lookup_table_size):
x = i * 1.0 / (lookup_table_size - 1)
value += pdf2(x)
CDFlookup_table[i] = value
CDFlookup_table /= value # normalize the cdf
if y < CDFlookup_table[0]:
t = y / CDFlookup_table[0]
return t / lookup_table_size
index = -1
for j in range(lookup_table_size):
if CDFlookup_table[j] >= y:
index = j
break
# linear interpolation
t = (y - CDFlookup_table[index - 1]) / \
(CDFlookup_table[index] - CDFlookup_table[index - 1])
fractional_index = index + t # 因?yàn)閕ndex從0開始,所以不是 (index-1)+t
return fractional_index / lookup_table_size
def inverseCDF():
"""
return the x value in PDF
"""
uniform_random = random.random()
return inverse_cdf_numerically(uniform_random)
def pdf2(x):
return (x * x * x - 10.0 * x * x + 5.0 * x + 11.0) / (10.417)
def draw_pdf(D):
global bin_count
D = collections.OrderedDict(sorted(D.items()))
plt.bar(range(len(D)), list(D.values()), align='center')
value_list = [(key + 0.5) / bin_count for key in D.keys()]
plt.xticks(range(len(D)), value_list)
plt.xlabel('x', fontsize=5)
plt.ylabel('counts', fontsize=5)
plt.title('counting bits')
plt.show()
for i in range(90000):
x = inverseCDF()
bin = math.floor(x * bin_count) # type(bin): int
count_dict[bin] = count_dict.get(bin, 0) + 1
draw_pdf(count_dict)
真實(shí)函數(shù)與模擬結(jié)果
擴(kuò)展：生成伯努利、正太分布
import numpy as np
import matplotlib.pyplot as plt
"""
reference:
https://blog.demofox.org/2017/07/25/counting-bits-the-normal-distribution/
"""
def plot_bar_x():
# this is for plotting purpose
index = np.arange(counting.shape[0])
plt.bar(index, counting)
plt.xlabel('x', fontsize=5)
plt.ylabel('counts', fontsize=5)
plt.title('counting bits')
plt.show()
# if dice_side=2, is binomial distribution
# if dice_side>2 , is multinomial distribution
dice_side = 2
# if N becomes larger, then multinomial distribution will more like normal distribution
N = 100
counting = np.zeros(((dice_side - 1) * N + 1))
for i in range(30000):
sum = 0
for j in range(N):
dice_result = np.random.randint(0, dice_side)
sum += dice_result
counting[sum] += 1
# normalization
counting /= np.sum(counting)
plot_bar_x()
以上這篇python 計(jì)算概率密度、累計(jì)分布、逆函數(shù)的例子就是小編分享給大家的全部?jī)?nèi)容了，希望能給大家一個(gè)參考，也希望大家多多支持腳本之家。
相關(guān)文章
python實(shí)現(xiàn)一般游戲的自動(dòng)點(diǎn)擊具體操作
這篇文章主要介紹了python實(shí)現(xiàn)一般游戲的自動(dòng)點(diǎn)擊,本文給大家分享具體操作代碼及需要的軟件，需要的朋友可以參考下2021-10-10
這篇文章主要為大家介紹了python中的變量，具有一定的參考價(jià)值，感興趣的小伙伴們可以參考一下，希望能夠給你帶來幫助2021-12-12
Python使用Socket(Https)Post登錄百度的實(shí)現(xiàn)代碼
以前都是用一些高級(jí)模塊，封裝的比較好，今天嘗試使用socket模塊登錄百度，弄了半天才弄好，主要由于百度在登陸頁使用了https，我們需要對(duì)socket進(jìn)行一定處理2012-05-05
Python的NLTK模塊詳細(xì)介紹與實(shí)戰(zhàn)案例
自然語言處理庫NLTK在Python中的應(yīng)用廣泛,提供了分詞、詞性標(biāo)注、句法分析等多種功能,本文介紹了NLTK的核心功能、基本概念以及通過具體實(shí)戰(zhàn)案例（如文本分詞、去除停用詞、詞干提取等）展示了其在NLP任務(wù)中的實(shí)際應(yīng)用2024-09-09
這篇文章主要介紹了python發(fā)送郵件接收郵件示例,大家參考使用吧2014-01-01
解決python os.mkdir創(chuàng)建目錄失敗的問題
今天小編就為大家分享一篇解決python os.mkdir創(chuàng)建目錄失敗的問題，具有很好的參考價(jià)值，希望對(duì)大家有所幫助。一起跟隨小編過來看看吧2018-10-10
Python利用keyboard模塊實(shí)現(xiàn)鍵盤記錄操作
模擬鍵盤操作執(zhí)行自動(dòng)化任務(wù)，我們常用的有pyautowin等自動(dòng)化操作模塊。今天介紹的這個(gè)模塊叫做keyboard,它是純Python原生開發(fā)，編譯時(shí)完全不需要依賴C語言模塊。一行命令就能完成安裝，非常方便，需要的可以了解一下2022-10-10
如何在keras中添加自己的優(yōu)化器(如adam等)
這篇文章主要介紹了在keras中實(shí)現(xiàn)添加自己的優(yōu)化器(如adam等)方式，具有很好的參考價(jià)值，希望對(duì)大家有所幫助。一起跟隨小編過來看看吧2020-06-06
Python實(shí)現(xiàn)爬蟲IP負(fù)載均衡和高可用集群的示例代碼
做大型爬蟲項(xiàng)目經(jīng)常遇到請(qǐng)求頻率過高的問題,這里需要說的是使用爬蟲IP可以提高抓取效率,本文主要介紹了Python實(shí)現(xiàn)爬蟲IP負(fù)載均衡和高可用集群的示例代碼,感興趣的可以了解一下2023-12-12
pytorch1.0中torch.nn.Conv2d用法詳解
今天小編就為大家分享一篇pytorch1.0中torch.nn.Conv2d用法詳解，具有很好的參考價(jià)值，希望對(duì)大家有所幫助。一起跟隨小編過來看看吧2020-01-01