       Re: Nonlinear fit to diff eqs.

• To: mathgroup at smc.vnet.net
• Subject: [mg5086] Re: Nonlinear fit to diff eqs.
• From: rubin at msu.edu (Paul A. Rubin)
• Date: Fri, 25 Oct 1996 22:48:51 -0400
• Organization: Michigan State University
• Sender: owner-wri-mathgroup at wolfram.com

```In article <53gpue\$5am at dragonfly.wolfram.com>,
mod <c530789 at showme.missouri.edu> wrote:
->Dear All,
->
->My problem is to fit a set of data to a system of coupled first order
->ordinary differential equations. The example is:
->
->   Experimental data of y1, y2 and y3 as a function of time
->             y1 = {a1,a2,a3,a4,.....,aN}
->             y2 = {b1,b2,b3,b4,.....,bN}
->             y3 = {c1,c2,c3,c4,.....,cN}
->             time = {t1,t2,t3,t4,.....,tN}
->
->    Equations to fit (for example):
->
->             y1'[t] = k1 y1[t] + k2 (1-y1[t]-y2[t]) + k4 y3[t]
->             y2'[t] = k2 y2[t] - k3 (1-y2[t]) + k4 y3[t]
->             y3'[t] = k3 (1-y1[t]-y2[t]-y3[t]) +k2 y3[t]
->
->I need to fit my data to get those k1, k2, k3, and k4. What is the
easiest
->
->
->SF

What I don't know about ODEs fills volumes (and in fact did - it filled my
textbooks from my undergrad DE courses :-), but since I don't see any
responses here, I'll try giving you a partial answer.

I assume you want to estimate k1 ... k4 to give the minimum
sum-of-squared-errors fit to the data.  For starters, we can put the
equations and functions in lists and define a function to solve the ODEs:

eqs = {y1'[t] == k1 y1[t] + k2 (1-y1[t]-y2[t]) + k4 y3[t],
y2'[t] == k2 y2[t] - k3 (1-y2[t]) + k4 y3[t],
y3'[t] == k3 (1-y1[t]-y2[t]-y3[t]) +k2 y3[t]};

funcs = {y1[t], y2[t], y3[t]};

f[ k_List ] :=
Flatten[
funcs /.
DSolve[ eqs /. Thread[ {k1, k2, k3, k4} -> N[ k ] ],
funcs, t ]
] /;
Length[ k ] == 4

Function f takes a list of four numbers as input, solves the diffeqs, and
spits up a solution.  I've checked to make sure the list has four entries,
but I haven't put in a check to make sure they're numbers (easy to add).
Also, I took the liberty of forcing them to machine precision (N[ k ]) to
avoid getting some seriously lengthy output as Mma tries to use "infinite
precision" where possible.

A typical application of f[] would look like this:

f[ {1, 2, -1, 3} ] /. t -> {-1., 1.}

1
Power::infy: Infinite expression -- encountered.
0.
<repeated a few times>
General::stop:
Further output of Power::infy
will be suppressed during this calculation.

{{2. - 0.00238455 C - 8.26016 C + 0.26013 C,
2. - 9.83365   C - 0.10936 C + 1.92212 C},
{-0.333333 - 0.010591 C -   1.85976 C -
0.26013 C, -0.333333 - 43.6763 C -
0.0246222 C - 1.92212 C},
{-0.222222 - 0.0111639 C + 1.96036 C -
-17
1.14325 10    C,
-0.222222 - 46.0389 C + 0.025954 C -
-17
8.44757 10    C}}

Dimensions[%]
{3, 2}

In this example, I've fed f[] estimates of the four parameters and
substituted into the ODE solution two choices of t (-1, +1).  The result is
a 3x2 matrix in which each row is the vector {y1[t], y2[t], y3[t]} for a
different value of t.

The basic game plan is to write a function that takes as input the list k
and spits back the sum of squared differences of observed vs. estimated y
values.  You would be summing over all time epochs and across the three
functions y1, y2, y3 to get a single total squared error.  That's the
function you minimize (with FindMinimum).  That function will use the
function f[] above to get teh "estimated" y values.  You'll need to make
those constants of integration C, C, C go away, either by
specifying some more equations (initial conditions?) or by adding them to
the parameter list (so that instead of solving for {k1 ... k4}, FindMinimum
solves for {k1, ..., k4, C, C, C}).  And you'll need either a fast
machine or some serious patience.

I'm not sure, but this sounds like it should work.

-- Paul

**************************************************************************
* Paul A. Rubin                                  Phone: (517) 432-3509   *
* Department of Management                       Fax:   (517) 432-1111   *
* Eli Broad Graduate School of Management        Net:   RUBIN at MSU.EDU    *
* Michigan State University                                              *
* East Lansing, MI  48824-1122  (USA)                                    *
**************************************************************************
Mathematicians are like Frenchmen:  whenever you say something to them,
they translate it into their own language, and at once it is something
entirely different.                                    J. W. v. GOETHE

```

• Prev by Date: Re: Fit in Excel Mathlink
• Next by Date: Re: Are two expressions equal?
• Previous by thread: Nonlinear fit to diff eqs.
• Next by thread: Gram matrix