In my work, I'm trying to fit a model to data using a likelihood based method. I've had this code work before for slightly different models but for some reason it has been throwing this error when used with this model: "RuntimeError: The map-like callable must be of the form f(func, iterable), returning a sequence of numbers the same length as 'iterable'". I'm not sure if it's the model or something wrong with the code, but if you could help me understand what this error message means and how to fix it I would appreciate it.
import pandas as pd
import numpy as np
from scipy.integrate import odeint
from scipy.optimize import curve_fit
from lmfit import minimize, Parameters, Parameter, report_fit
#==============================================================================
''' Make new data matrix, same as csv except infected cells are one total for convience '''
DataFrame = pd.read_csv('Cell_count_data_Tromas_2014.csv') # Read data from file
TROMAS_DATA = np.empty((DataFrame.shape[0], 5), int)
for i in range(DataFrame.shape[0]): # Number of rows
for j in range(5): # Number of columns desired
TROMAS_DATA[i][j] = DataFrame.iloc[i, j]
if j == 4:
TROMAS_DATA[i][j] = DataFrame.ix[i, 'Venus_only'] + DataFrame.ix[i, 'BFP_only'] + DataFrame.ix[i, 'Mixed']
'''
Col 0: Days post infection
Col 1: Leaf number
Col 2: Replicate plant number
Col 3: Number of unifected cells
Col 4: Number of total infected cells
'''
#==============================================================================
''' Make axis for negative log likelihood '''
ZERO_DAYS_AXIS = [0, 3, 5, 7, 10]
#==============================================================================
''' Parameter list for model '''
likeParams = Parameters()
likeParams.add('I0', value = .00372, min = .0000001, max = 1.0000)
likeParams.add('b', value = .5, min = .0001, max = 20.0000)
likeParams.add('x5', value = .5, min = .0001, max = 20.0000)
likeParams.add('x6', value = .5, min = .0001, max = 20.0000)
likeParams.add('x7', value = .5, min = .0001, max = 20.0000)
likeParams.add('psi3', value = .5, min = .0001, max = 1.0000)
likeParams.add('psi5', value = .5, min = .0001, max = 1.0000)
likeParams.add('psi6', value = .5, min = .0001, max = 1.0000)
likeParams.add('psi7', value = .5, min = .0001, max = 1.0000)
#==============================================================================
def model(Mk, t, parameters):
M3 = Mk[0]
M5 = Mk[1]
M6 = Mk[2]
M7 = Mk[3]
try: # Get parameters
b = parameters['b'].value
x3 = parameters['x3'].value
x5 = parameters['x5'].value
x6 = parameters['x6'].value
psi3 = parameters['psi3'].value
psi5 = parameters['psi5'].value
psi6 = parameters['psi6'].value
psi7 = parameters['psi7'].value
except KeyError:
b, x3, x5, x6, psi3, psi5, psi6, psi7 = parameters
if (M3 < psi3):
S3 = (1 - (M3 / psi3))
else:
S3 = 0
if (M5 < psi5):
S5 = (1 - (M5 / psi5))
else:
S5 = 0
if (M6 < psi6):
S6 = (1 - (M6 / psi6))
else:
S6 = 0
if (M7 < psi7):
S7 = (1 - (M7 / psi7))
else:
S7 = 0
dM3dt = b * M3 * S3 + x3 * S3 * M7
dM5dt = b * M5 * S5 + x5 * S5 * M7
dM6dt = b * M6 * S6 + x6 * S6 * M7
dM7dt = b * M7 * S7
return [dM3dt, dM5dt, dM6dt, dM7dt]
#==============================================================================
''' Compute negative log likelihood of Tromas' data given the model, see eq. (3) pg. 11 '''
def negLogLike(parameters):
# Solve ODE system to get model values; parameters are not yet fitted
Lk0 = [0, 0, 0, parameters['I0'].value]
MM = odeint(model, Lk0, ZERO_DAYS_AXIS, args=(parameters,))
nll = 0
epsilon = 10**-10
for t in range(4): # Iterate through days
for p in range(5): # Iterate through replicates
for k in range(4): # Iterate through leaves
Vktp = TROMAS_DATA[20 * t + 4 * p + k][4] # Number of infected cells
Aktp = TROMAS_DATA[20 * t + 4 * p + k][3] + Vktp # Total number of cells observed
Iktp = MM[t + 1][k] # Frequency of cellular infection
if (Iktp <= 0):
Iktp = epsilon
elif (Iktp >= 1):
Iktp = 1 - epsilon
nll += Vktp * np.log(Iktp) + (Aktp - Vktp) * np.log(1 - Iktp)
return [-nll]
#==============================================================================
''' Miminize negative log likelihood with differential evolution algorithm '''
result_likelihood = minimize(negLogLike, likeParams, method = 'differential_evolution')
report_fit(result_likelihood)
I figured it out. I didn't define "x3" in likeParams. I guess that error propagated and got spat out as the differential evolution error.