Scipy curve_fit指数拟合遇溢出警告，求正确拟合方法

问题背景 #

• 采用scipy.optimize.curve_fit拟合含参数a、k、b的指数函数：

  def exponential(x, a, k, b):
      return a*np.exp(x*k) + b
  

• 数据集情况：
  • x_data：int64类型，取值范围1~44098，共44098个数据点
  • y_data：float32类型，取值范围-0.44383433~1.0，共44098个数据点

遇到的问题 #

• 用初始参数p0=[1, -0.5, 1]拟合后，得到结果[1.28765636e+00, -3.27620187e-04, -8.91012481e-02]，但触发RuntimeWarning：exp和multiply操作溢出
• 尝试将数据类型改为complex256，无法解决溢出问题
• 调整参数边界后警告消失，但拟合结果完全依赖边界设置：
  • 设置bounds=([1, -1e100, 1], [1e100, 1, 1e100])时，结果为[5.e+99, -5.e+99, 5.e+99]
  • 设置bounds=([1, -1e1000, 1], [1e1000, 1, 1e1000])时，结果为[1.00448191, -5.15058477, 1.0]
• 尝试使用Decimal模块，同样无效果

原因分析 #

1. 数值溢出根源：x的最大值为44098，拟合过程中参数k可能跳转为正数，导致x*k成为极大的正数，np.exp(x*k)直接溢出为无穷大，触发警告。
2. 参数可识别性问题：原指数函数的参数组合存在退化情况：
   • 当k为绝对值极大的负数时，exp(x*k)趋近于0，函数近似为常数b，此时a的取值不影响结果，导致a和b无法同时确定
   • 当k为绝对值极大的正数时，大x对应的exp(x*k)溢出，若设置a的下界为1，拟合会被迫取极端参数值来满足边界，结果完全失控
3. 初始参数不合理：初始参数p0=[1, -0.5, 1]对应的函数在x=1时取值约1.606，超出y_data的最大值1.0，与数据分布不匹配，加剧了拟合过程中的参数震荡。

解决方法 #

1. 重新参数化指数函数 #

调整函数形式，避免数值溢出并解决可识别性问题，推荐使用衰减型指数函数（符合y_data有上限1.0的特点）：

def exponential_decay(x, a, k, b):
    # k为正数，x增大时指数项衰减至0，函数趋近于b
    return b + a*np.exp(-k*x)

这种形式下，k恒为正，大x对应的指数项不会溢出，且参数意义更明确：b是x→∞时的渐近值（接近y的最大值1.0），a是初始偏差，k是衰减速率。

2. 缩放x值预处理 #

将x_data缩放至小范围，避免大x导致的指数爆炸：

x_scaled = x_data / 10000  # 缩放后x范围为0.0001~4.4098

用缩放后的x拟合，得到参数后再转换回原尺度：原函数中的k等于拟合得到的k除以缩放因子（这里是10000）。

3. 设置合理的参数边界与初始值 #

结合数据分布设置边界和初始值，比如：

# 初始值：b接近y最大值1.0，a为负数（因为y最小值低于1），k为小正数
p0 = [-1.5, 1e-4, 1.0]
# 边界：a在-10~0之间，k在1e-6~1e-2之间（保证衰减速率合理），b在0.9~1.1之间
bounds = ([-10, 1e-6, 0.9], [0, 1e-2, 1.1])

4. 改用鲁棒拟合工具 #

若数据存在异常值，使用scipy.optimize.least_squares代替curve_fit，选择鲁棒损失函数减少异常值影响：

from scipy.optimize import least_squares

def residual(params, x, y):
    a, k, b = params
    return y - (b + a*np.exp(-k*x))

result = least_squares(residual, p0=[-1.5, 1e-4, 1.0], bounds=bounds, loss='huber', args=(x_data, y_data))
fit_params = result.x

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