Memorandum and Diary

2017年はザコい英語力・カスい筋力をどうにかする

逆関数法による任意の確率密度関数の生成

import numpy as np
from scipy import integrate
import math
import matplotlib.pyplot as plt

"""
python3
"""

#dense 未使用
def dense(z, thresh):
    if z <= thresh:
        return integrate.quad(lambda x: 2*x/thresh, 0, z)[0]
    else:
        return integrate.quad(lambda x: 2*x/thresh, 0, thresh)[0] + integrate.quad(lambda x: 2*(x-1)/(thresh-1), thresh, z)

#累積分布関数の逆関数
def reverse_dense(x):
    thresh = 0.7
    if x <= thresh:
        z = math.sqrt(thresh*x)
    else:
        z = -math.sqrt((thresh-1)*(x-1))+1
    return z

#一様分布からサンプリング
x = np.random.rand(100000)

#逆関数法により任意の分布を生成
z = np.array(list(map(reverse_dense, x)))

#ヒストグラム作成
plt.hist(z, bins=1000)   
plt.savefig('probability_density_function.png')

result

f:id:TheTree:20170613143710p:plain

参考

Rで学ぶ逆変換(逆関数)法