向SciPy curve_fit传入bounds触发TypeError的原因及解决方法

问题描述 #

使用Python脚本拟合实验数据时，SciPy的curve_fit函数在不传入bounds参数时能正常执行并返回拟合结果，但只要传入bounds（试过多种取值，包括bounds = (np.array([0, 0,0,0]), np.array([2500, 2500, 2500, 1000]))）就会触发错误：

TypeError: ufunc 'isfinite' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

实验数据范围为B=0.002到1（步长0.004，不含B=0），dSIGMA为负值，拟合函数来自H. Nakamura等人发表在《Nature Communications》11, 1161 (2020)的公式1。核心代码如下：

import pandas as pd
import numpy as np
from mpmath import zeta
from scipy.optimize import curve_fit
import os

def sigmaWAL(B, l, lfi, lso, N):
    q = 1.602*10**-19 
    h = 6.626*10**-34
    hbar = h/(2*np.pi)
    
    lb = np.sqrt(hbar/(4*q*np.abs(B)))
    
    l = l*10**-9
    lfi = lfi*10**-9
    lso = lso*10**-9
    
    # 将mpmath.zeta转为numpy兼容版本
    npzeta = np.frompyfunc(zeta,2,1)
    
    zeta1 = npzeta(0.5, 0.5 + (lb**2)/(l**2))
    zeta2 = npzeta(0.5, 0.5 + (lb**2)/(lfi**2))
    zeta3 = npzeta(0.5, 0.5 + 4*(lb**2)/(lso**2)+(lb**2)/(lfi**2))

    return N*((q**2)/(4*np.pi*h*lb))*(2*zeta1 + zeta2 - 3*zeta3)

def sigmaWAL_normalized(B, l, lfi, lso, N):
    return sigmaWAL(B, l, lfi, lso, N) - sigmaWAL(0.0001, l, lfi, lso, N)

def WAL_fit(data):
    Bexp = data['B']
    dS = data['dSIGMA']
    p0 = [10, 50, 10, 1]

    parameters, covariance = curve_fit(
        sigmaWAL_normalized, Bexp, dS, p0,
        bounds = ([0, 0, 0, 0], [2500, 2500, 2500, 1000]),
        maxfev = 100000
    )    
    l = parameters[0]
    lfi = parameters[1]
    lso = parameters[2]
    N = parameters[3]

    return parameters, covariance, sigmaWAL_normalized(Bexp, l, lfi, lso, N)

for file in os.listdir(os.fsencode('data')):
    name = 'data/' + os.fsdecode(file)
    panda = pd.read_csv(name)
    fitting = WAL_fit(panda)

错误原因 #

1. 数组类型不兼容：np.frompyfunc(zeta,2,1)返回的是object类型的numpy数组，数组元素是mpmath的数值对象，而非原生numpy数值类型。
2. 优化器行为差异：不指定bounds时，curve_fit默认使用Levenberg-Marquardt（lm）优化器，对输出类型兼容性较强；指定bounds时，默认切换到Trust Region Reflective（trf）或Dogbox优化器，这类优化器会对函数输出执行isfinite检查，而object数组无法被numpy的isfinite ufunc处理，因此触发TypeError。
3. 潜在除零风险：bounds中参数下限设为0（比如l=0），会导致(lb**2)/(l**2)出现除零错误，进一步加剧数值类型异常。

解决方案 #

1. 转换数组类型为原生numpy数值类型 #

在sigmaWAL函数中，将npzeta返回的object数组转换为np.float64类型：

# 替换原zeta计算部分
zeta1 = npzeta(0.5, 0.5 + (lb**2)/(l**2)).astype(np.float64)
zeta2 = npzeta(0.5, 0.5 + (lb**2)/(lfi**2)).astype(np.float64)
zeta3 = npzeta(0.5, 0.5 + 4*(lb**2)/(lso**2)+(lb**2)/(lfi**2)).astype(np.float64)

2. 调整bounds下限避免除零 #

由于l、lfi、lso是物理长度，不可能为0，将bounds下限设为极小的正数（比如1e-6）：

# 修改WAL_fit中的bounds参数
bounds = (
    [1e-6, 1e-6, 1e-6, 1e-6],  # 参数下限，避免为0
    [2500, 2500, 2500, 1000]   # 参数上限不变
)

3. 可选：显式指定优化器方法 #

如果上述修改后仍有问题，可显式指定'dogbox'优化器，它对边界处理的兼容性较好：

parameters, covariance = curve_fit(
    sigmaWAL_normalized, Bexp, dS, p0,
    bounds = bounds,
    method='dogbox',
    maxfev = 100000
)

验证修改 #

修改后重新运行代码，curve_fit可正常处理bounds参数，不会触发isfinite相关的TypeError，同时避免了除零导致的数值异常。

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