使用Scipy在Python中拟合函数至数据点失败的求助

我尝试在Python中将名为concave_func的自定义函数拟合到给定数据点，但未取得理想效果。该函数包含7个需估计的参数（C_1、C_2、alpha_one、alpha_two、I_x、nu_t、T_e），但仅有6个数据点。我已尝试两种拟合方法：1). scipy.optimize.minimize；2). scipy.optimize.curve_fit，但拟合曲线与数据点不匹配。此外，我需要约束alpha_one和alpha_two在2到3之间，且参数T_e不超过10000。

代码实现 #

import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt

frequency = np.array([22,45,150,408,1420,23000]) #x_values
b_temp = [2.55080863e+04, 4.90777800e+03, 2.28984753e+02, 2.10842949e+01, 3.58631166e+00, 5.68716056e-04] #y_values

#Defining the function that I want to fit
def concave_func(x, C_1, C_2, alpha_one, alpha_two, I_x, nu_t, T_e): 
    one = x**(-alpha_one)
    two = (C_2/C_1)*(x**(-alpha_two))
    three = I_x*(x**-2.1)
    expo = np.exp(-1*((nu_t/x)**2.1))
    eqn_one = C_1*(one + two + three)*expo
    eqn_two = T_e*(1 - expo)
    return eqn_one + eqn_two

#Defining chi_square function
def chisq(params, xobs, yobs):
    ynew = concave_func(xobs, *params)
    #yerr = np.sum((ynew- yobs)**2)
    yerr = np.sum(((yobs- ynew)/ynew)**2)
    print(yerr)
    return yerr

result = minimize(chisq, [1,2,2,2,1,1e6,8000], args = (frequency,b_temp),  method = 'Nelder-Mead', options = {'disp' : True, 'maxiter': 10000})
x = np.linspace(-300,24000,1000)
plt.yscale("log")
plt.xscale("log")
plt.plot(x,concave_func(x, *result.x))
print(result.x)
print(result)
plt.plot(frequency, b_temp, 'r*' )
plt.xlabel("log Frequency[MHz]")
plt.ylabel("log Temp[K]")
plt.title('log Temparature vs log Frequency')

plt.grid()
plt.savefig('the_plot_2060.png')

拟合结果 #

解决思路 #

• 参数与数据点数量不匹配：7个参数仅对应6个数据点，属于欠定问题，优化算法无法收敛到唯一合理解。要么补充更多数据点，要么简化函数（比如合并相似参数项，或根据物理背景固定部分参数）。
• 约束条件未生效：Nelder-Mead方法不支持参数约束，换用支持约束的优化器（如L-BFGS-B、TNC），通过bounds参数指定约束：

  bounds = [(None, None), (None, None), (2, 3), (2, 3), (None, None), (None, None), (None, 10000)]
  result = minimize(chisq, [1,2,2,2,1,1e6,8000], args=(frequency,b_temp), method='L-BFGS-B', bounds=bounds, options={'disp': True, 'maxiter':10000})
  

• 初始值设置不合理：nu_t初始值设为1e6，远大于最大频率23000，导致(nu_t/x)过大，expo趋近于0，函数行为被T_e主导。建议把nu_t初始值设为数据中的中间频率（如408或1420），其他参数也尽量贴合数据趋势设置。
• 损失函数存在缺陷：当前用的相对误差平方和在ynew接近0时会数值爆炸，导致优化不稳定。建议改用对数空间拟合（因为数据是对数刻度），计算对数后的误差平方和，或者使用加权均方误差。
• 绘图范围错误：x = np.linspace(-300,24000,1000)包含负数，而频率是正数，负数代入x**(-alpha)会产生复数，导致曲线不完整。应将x范围改为正数区间：x = np.linspace(10,24000,1000)。

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