采用二分法逼近寻找非标准正态分布函数中某概率对应的反函数值

        该方法用于寻找非标准正态分布函数的累积概率的反函数值，举例来说就是你给出一个概率比如0.05，该系统就会返回样本中的两个值，一个出现概率为0.05，一个出现概率为0.95，其作用可以用来求重现度或统计极值，其效果相当于标准正太分布中的 scipy. stats.norm.ppf(x,3,1)，但上述函数不支持非标准正态分布，同时在网上也未找到相关工具，因此编写集成代码以供使用。

        由于该系统使用二分法进行逼近寻找，所以在大样本数时其准确性和速度才具有优势。

from scipy import stats
from scipy.integrate import simps
import numpy as np

class pdf:
    # 采用二分法逼近寻找非标准正态分布函数的累积概率的反函数值, 对同一个值返回左右两个
    # 数值可用于求极值, 同时由于采用的是二分法, 因此对大量数据样本时有较好的计算效率。 
    # 但由于采用的是数值积分, 因此在样本数较少时误差会较大. 
    # self.x1 排列去重后的数据
    # self.gauss 为对应X1的高斯函数值
    # self.min 对应累计概率为prop的反函数值
    # self.max 对应累计概率为 (1-prop) 的反函数值
    # self.skew 偏度
    # self.kurtosis 偏度
    def __init__(self, x ,prob= 0.5):  #x为数据 prob为索求值对应的概率
        self.gauss = stats.norm.pdf(x, loc = np.mean(x), scale = np.std(x))
        data = []
        data.append(x)
        data.append(self.gauss)
        data = np.unique(data, axis=1)
        self.x1 = data[0]
        self.gauss = data[1]
        total = simps(self.gauss, self.x1)
        target_min = total * prob
        target_max = total * (1 - prob)
        self.skew = stats.skew(self.x1)
        self.kurtosis = stats.kurtosis(self.x1)
        
        def target(target):
            i1 = len(self.x1)
            if i1%2 == 0: i2 = int(i1/2)
            else: i2 = int((i1+1)/2)
            while True:
                a = simps(self.gauss[0:i1+1], self.x1[0:i1+1])
                b = simps(self.gauss[0:i2+1], self.x1[0:i2+1])
                if a  == target: return self.x1[i1]
                elif a > target:
                    if b > target:
                        i1 = i2
                        i2 = int(i1/2)
                        continue
                    elif b < target and i1 - i2 != 1:
                        i2 += int((i1-i2)/2)
                        continue
                    elif b < target and i1 - i2 == 1: 
                        return 0.5 * (self.x1[i1] + self.x1[i2])
                    elif b == target: return self.x1[i2]

        self.max = target(target_max)
        self.min = target(target_min)