Download Pharmacology Model Building Using Berkeley Madonna FlowChart Interface and more Slides Pharmacology in PDF only on Docsity!
MCB 137 WINTER 2008
PHARMACOLOGY
Introduction
This chapter is a short tutorial on model building using the Berkeley Madonna FlowChart interface.
We do this within the context of pharmacology because of the abundance of pharmacology models
that are both simple and practical. Pharmacology models are designed to predict the time
dependence of a drug's concentration in the body fluids following its administration, or to model the
action of the drug once it reaches its target organ. The complexity of pharmacological models varies
widely depending on the level of detail they encompass. We begin by representing the whole body
as a single c ompartment ; that is, a homogeneous volume. This simple model will provide us with an
easy introduction to modeling with Madonna, and suffices to introduce the general concepts of
steady state , exponential response , relaxation time , and the principle of superposition in linear
systems.
Although the single compartment model yields useful results in many cases, it is not a very realistic
representation of the body, especially if we are interested in the distribution of a drug into different
tissues. So the next step up in complexity is to model a body as two interconnected compartments: a
central compartment that includes the blood and tissues like the liver and kidney, that have a rapid
and profuse blood circulation, and a second peripheral compartment that includes the more slowly
perfused tissues. While more successful than a one compartment model, the abstract nature of the
two compartment system makes it difficult to assign anatomical locations to the compartments, and
to give physical interpretations to its parameters. These issues prompt us to move on to more precise
and more elaborate physiologically based pharmacokinetic models which take explicit account of the
blood supply to each represented tissue. Finally we conclude this chapter with an introduction to the
more complex and more challenging problem of the action of drugs at their target site.
One Compartment Pharmacokinetics
We begin with the simplest case where we infuse into the body the biologically inert polysaccharide,
inulin. Inulin is infused into animals or humans to trace fluids; it is used to estimate the volume of
extracellular fluids as well as the rate of filtration of blood plasma fluids into the kidney nephrons.
This filtration rate, a fundamental parameter for investigation and evaluation of kidney function, is
calculated from measurements of urinary flow and urinary and plasma concentrations. Interpretation
of these measurements requires some notion of how the inulin concentration in the plasma changes
following its administration by injection or infusion. Following the start of an inulin infusion, the
plasma concentration rises rapidly. When will the concentration stabilize?
Model Summary
The model computes the inulin concentration in the plasma at any time. It does this by comparing
the amounts entering the body with the amount leaving at each moment; any discrepancy must
represent a gain or loss of inulin within the body. Adding all these gains (or losses) up to a particular
time computes the amount and concentration of body inulin at that time.
Before starting, verify that the flowchart editor
is functioning properly. To do this, choose
About Berkeley Madonna from the Help
menu and look for text that reads “Flowchart
version 8.0.1”. If it says “Java not loaded”,
click the Load Java button. If Java fails to
load, you need to install (or reinstall) Java™
support on your system.
Figure 1. Changes in the inulin body compartment depend on the entry and exit flow rates.
The results enable us to study and evaluate protocols for making measurements of extracellular
volume and of renal filtration.
Model For The Accumulation Of Inulin In The Body Fluids
Building the Model --- Step by Step
In this section, we construct our inulin model
graphically using the flowchart editor. The
model consists of three dynamic elements: a
reservoir icon (compartment) representing
the mass [mg] of inulin in the body fluids,
and two flow icons representing inulin
entering and leaving the body [mg hr
]. In
addition we will employ supporting formula icons to represent the rate of rate of renal filtration [L
hr
], volume of body fluids L, and the concentration of inulin [mg L
]. Begin:
- Choose New Flowchart from the File menu. An empty flowchart window appears.
- Position the mouse over the reservoir tool , press the mouse button, drag the reservoir to
the flowchart, and release the mouse button. You’ve now placed a reservoir on the flowchart
with the default name R1. Note that the reservoir is colored red which means it is selected.
- Change the name of the reservoir to m by typing the new name. (It is not necessary to first click
the reservoir since it is already selected.) This single reservoir represents the total mass in
milligrams [mg] of inulin which we assume is uniformly dissolved in the extracellular fluid.
- Click the flow tool , move the mouse over the mass reservoir, press the mouse button, drag
the mouse to the right a few inches, and release the mouse button. This places a flow icon on the
flowchart which drains the mass reservoir into an infinite sink.
- Make sure the reservoir is connected to the flow. If the flow looks like this:
Figure 2. Inulin Model before
setting the model
parameter values.
- So far we have not supplied the model with any numerical values. Madonna indicates this by the
question marks on the icons. Double-clicking on any icon will open a dialog box asking for
numerical data. For example, double-click on the C icon. The dialog box says that m and V are
required inputs. This means that these variables must be used in the right-hand side of the
equation defining the icon’s value. Enter m/V in the data box. You can type in the names from
the keyboard or, more conveniently, simply double-click on the input’s name in the list. Click
OK (or press Return) when done, and the dialog box will close; the question mark on the icon
has disappeared indicating that Madonna is satisfied with your input to the flow icon. Next we
must assign numerical values to V and k. The basis for these values is discussed in the
Parameter Estimation box.
Parameter Estimation
The molecular weight of inulin is ~ 5000, so it is too large to pass through cell membrane channels,
but it is small enough to easily pass through porous capillaries. There are no other significant
transport pathways for inulin. Therefore, to a good approximation, we can assume inulin is
uniformly dissolved in the extracellular fluid that makes up about 20% of the body weight. The
extracellular volume for a 70 kg body is 14 liters [L] (20% body weight = 0.20 × 70 kg = 14 kg ≈ 14
L).
We calculate the rate constant k for excretion from the mass action relation Ex = k c C. Since inulin is
biologically inert, its only escape from the body fluids is by renal excretion so Ex is easy to measure:
simply collect the urine that flows over a given period of time, measure its inulin content and divide
by the collection time. When C = 100 mg L
, typical values of Ex are around 720 mg hr
these figures, we estimate k c = Ex/C = 720/100 = 7.2 L hr
- Double-click the k c icon and enter 7.2. Note that you can use the dialog’s calculator keypad to
insert digits and common arithmetic operators. Double click the V icon and enter 14.
- Double click the reservoir ( m icon). In this dialog, you specify the reservoir’s initial value. Since
there is normally no inulin in the extracellular fluids, we assume there is none at the beginning of
the experiment, so enter 0 and close the box.
- Now double click on the flow icon (click in the center circle with the question mark). k c and C
are required, so enter their product, k c
C, and click OK. (Flows are rates and have time in their
dimensions. The dimensions of k c
= [L⋅hr
] while C = [mg⋅L
] so that the dimensions of the
product k c C = mg⋅hr
, the rate that the reservoir is losing inulin.
- Double click on the last unfilled icon, Inj. We will simulate a constant infusion at the rate 600
mg⋅hr
- 1 . Enter 600 in the dialog and click OK
- The model flowchart is now complete, and all question marks should have disappeared. Now run
your model by choosing Run from the Compute menu. Berkeley Madonna runs your model and
displays the results in a graph window. Click on the button labeled L located on top of the graph
to insert a Legend on the graph. Click the button again to remove the legend. All buttons on the
graph are toggle (on-off) buttons.
Figure 3. The graph window.
- The graph shows the variables m , whose scale is on the left vertical axis, and Ex scaled on the
right. (Madonna automatically plots the first two variables that were defined in the construction
of the model window The first is scaled on the left axis, the second on the right.) There are 3
toggle buttons located on the left bottom margin of the graph window. Click on the C button.
This will place a plot of C on the graph. Click on the Ex button to remove Ex. Both m and C are
scaled on the same left hand axis. To scale one of them, say C , on the right, click on the C button
while holding down the shift key.
Figure 5. Equation Window obtained from the menu
Model>Equations.
The model predicts an exponential response
Many experiments begin with a system in a steady state where all flows are constant and the
reservoir levels are unchanging. The system is suddenly perturbed, for example by a sudden
concentration change, after which the system moves to a new steady state. We can use the model to
determine the new steady state and how long it took to get there.
The steady state
Whenever the flow into the body through injection is equal to the flow out through excretion the
system is in a steady state condition; that is, the flow, Inj , into the reservoir m is equal the flow out
Ex = k c
C. Then the steady state concentration of inulin will be given by
Equation 1 steady state: C = Inj / k c
[mg/L]
When you select an icon in the flowchart, Berkeley Madonna highlights the corresponding
equation(s) in the Equation window. This makes it easy to see an icon’s equations without
opening its icon dialog. Conversely, you can show which icon corresponds to an equation in
the equation window. Simply position the equation window’s insertion point (caret) within
one of the equations, then choose Show Icon from the Edit menu. This activates the
flowchart window and highlights the corresponding icon.
You can also open the icon dialog directly from the equation window by double-clicking an
equation. For example, double-clicking anywhere within the line C =m/V opens the icon
dialog for the C icon. Select the m reservoir, type Inulin and press return. The name of the
icon changes and Berkeley Madonna updates the name in all equations that depend on it. To
see this in action, change the name again and keep your eye on the equation window as you
press the return key. Change the name back to m
When Inj = 0, C = 0; this is our initial condition (a trivial result simply restating what we already
knew: inulin is not a normal body constituent). When Inj = 600 mg/hr, the final steady state
concentration is C = 600/7.2 = 83.3 mg L
- 1 . You can verify this on the graph by ‘zooming in’ to view
the last portion of your plot: click and drag the mouse to draw a rectangle over the relevant area.
When you release the mouse, Berkeley Madonna adjusts the axis limits so that the portion that was
in the rectangle now fills the window. To undo the effect of a zoom-in, click the Z button located in
the toolbar on top of the plot. Now click the Readout button, on the Graph window toolbar
and place, the ‘crosshair’ cursor over a point on the plot; the coordinates are displayed
in the information area at the top of the window. Alternatively, you could view the data in
Table form and scroll down to the last time row.
More complex models will have more reservoirs, but the same considerations hold: in a steady state
the inflow to each reservoir must equal its outflow. This can involve finding roots of a complicated
set of simultaneous algebraic equations, but Madonna has provisions for finding the required
numerical solutions.
Relaxation Time
Letting τ = V/k c , the mathematical solution for the inulin model as derived in the MATH NOTES box:
The exponential response is:
Equation 2
C ( t ) =
Inj
k c
1 − e
−
t
τ
C(t) approaches asymptotically (but never really reaches) the time independent steady state. So it
makes no sense to ask how long it takes the system to reach its final destination. Instead, we ask how
long it takes for C to reach a prescribed fraction of its final value, say 50%. In general, e
will
reach the 50% mark when t ≈ 0.69⋅τ ( e
≈ 1/2). The half-time, t 1/ , is given by
Equation 3 t 1/
= 0.69⋅τ
Rather than carry around the numerical factor 0.69, we use τ itself as a measure of the time scale of
an exponential process. τ is called the relaxation time or sometimes the time constant ; it is the time
required for the exponential to reach 63% of its final value.
In addition to the notational convenience, relaxation times are preferred because they have a physical
interpretation: they measure the mean residence time of a molecule in the reservoir. In our inulin
example each injected inulin molecule has a different fate. After an injection into the circulation
some molecules travel to the kidneys and are excreted almost immediately. Others escape into tissue
spaces and may meander around for long periods before entering the kidney. These molecules may
have ‘life times’ of many hours. The mean residence time is the average lifetime of all the injected
molecules, and this turns out to be equal to the relaxation time. This is derived in the MATH NOTES
box: The mean residence time. Similar conclusions hold for any exponential process; for example,
the relaxation time ( = t 1/ /0.69) for radioactive decay equals the average time that molecules persist
before undergoing a radioactive disintegration.
Note that the relaxation times for the uptake and decay processes are identical
MATH NOTES: MEAN RESIDENCE TIME
Assume that at any arbitrary time we set t =0 and label the inulin molecules present in the body. Let
there be N 0
of these molecules and let us note the time each labeled molecule leaves the body. If n i
the number of labeled molecules that leave at time t i ,
then the average time a molecule spends in the
body is given by
t mean
n 1
t 1
t 2
N
0
n i
t i ∑
N
0
We can calculate these times from our results. We begin with N o
molecules at t=0. Since the
number of molecules at any time N is proportional to m, we have
N = N
0
e
−
t
τ
and the number that leave in the time interval between t i and t i +dt is given by – dN i
, (The minus sign
occurs because we are counting the number of molecules that leave the body for the environment -
each time the environment gains one, the body loses one)
− dN i
dN i
dt
dt =
τ
N
0
e
−
t i
τ
Using this result together with our expression for the mean, we have
t mean
N
o
− dN i
( ) t i
i = 1
∞
∑
τ
e
−
t i
τ
t i
dt
i = 1
∞
∑
As dt gets smaller, the sum approaches the following integral that we evaluate using standard
integral tables:
t mean
t
τ
e
−
t
τ
0
∞
∫
dt = τ
Note that this result ( t mean
= τ ) was derived entirely from the rate of disappearance of molecules. As
long as this rate is independent of incoming molecules (e.g. when Ex is proportional to C ), the result
is also independent of incoming molecules.
Dosage regimens
In the inulin model the inulin plasma concentration is proportional to the amount of drug in the body
( C =m/V ) and its excretion rate is proportional to the plasma concentration ( Ex = k ⋅ C ) But unlike
inulin, drugs are not inert. They are absorbed on plasma proteins, dissolve in fatty tissue, and enter
cells where they may be metabolized. To model the fate of an injected drug we can use the same
approach to compute C and Ex, but t he constants of proportionality are now empirical quantities
that may not have an explicit physical interpretation. We write formally as before
Equation 4 C =
m
V
d
Equation 5 E x
= k c
C
where V d and k c are experimentally determined parameters. V d is called the volume of distribution,
and k c is called the clearance. Loosely speaking, V d is a volume ‘corrected’ for absorption and
solubility in different tissues, and k c
accounts for metabolic destruction as well as excretion (see
Empirically defined parameters box) These two parameters are very useful because their values
for virtually any common drug are listed in tables found in textbooks and handbooks of
pharmacology
Unlike inulin, drugs have biological effects, some desirable, some not. A common problem is to
determine a dosage regimen (schedule) to produce plasma concentrations that are high enough to
Empirically defined parameters
The use of an empirical constant to correct for absorption and solubility is valid only over the range
of drug concentrations where absorption and solubility are linear functions of concentration. For
example, if the amount of drug absorbed = a ⋅ C and the amount dissolved in fatty tissue = bC,
where a and b are proportionality constants. Then the total amount in the body = m = VC +aC +bC
= (V+a+b)C. Solving for C leaves
C =
m
V + a + b
m
V
d
so that V d = V +a + b.
Similarly, if drug removal by metabolism and excretion are both linear then
Ex = k metab
C + k excr
C = k metab
( ) C = k c
C
where k metab
and k excr
are constants. Thus, the clearance, k c
= k metab
Where q is the total amount that is delivered, first is the time of the first pulse, and interval is the
interval between pulses.
You now have a generic one-compartment model that has been applied to hundreds of drugs whose
parameters are tabulated in standard references.
To deliver a single pulse, use the PULSE function with an interval > STOPTIME. For example, set
interval = STOPTIME +1. Select the Parameter Window under the Parameters menu to be certain that
STOPTIME = 10 and run the model. You should get a single injection at the beginning.
Now change the STOPTIME in the Parameter Window to 240 (10 days), and run the model. Use the
Pulse function to simulate different practical regimens (i.e. 1, 2, 3, and 4 times per day. The
corresponding intervals are 24, 12, 8, and 6 hours.) Set the Interval and vary the Dose to find a
dosage that, after a few days, will keep the concentration, C, within the safe-but-effective range of
10 −20 mg L
Madonna’s sliders offer a convenient way do that.
Open the Define Sliders dialog (in the Parameters menu), double click on dose to add it to the
Sliders list, and then click OK. Move the slider back and forth: when you release the mouse,
Madonna automatically runs the program with the current parameter value shown on the slider. See
Help>How Do I >Define and Use Sliders for details on how to tailor the slider to specific needs.
Oral Ingestion
Figure 7. Theophylline Drug Model: Oral Ingestion
So far, we have used a single, well stirred, compartment to represent the entire body. Most drugs,
including theophylline, are given orally, not injected directly into the veins. Thus they pass first into
the GI tract and then into the circulatory system. To simulate oral administration, add an additional
compartment, labeled GI , between the Dosage inflow and the m compartment. You will also need a
new flow representing absorption (labeled E gm ) between the GI and m compartments. Choose the
flow icon, click on the GI reservoir and drag the mouse until it contacts the m reservoir. Set
Equation 6 Ex
gm
= k
gm
m
where k gm
is an absorption coefficient. Set k gm
= 0.1 and run the model. Experiments with
theophylline show that its plasma concentration reaches a peak value in about 2 hours following oral
ingestion. Reset stoptime back to 10 so that you can see the details of a single injection. Then set up
a slider to adjust the value of k abs
so that the peak value in the model will appear at an appropriate
time. Change stoptime to 240 and find the best dose to use at different dosage regimens (1, 2, 3,and 4
times per day).
Two Compartments
Salicylic Acid
Although the last model includes an extra reservoir to accommodate oral ingestion, the lumen of the
GI tract is essentially part of the delivery system; it is continuous with the external environment and
is not considered part of the body. Following an intravenous drug injection, a single compartment
model predicts that the drug plasma concentration will decay exponentially. In a number of cases
Tissue Region Blood Flow mL/100g /min
Liver 57.
Kidneys 420.
Brain 54.
Skin 12.
Skeletal Muscle 2.
Heart Muscle 84.
Rest of Body 1.
Whole Body 8.
Table 1. Regional Blood Flow Density (ml/100g tissue/min)
The fast and slow phase for salicylic acid can be modeled by defining two compartments (reservoirs)
in series as shown in Figure 9. The first, central compartment , consists of the plasma together with
all rapidly exchanging tissues that are assumed to be equilibrated at all times. The second, peripheral
compartment is composed of the slowly exchanging compartments. Since drugs are generally
eliminated by kidney excretion and by metabolic inactivation in the liver, the route of elimination
emanates from the central compartment.
Figure 9. Two compartment
model for salicylic acid ( SA ).
The formula icon SA 1 norm
represents the value of SA 1
normalized by the Dose , i.e.
its value equals the fraction of
the original dose that appears
in the central compartment at
any time.
Although it is tempting to assign specific organs to each compartment, it is usually not feasible and
the anatomical identity of the two compartments as well as their associated volumes are generally
left as abstract elements of the model. The model shown in Figure 9 circumvents the ambiguity in
the volumes by dealing explicitly with amounts (e.g. mg) rather than concentrations. There is no
problem identifying quantities in the reservoirs as amounts; this is consistent with all of the models
we have considered. However transfer between compartments as well as rates of elimination are
determined by concentration, not by amount. This is easily resolved as follows:
We expect the transfer between compartments 1 and 2 to be proportional to the concentration
difference, ( C 1 – C 2 ) between the compartments. If the proportionality constant is denoted by k , and
the distribution volumes of the two compartments are V 1 and V 2 , then the flow between 1 and 2,
denoted jSA , will be given in terms of the drug masses SA as
Equation 7 jSA = k C 1
− C
2
( )
k
V
1
SA
1
k
V
2
SA
2
= k 12
SA
1
− k 21
SA
2
where we have replaced k by two rate constants defined by k 12
= k/V 1
and k 21
= k/V 2
. Similarly, we
assume Ex = k ex
C
1 to arrive at
Equation 8
Ex = k 1 x
SA
1
where k 1x = k ex
/ V
1
Eliminating the concentrations has the advantage of reducing the number of unknown parameters;
we have replaced four unknown parameters ( k, k ex
, V
1
, and V 2
) with three ( k 12
, k 21
, and k 13
Guess the Parameter Values
Using Equation 7 and Equation 8 together with initial conditions SA 1 = SA 2 = 0 , set up the model
illustrated in Figure 9. We will use the curve fitting routine in Madonna to find values for the rate
constants from the experimental data shown in Figure 8. To do this, we will need initial guesses for
these parameters. The time scale in the figure runs from 0 to 500 minutes suggesting that residence
times of the order of 250 min ought to produce results that change within the required range, Set k 12
= k 21
= k 13
= 1/250 = 0.004 min
, and run the model. Your results should look like the blue curve in
Figure 10.
Import the Data:
- Select Import Dataset from the File menu. The
standard Open dialog will appear.
- Open the enclosed data file called ‘Salicylic Acid’.
This is a simple 2 column text file; the first column
contains values for the x (time) axis, the second
contains y ( SA 1 ) values. A dialog will appear. It
describes the data file, allows you to rename it, and
allows you to choose which column represents the
independent (x axis) variable and which represents y.
Numerical data files in Berkeley Madonna are
distinguished by the prefix #.
- Click OK
The data is now superimposed on the plot. The 13 data points should be represented by open
unconnected circles. If your results show data points connected by lines, then open Edit >
Preferences. Click on the Graph Windows tab and check Plot Imported Data with Circles.
Notice the large discrepancy between the data and the model run with the guessed parameters.
Find the Parameters by Curve Fitting the Data:
- Click Curve fit in the Parameters menu to bring up the dialog shown below.
Figure 10
Aspirin
Our model for salicylic acid provides one part of the more complex model for aspirin (acetyl
salicylic acid, ASA ). When ASA enters the blood plasma it is rapidly hydrolyzed into salicylic ( SA )
and acetic acid by esterases contained in both blood and tissues. Both SA and ASA are effective
analgesic, antipyretic, and anti-inflammatory drugs (see Error! Reference source not found. )
Figure 11 Pharmacological Effects of SA and ASA
Our model for an iv injection of aspirin is shown in Figure 12. If ASA was not transformed into SA ,
their models would be almost identical. However, the transformation is rapid and irreversible; it
provides the principal route J 31 for elimination of ASA (elimination by other means, e.g. excretion,
are negligible compared to J 31
so they are left out of the model). Save a copy of your SA model and
then use the original to build the aspirin model below.
