如何在Python中拟合非对称高斯曲线？

我有一条非对称高斯曲线需要拟合，但在Python中无法完成拟合。原高斯函数可以正常运行，但非对称高斯曲线的拟合存在问题。数据从Excel文件读取，x（结合能）和y（强度）值位于不同工作表中。

代码中的核心问题 #

1. 非对称高斯函数参数结构不合理：原函数将所有参数打包为一个列表p，不利于优化函数识别和迭代，也增加了初始值设置的复杂度。
2. 初始值设置错误：原代码中initials的第一个参数使用了df.index（数组），这完全不符合振幅参数的数值要求，直接导致拟合不收敛。
3. 语法与调用错误：注释掉的curve_fit部分存在语法错误，参数传递逻辑混乱。

修正后的完整代码 #

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.special import erf
from scipy.optimize import leastsq
from scipy.optimize import curve_fit

# 定义非对称高斯函数：参数为x，振幅A，中心位置x0，宽度sigma，偏斜系数alpha
def asymGaussian(x, A, x0, sigma, alpha):
    amp = A / (sigma * np.sqrt(2 * np.pi))
    spread = np.exp(-(x - x0)**2 / (2 * sigma**2))
    skew = 1 + erf(alpha * (x - x0) / (sigma * np.sqrt(2)))
    return amp * spread * skew

# residuals函数适配新的参数格式
def residuals(p, y, x):
    return y - asymGaussian(x, *p)

# 设置绘图字体参数
params = {'legend.fontsize': 'large',
          'axes.labelsize': 'x-large',
          'axes.titlesize': 'x-large',
          'xtick.labelsize': 'x-large',
          'ytick.labelsize': 'x-large'}
plt.rcParams.update(params)
filename = "Au intensity.xlsx"

# 读取数据
var = pd.read_excel(filename, sheet_name=0).dropna()

# 读取结合能（Sheet2的E列）作为行索引并转换为浮点数组
BEs = pd.read_excel(filename, sheet_name=1, header=None)
var.index = BEs.values.reshape(var.shape[0])
x_data = var.index.values.astype(float)

# 读取角度（Sheet3的A列）作为列名
angles = pd.read_excel(filename, sheet_name=2, header=None)
var.columns = angles.values.reshape(var.shape[1])

angle_range = 6
vertical_shift = 250

fig, ax = plt.subplots(figsize=(5, 6))
for i in range(0, int(var.columns.max())+1, angle_range):
    df = var.loc[:, (var.columns >= i) & (var.columns < i+angle_range)]
    y_data = df.sum(axis=1).values.astype(float)
    
    # 设置合理的初始参数猜测
    peak_y = np.max(y_data)
    peak_x = x_data[np.argmax(y_data)]
    initial_sigma = 1.0  # 可根据数据调整的初始宽度
    initial_alpha = 1.0  # 初始偏斜系数
    # 振幅A的初始值由峰值反推：A ≈ peak_y * sigma * sqrt(2pi)
    initials = [peak_y * initial_sigma * np.sqrt(2*np.pi), peak_x, initial_sigma, initial_alpha]
    
    # 使用leastsq进行拟合
    cnsts, _ = leastsq(residuals, initials, args=(y_data, x_data))
    
    # 生成拟合曲线
    y_fit = asymGaussian(x_data, *cnsts)
    
    # 绘制原始数据与拟合曲线
    ax.plot(x_data, y_data + i * vertical_shift, 
            label=f"{i}°-{i+angle_range}°")
    ax.plot(x_data, y_fit + i * vertical_shift, linestyle='--', color='red')
    
    print(f"角度范围: {i}-{i+angle_range}° | 列数: {df.shape[1]}")
    print(f"拟合参数: 振幅A={cnsts[0]:.2f}, 中心x0={cnsts[1]:.2f}, 宽度sigma={cnsts[2]:.2f}, 偏斜alpha={cnsts[3]:.2f}\n")

ax.set_xlim(max(x_data), min(x_data))
ax.set_yticks([])
ax.set_xlabel("结合能 (eV)")
ax.set_ylabel("强度 (a.u.)")
ax.legend()
plt.show()

关键修正说明 #

• 函数重构：将非对称高斯函数的参数拆分为独立的A, x0, sigma, alpha，优化函数可以更高效地迭代每个参数。
• 初始值优化：基于数据的峰值位置和高度设置初始值，确保拟合从接近真实解的起点开始，避免迭代发散。
• 数据类型明确：将x、y数据转换为浮点数组，避免索引类型导致的计算错误。
• 可视化区分：给拟合曲线添加虚线和红色标识，便于对比原始数据与拟合结果。
• 参数输出：打印拟合得到的参数，方便验证拟合效果。

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