Barnes D.J., Chu D. Introduction to Modeling for Biosciences

Подождите немного. Документ загружается.

Chapter 5

Mathematical Tools

Nowadays computer algebra systems (CAS) are an essential tool for the mathemat-

ical modeler. In a strict sense, a CAS is a computer program that can perform sym-

bolic computation, that is manipulate mathematical expressions. Such programs can

be an important aid to modelers of all skill levels. Apart from CAS as tools for per-

forming symbolic mathematics, we will also include a second type of computational

mathematics tool into the deﬁnition of CAS: Programs that can aid the modeler in

producing numerical solutions to mathematical models.

We start by giving a brief overview of the range of available software, and then

the bulk of the chapter will introduce one particular system—Maxima—in some

detail. We include some case studies that explicitly work through models encoun-

tered previously in this book. At the end of the chapter, the reader should be able to

formulate and solve simple models in Maxima.

5.1 A Word of Warning: Pitfalls of CAS

Hours and days and weeks of the life of a mathematical modeler are spent looking at

formulas, deriving expressions and proving theorems. Many of us enjoy this activity

and feel a sense of satisfaction once the task is accomplished. Yet, particularly when

the mathematics becomes more involved and the formulas longer, it is often helpful

to hand over parts of the task to a CAS. Unlike us, these systems can perform the

chores of mathematical modeling in an instant, without complaint and, most signif-

icant of all, without the typical errors of negligence that humans are liable to make

when manipulating strings of Greek letters and strange symbols. However, while ex-

tremely useful, one must always keep a watchful eye to avoid the potential pitfalls

of handing over large portions of one’s own work to a machine.

Machines lack the innate laziness we have and do not mind producing lines and

lines of results when the same statement can be expressed much more concisely.

While most computer algebra systems have a simpliﬁcation command to reduce the

size of expressions, in practice this sometimes fails. Blind trust in the capabilities of

the machine then results in long incomprehensible equations, which only serve to

D.J. Barnes, D. Chu, Introduction to Modeling for Biosciences,

DOI 10.1007/978-1-84996-326-8_5, © Springer-Verlag London Limited 2010

183

184 5 Mathematical Tools

obscure insight rather than providing it. For this reason, it is essential to check the

results of computer generated computations by hand, wherever possible.

In the practice of mathematical modeling, it is often the very process of wrestling

with the equations that leads to insight. The ﬁnal result is then only the distilled

version of weeks of struggle with the model. It is the secret shortcut that can only

be discovered by walking backwards from the destination, but could not have been

discovered from the starting point. If we hand over the struggle to the machine

then we risk losing or never gaining some of that insight and perhaps miss out on

the crucial formula that will then be published by our academic arch-enemy the

following year (in Nature, even). It seems that the arch enemy is still computing by

hand.

Undoubtedly, computer algebra systems are useful, and even more so for non-

mathematicians. Built into these systems are libraries of mathematical knowledge

that few non-mathematicians possess (and indeed maybe few mathematicians). Of-

ten a quick query of a CAS will deliver an instant answer that would otherwise re-

quire weeks of pestering mathematically more versed colleagues for help; or, even

more frustrating, hour upon hour spent in the library on fruitless self-study. There is

signiﬁcant value in computer algebra systems, but only as long as the modeler re-

mains the master and does not succumb to a naive trust in the power of the system.

While powerful, CAS still need to be systematically doubted and their results

checked—not because the results are likely to be wrong but because they often

are presented in a form that is not useful or hides what we are seeking. Whenever

something can be calculated by hand then it should be done by hand; the computer

then serves as a useful tool in conﬁrming the result.

These cautions mainly apply to computer algebra systems whose primary pur-

pose is to assist the user in symbolic manipulation of equations. When it comes to

using computers to generate numerical solutions the situation is different. They are

more of a necessity. Sooner or later, most practical modeling projects need to per-

form some sort of numerical computation, if only to plot an example curve for a

paper to be published. Many mathematical models are so complicated that symbolic

solutions cannot be obtained. In these cases there is an intense need to solve equa-

tions numerically. Computer systems that have the available algorithms to do this

are then an indispensable tool for modelers.

Designing efﬁcient algorithms to perform numerical computations is a ﬁeld of

research in its own right. In rare cases a modeler may need to design such a algo-

rithms herself but, in the vast majority of cases, it would be a waste of time not to

use the ready made, off-the shelf systems that are available. CAS are a library of

implementations of such algorithms and can save the user months of preliminary

studies. It is worth bearing in mind that a great deal of sophistication has gone into

the design of these algorithms and the particular implementations that underly the

computer programs are complex pieces of software engineering. Still, we should

always be on guard and critically examine the sanity of the results we obtain. Like

all software, some implementations will contain inefﬁciencies, bugs and errors that

may lead to wrong results. If in doubt and, in the case of very complex models as a

matter of routine, it is also worthwhile cross-checking the obtained results by using

5.2 Existing Tools and Types of Systems 185

two independent pieces of software. This is time-consuming, but an essential insur-

ance against falling prey to errors in the software, artifacts produced by a quirky

algorithm, but also our own mistakes in formulating the model.

This section will introduce the reader to a freely available computer algebra

system and, in addition, provide a general introduction, a few practical examples,

mostly of problems that have been encountered before in the book.

5.2 Existing Tools and Types of Systems

There are a number of systems available that can be useful for performing math-

ematics on a computer. The three dominant commercial systems are M

APLE,

MATHE MATICA,

and MATLAB.

The main strengths of MAPLE and MATHEM AT-

ICA are their powers of performing symbolic computations. In addition, both have

evolved from bare-bones computer algebra systems to programs that also provide

a visually appealing interface capable of producing potentially publishable docu-

ments.

The main strengths of M

ATLAB is its numerical capability. It is designed to facil-

itate rapid proto-typing of complex numerical modeling and programming. It comes

with a number of toolboxes that allow one to build on the work of previous users.

In many areas of science, including biology, M

ATLAB has evolved to a de facto

industry standard.

A potential drawback common to all three systems is their high cost. Universities

and companies often have site licenses that cover various types of use, and may even

allow installation of these programs on private machines. Student licenses may also

be available relatively cheap. At the full price, these systems are very expensive

which motivates the search for cheaper alternatives . . . and these exist.

A well known system that allows much of the capabilities of M

ATLAB is Octave.

There is a near cross-compatibility between Octave and MATLAB, in the sense that

most code for Matlab will work on Octave and vice versa. At the time of writ-

ing, however, Octave does not have a graphical user-interface like M

ATLAB, which

might deter some from using it. On the plus side, Octave is freely available. For

readers familiar with Matlab who would like to undertake primarily numerical com-

putation, Octave is certainly a possible alternative to its commercial counterpart.

Another freely available system that is closely modeled on Matlab, but does provide

a graphical user interface, is Scilab.

For heavy duty mathematical modeling the

ASCEND software

is also worth a look. It is also free but, in our experience, not

as user-friendly as Octave and Scilab.

http://www.maplesoft.com.

http://www.wolfram.com/products/mathematica.

http://www.mathworks.com/.

http://www.gnu.org/software/octave.

http://www.scilab.org.

http://ascend.cheme.cmu.edu.

186 5 Mathematical Tools

A system that is useful for some numerical tasks is R.

It is primarily intended

for statistical computing and has a number of statistical functions built into it. Like

Octave, it is freely available, can be used completely from the command line and

supports a simple scripting language. In addition, there are also freely available

graphical user-interfaces available for it. Unlike Octave, it is mainly intended as a

tool to manipulate and analyze data, but has limited use in computationally more

intense mathematical computing, such as differential equations. In our view it also

suffers from a rather idiosyncratic syntax which certainly makes sense, but not to

the uninitiated.

A high quality freely available system that allows a quick and relatively reliable

exploration of systems of differential equations is XPPAUT.

It has no analytical

capabilities and is speciﬁcally meant to explore systems of differential equations. It

can do nothing that cannot be done with other more generalist systems like Maple

and Octave. Its main value, in our opinion, lies in the efﬁcient way in which it allows

the user to explore systems of differential equations. XPPAUT is designed to enable

the user to understand rapidly how varying initial conditions impact the result. It

can quickly display trajectories in phase space and even plot bifurcation diagrams.

The downside is its idiosyncratic user-interface. It may take the uninitiated some

time to ﬁnd her way around the system. Initially, this may be frustrating, but it

is well worth persevering. Should the reader anticipate using differential equation

models then we would urge her to set aside an afternoon investigating XPPAUT.

We are conﬁdent that the initial investment of patience will pay off in dividends of

improved efﬁciency.

There are also a few freely available software packages that allow symbolic

computation. A well known system is Sage

which can be used via a web-based

interface. Another system worth mentioning is Mathomatic

which can perform

symbolic computations. It runs on very small systems as well, including handheld

devices. One of its particularly nice features is that it can be accessed via the tel-

net program, available on most operating systems. It needs no installation or fancy

graphical user interface. Here is a little demonstration.

Assuming connection to the Internet, the mathomatic server for interactive ses-

sions can be reached by typing into a shell the following command:

telnet mathomatic.org 63011

The reader can then have a session such as the following to compute the expres-

sion



−bx

−a



and to evaluate it for x =1, b =2, a =3. (Text after -> represents user input.)

http://www.r-project.org.

http://www.math.pitt.edu/~bard/xpp/xpp.html.

http://www.sagemath.org.

http://mathomatic.org.

5.3 Maxima: Preliminaries 187

> telnet mathomatic.org 63011

Trying 213.133.98.206...

Connected to mathomatic.org.

...

1-> x^b - b*x^(-a)

#1: (x^b) - (b*(x^(-1*a)))

1-> differentiate

Enter variable: x

Differentiating with respect to (x) and simplifying

b*((x^b) + -----)

(x^a)

#2: -----------------

2-> calculate

Enter x: 1

Enter b: 2

Enter a: 3

2-> exit

We will leave it at this one small example but encourage the interested reader to

explore this system in more detail. The remote availability of Mathomatic is a po-

tentially attractive feature.

5.3 Maxima: Preliminaries

In this book we will introduce the reader to another computer algebra system, Max-

ima.

It is free and can perform both symbolic and some numerical mathematics. It

can be compiled for the most common operating systems. Unlike Maple or Matlab,

Maxima does not come with its own graphical user interface, although there is a

choice of interfaces within which Maxima can be used. The simplest approach is

to run Maxima in a command-line environment. This does not provide any fancy

interface, but does the job with a minimum of hassle. A simple, yet useful graphical

user interface for Maxima is xMaxima which comes with most distributions. The

http://maxima.sourceforge.net.

188 5 Mathematical Tools

reader may wish to check her installation for xMaxima. Additionally there is TeX-

macs

as a possible interface. TeXmacs is also a document preparation system in

its own right. A somewhat simpler alternative interface is wxMaxima,

which is

purely a graphical user interface for Maxima. The choice between these is left to the

reader. Each has its own merits and each is certainly very spartan compared to what

one has become used to from commercial interfaces. While wxMaxima has a con-

venient, well thought out graphical interface, producing, manipulating and saving

plots is perhaps easiest in xMaxima. The best mathematical output is provided by

TeXmacs, although its general purpose interface is not as user-friendly as the other

two. In the remainder of this chapter we will remain agnostic about the interface to

be used.

To reap maximum beneﬁt from the rest of the section, we strongly urge the reader

to choose one installation of Maxima, to actually install it and familiarize herself

with the system as we progress through the examples to follow.

Before starting the description of how to use Maxima, it is worthwhile providing

a general overview of the software and how it is organized. Unlike the very polished

commercial products that the reader may be acquainted with, Maxima is the product

of an academic open source project. Its main focus is on substance and quality. Also,

it is born out of the more spartan world of UNIX where things are, to a large degree,

text-based and centered on the use of a command line interface. This might sound

like an antiquated approach to software but, in reality, it has its advantages, partic-

ularly once one is used to the concept and has overcome any fear of command-line

windows. In fact, the minimalism of simple interfaces has the advantage of focus-

ing the user on the task at hand, rather than distracting with unnecessary dialogues,

bells and whistles. More importantly, programs of this sort tend to be more stable

and faster than their more visually-appealing cousins.

An apparent drawback is that it appears to be more contrived to extract graphical

output from these programs. Maxima itself does not have its own graphical routines,

yet it is still possible to use it to produce plots and to use the plots in high-quality

documents such as scientiﬁc papers. To achieve this Maxima uses a third-party pro-

gram (also free) called Gnuplot.

This is a very powerful program in its own right

(and has been used to prepare many of the graphs in this book). In practice, when-

ever a curve is plotted from Maxima, the program produces a set of points that are

then (automatically) passed to Gnuplot, which in turn produces the plot window that

is shown. Gnuplot can also export images in different ﬁle formats. Depending on the

installation, this can include PDF, JPEG, GIF, and many others. In what follows, we

will always assume that the ﬁle output format is Postscript, which is the simplest for-

mat to use and is always available with Gnuplot. On MacOS X and Linux, display

of Postscript ﬁles is unproblematic. On other platforms this may be more difﬁcult

so, as an alternative, the reader may like to choose PDF as the output format. The

Maxima documentation provides details.

http://www.texmacs.org.

http://wxmaxima.sourceforge.net.

http://www.gnuplot.info.

5.4 Maxima: Simple Sample Sessions 189

Unlike many of the programs the reader may be familiar with, Gnuplot does not

allow the saving of graphs by direct mouse-interaction. Instead, commands must

be issued within the Maxima environment to write a graph to a ﬁle that can subse-

quently be integrated into a document, say. This may sound complicated but it turns

out to be a very effective process for producing high quality graphs.

The xMaxima program provides a viable alternative for those readers who do

not wish to embark on this adventure into command-line computing. The xMax-

ima interface allows direct plotting within the interface, using openmath rather than

Gnuplot. What is more, every graph comes with an interactive menu that lets the

user manipulate the graph and save it to ﬁle.

5.4 Maxima: Simple Sample Sessions

The best way to quickly get to work with any system is to follow through a num-

ber of examples. In this section we will present some sample sessions designed to

illustrate the main features of the Maxima computer algebra system. The examples

will be sufﬁcient to reproduce the mathematical models presented in this book, and

to get going with some initial mathematical modeling. However, in the long run,

the serious user will need to dig deeper herself. In our opinion, the best approach

for this is a “learn-as-you-go” approach, whereby the user only familiarizes herself

with new concepts as and when she actually requires them.

5.4.1 The Basics

The most efﬁcient way to acquire the basic skills necessary to become at least a

novice user of Maxima is to follow the examples given here, modify them, and thus

develop an intuition for the program and how to use it. The results of any compu-

tations should be reproducible independently of the particular version of Maxima

used. However, the formatting of the output may vary from interface to interface,

and may not be as visually appealing as the examples given here, which were pro-

duced using TeXmacs.

A few key principles of Maxima ﬁrst:

• The operators for multiplication, division, addition, subtraction and exponentia-

tion are: *,/,+,-,ˆ, respectively.

• The end of every statement and expression must be indicated with either a semi-

colon or with a dollar sign. If the former is used then the input will be echoed back

to the user; the dollar sign suppresses any output. This can be desirable when the

output is expected to be long. A common error of novice users is to forget to end

the line with a semicolon or a dollar sign.

• Equations, variables and statements can be named. This is done by preceding the

expression with the name and a colon sign (e.g., name : expression).

190 5 Mathematical Tools

• Every input and output is assigned a line number of the form %ix or %ox,for

input and output respectively. In this case “x” will be a number. These numbers

can be used to recall previous statements.

What follows is a full Maxima session that illustrates its basic use. A look at

these example interactions demonstrates that input of expressions and equations is

very easy.

Maxima 5.18.1

http://maxima.sourceforge.net

Using Lisp SBCL 1.0.29

...

(%i1) x;

(%o1) x

(%i2) x$

(%i3) eq1: x*x*x - x^3 + y - x^4;

(%o3) y −x

(%i4) eq1;

(%o4) y −x

(%i5) %o3

(%o5) y −x

(%i6) eq2: a*x - z=0;

(%o6) ax −z =0

(%i7) a: 3$

(%i8) eq2;

(%o8) ax −z =0

(%i9) ’’eq2;

(%o9) 3x −z =0

(%i10) eq3: a*x^a = 3;

(%o10) 3x

(%i11) quit();

The end

Most of the above is self-explanatory, however a short comment is necessary

regarding the line marked %i3. The user will notice that the system has performed

a certain amount of processing between the input and the output. Speciﬁcally, it has

recognized that x ∗x ∗x is the same as x

. Often output from Maxima needs to be

simpliﬁed by issuing a speciﬁc command. We will discuss this in more detail later

on.

The reader should also pay particular attention to input line %i7, where we assign

the value 3 to the symbol a. Recalling in %i8 eq2—an equation that contains

symbol a—we ﬁnd that the symbol has not been replaced by its numerical value.

We can recover eq2 with the speciﬁc value of a by preceding it with double quotes

(%i9). However, if we deﬁne a new equation containing a, in this case eq3, then

5.4 Maxima: Simple Sample Sessions 191

we ﬁnd that it has been replaced by the numerical value. Maxima sessions are ended

by issuing the quit() command.

One often wants to write a formula with parameters that are initially unspeciﬁed,

and only later commit to speciﬁc values for these variables—for instance, to gener-

ate examples for plotting. To be able to do this, Maxima provides a way for the user

to assign values to variables.

(%i1) eq1: x^2 - u;

(%o1) x

−u

(%i2) eq2: ev(eq1,u=2);

(%o2) x

−2

(%i3) eq3: ev(eq2, x= sin(y));

(%o3) sin(y)

−2

(%i4) eq4: ev(eq3,y=1);

(%o4) sin(1)

−2

(%i5) eq5: ev(eq4,numer);

(%o5) −1.291926581726429

(%i6) eq6: ev(eq1,u=2,x=sin(y));

(%o6) sin(y)

−2

(%i7) eq7: ev(eq1,u=2,x=sin(y),y=2,numer);

(%o7) sin

(

)

−2

The command ev(...) is used both to replace parameters or symbols by num-

bers and to replace them by other expressions. Line %i3 illustrates this latter capa-

bility, where we replace x by the expression sin(y). Since Maxima is primarily a

computer algebra system (rather than a fancy calculator) it tends to prefer to keep

things as general and precise as possible. As a result, in %i7, even when we have

speciﬁed a value for y, the system still returns the “symbolic” form of the equa-

tion. Yet, sometimes we want to know the actual number, rather than just being told

that the result is sin(1) −2. In this case the command ev is used with the keyword

numer to force evaluation of the expression to a numerical value.

A particularly useful feature of ev is that it allows a sequence of replacements

to be speciﬁed, simply by stringing them together in a single expression, as in %i6.

Note, however, that due to the way ev is processed internally by Maxima, the value

2 is not actually assigned permanently to y in %i6, whereas the variable u is re-

placed. The difference is that y is substituted into the expression, whereas u is not.

Naively, one would have expected that output lines 5 and 7 produce the same result,

but they do not. In order to rectify this, one needs to tell ev to add an additional

evaluation step. This can be done by adding the keyword eval after the expression

to be evaluated, as so:

(%i1) eq1: x^2 - u$

(%i2) ev(eq1,eval,u=2,x=sin(y),y=1,numer);

(%o2) −1.291926581726429

In the remaining example sessions we will not explicitly write the quit() command any longer.

192 5 Mathematical Tools

We have seen that Maxima allows naming of expressions so that they can be

recalled later on (for example eq1 in the last session). Often we need to deﬁne

customized functions. Not only does Maxima allow us to deﬁne our own simple

functions, such as f(x)=x

−y, but also more complicated compound functions,

which may perform very complex evaluations. The following session illustrates how

a simple function can be deﬁned.

(%i1) f(x,y):=x^2 - sqrt(y);

(%o1) f

(

x, y

)

:=x

−

√

(%i2) f(d,f);

(%o2) d

−

√

(%i3) f(1,2);

(%o3) 1 −

√

This principle immediately generalizes. So, for example, the function f(x) =

can be deﬁned by writing f(x) := xˆ2; note that “:=” is used for function

deﬁnitions and not the simple colon assignment operator. Indeed, entire procedures

can also be deﬁned with this notation by enclosing them into a block environment.

(%i1) f(x,y):= block(

fv1: x^2 - 6,

fv2: y^5 - 7,

fv1/fv2);

(%o1) f(x, y) := block



fv 1 :x

−6, fv2 :y

−7,

fv 1

fv 2



(%i2) f(a,b);

(%o2)

−6

−7

Within a block environment we can specify any sequence of steps to be evaluated

one after the other. The last statement is returned as output. Note that statements

are separated by commas rather than by dollar signs or semi-colons. (The special

formatting we used in input line 1 is not essential, but chosen for readability. Blocks

can be deﬁned on a single line as well.)

Variables and names in functions are globally valid. In the next example, pay

attention to the sequence of assignments.

(%i1) x: 1$

(%i2) u: s^2 - o$

(%i3) f(x,y):= block(

u: x^2 - 1,

v: y,

sqrt(u + v))

(%o3) f(x, y) := block



u :x

−1, v :y,

√

u +v



(%i4) u;

(%o4) s

−o

(%i5) f(1,2);

(%o5)

√

(%i6) f(2,2);

(%o6)

√