pythonnumpyscipysage

"only length-1 arrays can be converted to Python scalars" using scipy.optimize in Sage


I want to adjust the parameters of a model to a given set of data.

I'm trying to use scipy's function curve_fit in Sage, but I keep getting

TypeError: only length-1 arrays can be converted to Python scalars

Here´s my code:

from numpy import cos,exp,pi
f = lambda x: exp(  - 1 / cos(x)  )
import numpy as np

def ang(time): return (time-12)*pi/12

def temp(x,maxtemp):
    cte=(273+maxtemp)/f(0)**(1/4)
    if 6<x and x<18:
        return float(cte*f(ang(x))**(1/4)-273)
    else:
        return -273

lT=list(np.linspace(15,40,1+24*2))
lT=[float(num) for num in lT]             #list of y data

ltimes=np.linspace(0,24,6*24+1)[1:]
ltimes=list(ltimes)                     #list of x data
u0=lT[0]

def u(time,maxtemp,k):                    #the function I want to fit to the data
    def integ(t): return k*exp(k*t)*temp(t,maxtemp)
    return exp(-k*time)*( numerical_integral(integ, 0, time)[0] + u0 )

import scipy.optimize as optimization

print optimization.curve_fit(u, ltimes, lT,[1000,0.0003])

Solution

  • scipy.optimize.curve_fit expects the model function to be vectorized: that is, it must be able to receive an array (ndarray, to be precise), and return an array of values. You can see the problem right away by adding a debug printout

    def u(time,maxtemp,k):                    
        print time                % for debugging
        def integ(t): return k*exp(k*t)*temp(t,maxtemp)
        return exp(-k*time)*( numerical_integral(integ, 0, time)[0] + u0 )
    

    The output of print will be the entire array ltimes you are passing to curve_fit. This is something numerical_integral is not designed to handle. You need to give it values one by one.

    Like this:

    def u(time,maxtemp,k):                    
        def integ(t): return k*exp(k*t)*temp(t,maxtemp)
        return [exp(-k*time_i)*( numerical_integral(integ, 0, time_i)[0] + u0) for time_i in time]
    

    This will take care of the “only length-1 arrays can be converted" error. You will then have another one, because your lists ltimes and lT are of different length, which doesn't make sense since lT is supposed to be the target outputs for inputs ltimes. You should revise the definitions of these arrays to figure out what size you want.