0
$\begingroup$

I'm trying to estimate the parameters for the following pair of differential equations

meq = gamma*v[t]*m[t] - mu*m[t]; veq = v[t]*(k - epsilon - mu) - v[t]^2*k + m[t]*l - m[t]^2*l - m[t]*v[t]*(k + l + gamma); equations = {m'[t] == meq, v'[t] == veq};

where gamma, k, l and epsilon are parameters to be estimated, and mu is a known parameter with value 4.58*10^(-5). Also I have initial conditions given by

initials = {m[0] == 2.03, v[0] == 3.09}; mu = 4.58*10^(-5)

and the data for the function m is the following

data = {{60, 2.13597}, {300, 2.27247}, {390, 2.29472}, {420, 2.40096}, {630, 2.59312}, {660, 2.84918}, {780, 2.93677}, {1020, 3.02945}, {1110, 3.04794}, {1140, 3.05796}, {1140, 3.08739}, {1380, 3.21218}, {1380, 3.2873}, {1500, 3.3347}, {1680, 3.44467}, {1710, 3.47574}, {1710, 3.48421}, {1830, 3.50433}}

I defined a parametric solution for the pair of equations with 

pfun = ParametricNDSolveValue[Join[equations, initials], m, {t, 0, 2000}, {epsilon, k, l, gamma}]

And then I tried to use FindFit as

fit = FindFit[data, pfun[epsilon, k, l, gamma][t], {epsilon, k, l, gamma}, t]

However this method seems to be super inefficient as it seems to take unreasonable amount of time to complete the fit. Is there a way maybe to form the parametric function in a way that the FindFit works faster or use the FindFit more efficiently? Of course if you have another faster method in mind it's very welcome.

I also tried to solve this with some different programs than Mathematica and they seems to be working much faster but in those I have no idea what algorithms they use and also in them I have serious data issues.

Improve this question
$\endgroup$
2
  • 1
    $\begingroup$ "seems to be super inefficient" - note that you had not even attempted to provide starting guesses for epsilon, k, l, gamma; thus, FindFit[] just assigns them all a value of 1 as a guess and iterates accordingly, which you have seen to be intolerably slow. $\endgroup$ Commented Oct 2, 2018 at 11:05
  • $\begingroup$ @J.M.issomewhatokay. I also tried that but it seemed to be no help to start them from some values that should be closer to their "real" values. $\endgroup$ Commented Oct 2, 2018 at 11:11

1 Answer 1

Reset to default
2
$\begingroup$

Here you need to use a different model. In this model, too, not everything is smooth, but the result is clear and fast

meq = gamma*v[t]*m[t] - mu*m[t]; veq = v[t]*(k - epsilon - mu) - v[t]^2*k + m[t]*l - m[t]^2*l - m[t]*v[t]*(k + l + gamma); equations = {m'[t] == meq, v'[t] == veq}; initials = {m[0] == 2.03, v[0] == 3.09}; mu = 4.58*10^(-5); data = {{60, 2.13597}, {300, 2.27247}, {390, 2.29472}, {420, 2.40096}, {630, 2.59312}, {660, 2.84918}, {780, 2.93677}, {1020, 3.02945}, {1110, 3.04794}, {1140, 3.05796}, {1140, 3.08739}, {1380, 3.21218}, {1380, 3.2873}, {1500, 3.3347}, {1680, 3.44467}, {1710, 3.47574}, {1710, 3.48421}, {1830, 3.50433}}; model[epsilon_?NumberQ, k_?NumberQ, l_?NumberQ, gamma_?NumberQ] := (model[epsilon, k, l, gamma] = First[m /. NDSolve[{v'[t] == v[t]*(k - epsilon - mu) - v[t]^2*k + m[t]*l - m[t]^2*l - m[t]*v[t]*(k + l + gamma), m'[t] == gamma*v[t]*m[t] - mu*m[t], m[0] == 2.03, v[0] == 3.09}, {m, v}, {t, 2000}]]); fit = FindFit[data, model[epsilon, k, l, gamma][ t], {{epsilon, -4.8}, {k, 4.9}, {l, 0.}, {gamma, .1}}, t, PrecisionGoal -> 4, AccuracyGoal -> 4] (*Out[]= {epsilon -> -4.78195, k -> 4.89863, l -> -0.00336569, gamma -> 0.0769293}*) Show[ Plot[model[epsilon, k, l, gamma][t] /. fit, {t, 0, 2000}], ListPlot[data, PlotStyle -> Orange]]

fig1

Improve this answer
$\endgroup$
7
  • $\begingroup$ What do you mean not everything is smooth? And how is the model different? I see thar its fast but how is it different. $\endgroup$ Commented Oct 3, 2018 at 5:24
  • $\begingroup$ Not smooth because there is a message system that not everything is smooth. My code does not contain a parametric function that slows down the solution. $\endgroup$ Commented Oct 3, 2018 at 5:33
  • $\begingroup$ Okay, thank you. Still testing this trough more. $\endgroup$ Commented Oct 3, 2018 at 5:59
  • $\begingroup$ How did you choose the initial quesses {{epsilon, -4.8}, {k, 4.9}, {l, 0.}, {gamma, .1}} for the four parameters? Because the ones you selected seem to be much better than other ones that I'm trying. $\endgroup$ Commented Oct 3, 2018 at 13:41
  • $\begingroup$ @Kplusn I found them in two steps. First, I took all the parameters equal to 1, except for l = 0. Got a solution and took it as the next step. $\endgroup$ Commented Oct 3, 2018 at 14:18

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.