ODE fitting with Python - least squares takes to long to compute

164 views Asked by At

I'm trying to replicate the results of "Modeling Plasma Virus Concentration during Primary HIV Infection", by fitting the data sets to my ODE. Essencially I need to fit 4 parameters (d, k, delta, p) to a dif. equation given a data set.

Up until now I got the following code done (I'm a begginer):

c = 3 #global constant 
v0= 10e-6

#dados
timeDataPatient1=np.array([0,22,43,78,106,146,183,230,268,358,435,489,519,534,584,610,687,778])
VirDataPatient1=10e3*np.array([v0*10e-3,27.2,210,85.9,81.1,46.2,60.1,82.8,103,72.1,79.4,70.4,207,42.6,10.8,54.2,22.3,40.8])


xdata = timeDataPatient1
ydata = VirDataPatient1

#inicial condition and max time of data set
IC = [10,0,10e-6]
t_max = max(xdata)

def diffeqs(t,y,d,k,delta,p):
    '''Model of the differencial equation'''
    global lam, c
    
    T=y[0]
    I=y[1]
    V=y[2]
    
    dT = 10*d - d*T-k*T*V
    dI = k*T*V - delta*I
    dV = p*I-c*V
    
    return [dT,dI,dV]

def Vtime(theta):
    '''solves the differencial equation, and returns the values of the model on the time points corresponding to the data set'''

    global t_max, IC, xdata
    
    sol = solve_ivp(diffeqs, (0,t_max), IC, args = (theta),t_eval=xdata)
    
    return sol.y[2]

#now I define the objetive function to minimize. For a parameter theta, that corresponds to d,k,delta,p in the ODE model,it calculates the difference between the log of the data given and the predicts by the ODE model.
  
def calcErrorHIV(theta):
    '''objetive function given in the paper'''
    global ydata
    
    dif = np.log(ydata)-np.log(Vtime(theta))
    
    return sum(dif**2)

#should return theta parameters that best fits the data, but It doesn't compute.
sp.optimize.least_squares(calcErrorHIV,[0.013,0.46e-2,0.40,980])

but when I try to run the otimize or least_squares function I get no result.

Does my code have a problem or the 4 parameters problem takes a long time to solve?

Thanks for the help.

Edit: I tried the inicial guess [0.013,0.46 10*(-3),0.40,980]

Edit2: Added extra info about the problem.

0

There are 0 answers