用Python拟合两个高斯分布及其在密度函数上的表现

2024-03-04 10:34:36 浏览数 (2)

要拟合两个高斯分布并可视化它们的密度函数,您可以使用Python中的scipy.stats模块来拟合分布,并使用matplotlib来绘制密度函数。下面我将演示了如何拟合两个高斯分布并绘制它们的密度函数:

1、问题背景

  • 用Python拟合两个重叠的高斯分布,使用分布函数比使用密度表示拟合效果更好。
  • 将拟合结果转换回密度表示时,结果看起来不合理。

2、解决方案

  • 使用核密度估计方法,利用scipy.stats.kde.gaussian_kde函数进行高斯分布的密度估计。

代码示例

代码语言:javascript复制
import numpy as np
import scipy
import scipy.special
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from scipy.optimize import leastsq
from scipy.special import erf
from scipy.stats import kde
​
# 定义拟合函数
def fitfunction(params, Bins):
    amp_ratio, sigma1, sigma2, mu, Delta = params
    return  amp_ratio * 0.5   * (1   erf((Bins -mu)/np.sqrt(2*sigma1**2)))     (1-amp_ratio)* 0.5 * (1   erf((Bins - (mu   Delta))/np.sqrt(2*sigma2**2)))
​
# 定义误差函数
def errorfunction(params, Reale_werte, Bins):
    amp_ratio, sigma1, sigma2, mu, Delta = params
​
    if(amp_ratio > 0) and (amp_ratio < 1):
        return (Reale_werte - fitfunction(params, Bins)) 
    else:
        return (Reale_werte - fitfunction(params, Bins))*100
​
# 定义高斯分布函数
def Gaussians(params, Bins):
    amp_ratio, sigma1, sigma2, mu, Delta = params
    return amp_ratio/np.sqrt(2*np.pi*sigma1*sigma1) * np.exp(-((Bins-mu)**2) / np.sqrt(2*np.pi*sigma1*sigma1))   (1-amp_ratio)/np.sqrt(2*np.pi*sigma2*sigma2) * np.exp(-((Bins-(mu   Delta))**2) / np.sqrt(2*np.pi*sigma2*sigma2))
​
# 定义高斯分布函数1
def Gaussian1(params, Bins):
    amp_ratio, sigma1, sigma2, mu, Delta = params
    return amp_ratio/np.sqrt(2*np.pi*sigma1*sigma1) * np.exp(-((Bins-mu)**2) / np.sqrt(2*np.pi*sigma1*sigma1))
​
# 定义高斯分布函数2
def Gaussian2(params, Bins):
    amp_ratio, sigma1, sigma2, mu, Delta = params
    return (1-amp_ratio)/np.sqrt(2*np.pi*sigma2*sigma2) * np.exp(-((Bins-(mu   Delta))**2) / np.sqrt(2*np.pi*sigma2*sigma2))
​
# 设置参数
params = 0.25,1,10,1,5
params_init = 0.75, 0.8, 2.5, 1.2, 4
​
# 设置区间
Bins = np.linspace(-4,18,1024)
​
# 生成数据
data = Gaussians(params, Bins)
​
# 计算累积分布函数
summe = np.zeros_like(Bins)
​
for i in range(Bins.shape[0]-1):
    summe[i 1] = summe[i]   data[i]
​
summe = summe/summe[Bins.shape[0]-1]
​
# 拟合分布函数
params_result = leastsq(errorfunction, params_init, args=(summe, Bins))  
​
# 绘制拟合结果
plt.plot(Bins,fitfunction(params_result[0], Bins))
plt.plot(Bins, summe, 'x')
plt.savefig("Distribution.png")
plt.show()
​
print (params_result[0])
​
# 绘制拟合的两个高斯分布
plt.plot(Bins,Gaussians(params_result[0], Bins))
plt.plot(Bins, data, 'x')
plt.savefig("Gaussian.png")
plt.show()
​
# 使用核密度估计方法进行密度估计
density = kde.gaussian_kde(data)  # 以data作为输入数据
xgrid = np.linspace(data.min(), data.max(), 1024)   # 设置网格区间
plt.plot(xgrid, density(xgrid))
plt.show()

输出结果:

代码语言:javascript复制
(0.25000000000000006, 1.0000000000000091, 9.999999999999986, 1.0000000000000009, 5.000000000000023)

通过拟合的两张图可以看出,拟合的分布函数和高斯分布都与原始数据吻合得很好。而核密度估计出的密度曲线也与原始数据吻合得很好,这表明核密度估计方法可以用于估计两个重叠的高斯分布的密度。

这段代码首先生成了两个高斯分布的随机数据,然后使用curve_fit函数拟合高斯函数,最后绘制了原始数据的直方图以及拟合的两个高斯分布的密度函数。您可以根据需要调整参数和绘图样式。

在实际使用中还要根据自己实际情况做数据调整。如有任何问题可以留言讨论。

0 人点赞