Оптимизация с помощью scipy_optimize

Question

Я пытаюсь оптимизировать функцию, используя curve_fit для scipy.optimize. Вот мой код.

import pandas as pd 
import numpy as np 
from scipy.optimize import curve_fit 

xdata = [row[0] for row in pd.read_excel("C:\\Users\\310967\\Desktop\\Scholar\\Wound Chelation Draft\\ChelationFiles.xlsx", sheetname="Case2Data",skiprows=0).as_matrix()] 
ydata = [row[1] for row in pd.read_excel("C:\\Users\\310967\\Desktop\\Scholar\\Wound Chelation Draft\\ChelationFiles.xlsx", sheetname="Case2Data",skiprows=0).as_matrix()] 
SF = [row[4] for row in pd.read_excel("C:\\Users\\310967\\Desktop\\Scholar\\Wound Chelation Draft\\ChelationFiles.xlsx", sheetname="Case2Data",skiprows=0).as_matrix()] 
uncertainty = [(np.sqrt(np.exp(np.log(a)**2)-1))*b for a,b in zip(SF, ydata)] 

Tau = [0,1,5,7] 



def func(x, I, E, ic1, ic2, ih1, ih2): 

    def iu(t): 
     return ((0.01295*np.exp(-0.645974*t))+(4.3688e-4*np.exp(-0.04251*t))+(5.642452e-5*np.exp(-0.00160863*t))) 

    def ic(t,tj): 
     if t > tj: 
      return ic1*np.exp(-0.693/ih1*(t-tj))+ic1*np.exp(-0.693/ih1*(t-tj)) 
     else: 
      return 0 

    def listofic(t): 
     list1 = [] 
     for tj in Tau: 
      list1.append(ic(t,tj)) 
     return list1 

    def Kj(tj): 
     return iu(tj+1)*(E-1)/(ic(1,0)-iu(tj+1)) 

    def listofKj(): 
     list2 = [] 
     for tj in Tau: 
      list2.append(Kj(tj)) 
     return list2 

    Kjs = listofKj() 

    def listofOneMinusKj(t): 
     list3 = [] 
     for a in Tau: 
      if t > a: 
       value = 1-Kj(a) 
      else: 
       value = 1 
      list3.append(value) 
     return list3 

    return (iu(x)*np.prod(listofOneMinusKj(x))+sum([a*b for a,b in zip(Kjs,listofic(x))]))*I 


popt, pcov = curve_fit(func, xdata, ydata, sigma=uncertainty) 
print(popt) 

Когда я запускаю приведенный выше код, я получаю ошибку о том, что «значение истинности массива с более чем одним элементом является неоднозначным. Используйте a.any() или a.All()». Это относится к части «if t> a» в списке функцийOneMinusKj (t).

Однако, если я запустил следующий код, несмотря на наличие «if t> a», код работает, как я ожидаю. Мне было интересно, что проблема с кодом выше.

import numpy as np 
Tau = [0,1,5,7] 

def func(x, I, E, ic1, ic2, ih1, ih2): 

    def iu(t): 
     return ((0.01295*np.exp(-0.645974*t))+(4.3688e-4*np.exp(-0.04251*t))+(5.642452e-5*np.exp(-0.00160863*t))) 

    def ic(t,tj): 
     if t > tj: 
      return ic1*np.exp(-0.693/ih1*(t-tj))+ic1*np.exp(-0.693/ih1*(t-tj)) 
     else: 
      return 0 

    def listofic(t): 
     list1 = [] 
     for tj in Tau: 
      list1.append(ic(t,tj)) 
     return list1 

    def Kj(tj): 
     return iu(tj+1)*(E-1)/(ic(1,0)-iu(tj+1)) 

    def listofKj(): 
     list2 = [] 
     for tj in Tau: 
      list2.append(Kj(tj)) 
     return list2 

    Kjs = listofKj() 

    def listofOneMinusKj(t): 
     list3 = [] 
     for a in Tau: 
      if t > a: 
       value = 1-Kj(a) 
      else: 
       value = 1 
      list3.append(value) 
     return list3 

    return (iu(x)*np.prod(listofOneMinusKj(x))+sum([a*b for a,b in zip(Kjs,listofic(x))]))*I 

print(func(1,400,12.5,0.99,0.01,0.55,10)) 

Answer 1

Один из способов отладки такого рода проблем является бросить import pdb; pdb.set_trace() колдовство чуть выше линии точки отладочных в. Затем запустите свой код, и он остановится в точке останова --- и вы сможете интерактивно исследовать различные объекты и пройти через код. Здесь вы, скорее всего, обнаружите, что либо t, либо a - это массив numpy при вызове curve_fit, и вы не можете просто сделать if array > another_array.

Answer 2

Scipy Optimizes Curve Fitting Procedure пытается передать полный вектор xdata в вашу функцию func. Вы передаете это listofOneMinusKj.

Поэтому t > a (или прошел как x > a) создает вектор bools. Чем срабатывает следующее сообщение об ошибке:

«Значение истинности массива с более чем одним элементом является неоднозначным Используйте a.any() или a.All()»

Эта ошибка возникает потому что вы не можете проверить, соответствует ли вектор length > 1. Как и предлагалось, вы можете использовать (t > a).any() для проверки, если любое значение t больше, чем a или (t > a).all(), чтобы проверить, все ли значения t больше.