如何在Python中使用真实高斯分布拟合曲线？

我有一组近似高斯正态分布的散点，网上多数Python拟合示例用含C参数的Gauss1函数配合curve_fit，拟合效果不错，但Gauss1不是标准高斯分布。用标准高斯分布函数Gauss2调用curve_fit时会报错：OptimizeWarning: Covariance of the parameters could not be estimated，拿不到拟合参数。

相关代码片段 #

含缩放参数的高斯函数（Gauss1） #

def Gauss1(X, C, mu, sigma):
    return C * np.exp(-(X-mu) ** 2 / (2 * sigma ** 2))

导入curve_fit #

from scipy.optimize import curve_fit

标准高斯PDF（Gauss2） #

def Gauss2(X, mu, sigma):
    return (1 / (sigma * np.sqrt(2 * np.pi))) * np.exp(-(X-mu) ** 2 / (2 * sigma ** 2))

报错的拟合代码 #

popt2, pcov2 = curve_fit(Gauss2, x, y)

OptimizeWarning: Covariance of the parameters could not be estimated

完整测试代码（基于FIDE 2024年10月棋手评级数据） #

import matplotlib.pyplot as plt # Follow the convention
import numpy as np              # Library for working with arrays
from scipy.optimize import curve_fit

######################################## G L O B A L S
# X the ratings after grouping by 23 the chess players from FIDE 2024 Oct Standard
# Y is the average of players in the averaged ratings
XY = (
    (1411, 231), (1434, 271), (1457, 281), (1480, 287), 
    (1503, 292), (1526, 298), (1549, 293), (1572, 299), 
    (1595, 300), (1618, 303), (1641, 304), (1664, 308), 
    (1687, 301), (1710, 311), (1733, 304), (1756, 291), 
    (1779, 286), (1802, 283), (1825, 279), (1848, 268), 
    (1871, 260), (1894, 243), (1917, 229), (1940, 221), 
    (1963, 203), (1986, 169), (2009, 134), (2032, 112), 
    (2055, 102), (2078,  94), (2101,  81), (2124,  76), 
    (2147,  70), (2170,  61), (2193,  54), (2216,  45), 
    (2239,  42), (2262,  35), (2285,  31), (2308,  27), 
    (2331,  26), (2354,  22), (2377,  20), (2400,  18), 
    (2423,  15), (2446,  10), (2469,  10), (2492,   7), 
    (2515,   6), (2538,   5), (2561,   4), (2584,   3), 
    (2607,   3), (2630,   2), (2653,   2), (2676,   2), 
    (2699,   1), (2722,   2), (2745,   1), (2768,   1), 
    (2791,   1), (2837,   1))

class c:  # Constants
    X_LABEL_STEP = 50
    Y_LABEL_STEP = 10
    A4 = (11.69, 8.27)  # W x H in inches

###################################### F U N C T I O N S

def Graph_Gauss():
    global c, XY

    X = [] ; Y = []
    for e in XY:
        X += [e[0]]
        Y += [e[1]]

    # Set up the limits for X, ratings
    minX = min(X) ; maxX = max(X)
    # Set up the limits for Y, number of players with that rating
    minY = min(Y) ; maxY = max(Y)
        
    fig, ax = plt.subplots()
    fig.set_size_inches(c.A4)
    fig.suptitle("Rating Distribution", fontsize=18, y=0.935)

    x1  = []  # Make the x axle
    xlb = []  # Make the labels for x
    stp = c.X_LABEL_STEP
    for k in range(stp*(minX//stp), stp*(maxX//stp + 2), stp):
        x1 += [k]
        xlb.append(f"{k}")
    ax.set_xticks(ticks=x1, labels=xlb, rotation=270)
    ax.set_xlabel("Rating", fontsize=15, labelpad=7)

    stp = c.Y_LABEL_STEP
    yticks = np.arange(stp*(minY//stp), stp*(maxY//stp + 2), stp)
    ax.set_ylabel("Number of Players", fontsize=15, rotation=270, labelpad=18)
    ax.set_yticks(ticks=yticks)
    ax.grid(which="major", axis="both")

    ax.scatter(X, Y, color="#000FFF", marker='o', s=14)

    # plt.show()

    # Fit a normal distribution, aka Gaussian fitting curve
    # mu = mean = sum(x) / len(x) ; In our case sum(x * y) / sum(y)
    # sigma = standard deviation = sqrt((sum((x - mean)**2) / len(x))
    #
    # 1/(sigma*sqrt(2*pi)) * e**(-(x - mean)**2 / (2 * sigma**2))
    # Calculating the Gaussian PDF values given Gaussian parameters and random variable X
    def Gauss1(X, C, mu, sigma):
        return C * np.exp(-(X-mu)**2 / (2 * sigma**2))
    def Gauss2(X, mu, sigma):
        return (1/(sigma*np.sqrt(2*np.pi))) * np.exp(-(X-mu)**2 / (2 * sigma**2))

    x = np.array(X)
    y = np.array(Y)
    mu    = sum(x * y) / sum(y)                  
    sigma = np.sqrt(sum(y*(x - mu)**2)/sum(y))

    print(f"{1/(sigma*np.sqrt(2*np.pi))=:.2f} {mu=:.2f} {sigma=:.2f}")

    popt1, pcov1 = curve_fit(Gauss1, x, y, p0=[max(y), mu, sigma], maxfev=5000)
    popt2, pcov2 = curve_fit(Gauss2, x, y)  # This generates:
    # OptimizeWarning: Covariance of the parameters could not be estimated

    yg1 = Gauss1(x, *popt1)
    yg2 = Gauss2(x, *popt2)
    
    ax.plot(x, yg1, color="#FF0F00", linewidth=3,
            label=f"Normal: mu={popt1[1]:.2f}, sigma={popt1[2]:.2f}")
    plt.legend(fontsize=12)

    breakpoint()  # DEBUG

    plt.show()
    pass  # To set a breakpoint

#######################################################################
if __name__ == '__main__':
    # breakpoint()  # DEBUG, to set other breakpoints
    Graph_Gauss()

补充疑问 #

1. 高斯分布函数的准确定义是什么，其中sigma的位置在哪里？
2. 使用Gauss2时，将y值归一化（如除以230000）能得到较好拟合，但不清楚这类数值的正确推导方法，求解决方案或参考资料。

1. 高斯分布的准确定义及sigma的位置 #

高斯分布（正态分布）的**概率密度函数（PDF）**标准定义为：
$$f(x) = \frac{1}{\sigma\sqrt{2\pi}} e^{{-\frac{(x-\mu)}2}{2\sigma^2}}$$
其中：

• $\mu$：分布的均值，决定曲线的中心位置
• $\sigma$：分布的标准差，决定曲线的“胖瘦”——$\sigma$越大，曲线越扁平；$\sigma$越小，曲线越尖锐
• 分母的$\sigma$是标准差，和$\sqrt{2\pi}$共同保证整个PDF在$(-\infty,+\infty)$上的积分等于1（即总概率为1）

你写的Gauss1是缩放后的高斯函数，相当于把标准PDF乘以一个常数$C$，这个$C$用来匹配数据的幅值（比如你的数据是玩家数量，不是概率密度），所以它不是标准PDF，但适合拟合非归一化的观测数据。

2. 用标准高斯PDF拟合观测数据的解决方案 #

你遇到的报错，核心原因是：标准高斯PDF的输出值是概率密度（范围通常很小，比如你的数据里计算出的峰值约0.001），但你的原始y值是玩家数量（峰值300+），两者量级差了5个数量级，curve_fit无法找到合适的参数收敛路径，所以报协方差无法估计的错误。

正确的归一化方法 #

要让标准高斯PDF能拟合你的数据，需要把y值转换成概率密度，步骤如下：

1. 计算数据的总样本量：你的XY数据是按23分组后的平均玩家数，总样本量等于所有Y值乘以分组宽度（23）的总和，即：

   bin_width = 23
   total_players = sum(y * bin_width for y in Y)
   

2. 将原始y值转换为概率密度：概率密度 = （分组内玩家数）/（总样本量 × 分组宽度）

   y_pdf = np.array(Y) / (total_players * bin_width)
   

   这样转换后的y_pdf就和标准高斯PDF的输出量级匹配了，此时用curve_fit拟合Gauss2就不会报错。

完整修正后的拟合代码片段 #

# 计算总玩家数和分组宽度
bin_width = 23
total_players = sum(Y) * bin_width  # 每个Y是分组内的平均玩家数，乘23得分组总人数，再求和
y_pdf = np.array(Y) / (total_players * bin_width)

# 拟合标准高斯PDF，传入初始值加速收敛
popt2, pcov2 = curve_fit(Gauss2, x, y_pdf, p0=[mu, sigma])

# 如果要把拟合后的PDF转回玩家数量，只需反向计算：
yg2_players = Gauss2(x, *popt2) * total_players * bin_width
ax.plot(x, yg2_players, color="#00FF00", linewidth=2, label=f"Standard Gaussian: mu={popt2[0]:.2f}, sigma={popt2[1]:.2f}")

为什么初始值很重要 #

即使归一化后，最好给curve_fit传入参数初始值p0=[mu, sigma]（你已经提前计算了这两个值），这样能大幅提升拟合的成功率和速度，避免因参数搜索范围过大导致的收敛问题。

内容的提问来源于stack exchange，提问作者Vilmos Foltenyi