After starting Boomer or MultiForte (by double-clicking on the Boomer or MultiFort icon-Mac, typing Boomer or Multif - IBM, or run Boomer or MultiF - VAX) you should see the title, copyright information, and references to the original MULTI.
The method of data entry is then chosen from the menu option.
DATA ENTRY 0) From KEYBOARD 1) From .BAT file -1) From .BAT file (with restart) 2) From KEYBOARD creating .BAT file 3) From .BAT file (quiet mode) -3) From .BAT file (quiet mode-with restart) 4) to enter data only 5) to calculate AUC from a .DAT file -9) to quit -8) Registration Information
Data can be entered from the keyboard or from a batch file (verbose (1) or quiet (3) modes). Alternately a batch file can be prepared while entering the data from the keyboard.
Boomer can be started (Mac OS X and DOS version) from the command line using the command
boomer filename quiet mrun where filename is the BAT filename, quiet (or 3, q or Q) indicates quiet mode, and mrun is the maximum number of runs in multi run mode.
0) From keyboard
All subsequent data will be entered from the keyboard.
1) From .BAT file
Subsequent data will be read from the batch file (named xxx.BAT). Only the xxx part of the name should be entered as the program will add the .BAT extension. As with later filename entries the extension is not to be entered.
Boomer is looking for -> xxxx.BAT <- Do you have the right file name? PAUSE Filename incorrect - try again? statement executed To resume execution, type go. Other input will terminate the job.
If a batch file (named xxx.BAT) has been prepared the data can be read from this file instead of from the keyboard. Only the xxx part of the name should be entered as the program will add the .BAT extension. As with later filename entries the extension is not to be entered. DOS version. If the file is not found (or spelled incorrectly) a PAUSE message is displayed which allows the user to enter a DIR or other DOS command.
For example, DIR.*.BAT would display all the .BAT files in the current directory. The directory or drive could be changed with a CD xxx or A: type command. This is then followed by another request for the file name (less extension).
-1) From .BAT file with Restart Option
If a numerical integration method (especially RKF or Gear) is used with fitting problems it is possible that parameter values may be tried which cause 'unstable' integration with the program stopping. With Restart ON (checked) it is possible for the program to restart automatically using parameter values which are stored to disk during the iteration process (in a file called Scratch.PAR). This option is especially useful with unattended operation using large systems of differential equations. The current best parameter values are stored periodically as Scratch.PAR.
2) From keyboard creating .BAT file
Again a .BAT file name is requested. Data will be entered from the keyboard and used both in the analysis and to create a .BAT batch file. Be careful not to use the name of a .BAT file previously created and needed later. The program will check to confirm that you want to over-write files having the same name.
3) From .BAT file (quiet mode)
Subsequent data will be read from the batch file (named xxx.BAT). In this mode subsequent data entry is not displayed on the screen. This can result in faster operation.
-3) From .BAT file with Restart Option (quiet mode)
4) to enter data only
This choice allows the user to enter data into a xxx.DAT file for subsequent analysis.
5) to calculate AUC from a .DAT file
This choice allows the user to calculate AUC from data in a xxx.DAT file.
METHOD OF ANALYSIS 0) Normal fitting 1) Bayesian 2) Simulation only 3) Iterative Reweighted Least Squares 4) Simulation with random error or sensitivity analysis 5) Grid Search -5) To perform Monte Carlo run (Start of BAT file) -4) To perform multi-run (End of BAT file) -3) To run random number test subroutine -2) To close (or open) .BAT file -1) To finish Enter choice (-3 to 5)
0) Normal
This option allows a normal weighted non-linear regression analysis of the data set.
1) Bayesian
A Bayesian analysis is also possible. The observed data may be weighted. During data entry population mean and standard deviations are required for each adjustable parameter.
2) Simulation only
This option is useful for models consisting of large systems of differential equations. That is, systems taking a long time to solve or for which fitting is not required. This option could be useful in the analysis of systems containing fast and slow rate constants (so-called 'stiff equations'). These are commonly observed with physiological models.
3) Iteratively Reweighted Least Squares
With this option the weighting functions use the calculated 'y' value at each stage of the fitting rather than the observed 'y' values. For example using a 1/value weighting function with option results in a weight calculated as 1/ycalc, not 1/yobs. This can reduce some of the strange fitting that results when very small 'y' observed values are included in a data set. These values may be given a weight much larger that they really deserve [Peck, 1984a; Peck, 1984b].
4) Simulation with random error or sensitivity analysis
This option allows the simulation of data sets with added error. The error can be calculated as proportional or independent of the y-value and can be calculated as uniformly distributed or normally distributed. Log normal or exponential error can also be added to data. This same uniform, normal, log normal and exponential error can be randomly applied to any of the parameter values. See section Q, below, for more details. A repeat (Monte Carlo) option will allow a large number of simulations to be run, creating different data sets with each run.
In the sensivity (or optimal sampling) mode the program generates graphs of y-value sensitivity to each of the adjustable parameters. Each point of maximum sensitivity (in a 100 step increment) is also reported. The program calculates the sensitivity value, calculated using the equation below for each parameter.
5) Grid search method
This option allows a systematic search of a specified parameter space. The space is defined by specifying upper and lower limits and the number of steps in the grid. Output can be in the normal form or as disk files of parameter values and weighted sum of squares values for input into a graphing program for 3-D surface plotting. The Macintosh version will automatically create gnuplot control files and create a QuickTime™ movie of the surfaces.
-5) Monte Carlo run
This is similar to the multi run option below with one significant difference. If data files are input or output they will be created with unique (incremental) names. They can be useful in the generation of many multiple data sets that might be useful as a clinical trial simulation. A Monte Carlo run .BAT file can also be created from a 'ordinary' .BAT file by adding the number of runs (>10) as the second line in the .BAT file.
Monte Carlo run Option
It is possible to run repeated analyzes of 10 to 99999 data sets (memory and time permitting). This could be be useful in a Monte Carlo type analysis. I have found this useful for generating homework data sets. It could also be used to simulate data sets for other purposes such as clinical trial simulations. An additional output data file (extension .NDT) is generated with the same name as the 'normal' output (.OUT) file but in NONMEM input data file format.
These runs can be set-up during the creation of the .BAT file as above or by creating a single run batch files as normal with -2 at the end. Open the .BAT file in a text editor and insert above the first data line the number of repeats required (must be greater than 10).
Boomer Batch File
4 wls,bayes,sim,irwls,sim+error,grid
Top of the .BAT file before editing and
Boomer Batch File
25
4 wls,bayes,sim,irwls,sim+error,grid
after editing to simulate 25 data sets.
Save the file and select batch file run. When simulating data sets use the simulate with error option to make different data sets on each run. When saving the calculated data set (the first time through when creating the batch file) use any reasonable, short file name (e.g. test: A maximum of 4 characters are used with the DOS version). In repeat mode the program will create file names with a numbered extension. For example test0001.DAT, test0002.DAT, test0003.DAT ... in sequence. When reading the data sets in using the repeat option the same naming convention is used. That is specify 'test' and the program will expect to read files test0001.DAT, test0002.DAT, test0003.DAT etc. [There are four columns in the number field to accommodate file names from test0001.DAT to test9999.DAT].
The auxiliary program cDAT can read these multiple DAT files and combine them into a tab delimited summary file suitable for importing into Excel(R) or OpenOffice. The cDAT program is installed by the Mac OS X boomer installation package in /usr/bin. The program can be invoked by typing cDAT and then providing the DAT file name and the number of files. It can also be invoked as one command.
cDAT test nnn
where 'test' is the first part of the DAT file name (can be entered as test0001 or the complete name test0001.DAT) and 'nnn' is the number of DAT files.
The auxiliary program eOUT can extract specified lines from .OUT files (especially useful for Monte Carlo or Multi runs), optionally converting spaces to tabs. The output file is suitable for importing into Excel(R) or OpenOffice. The eOUT program is installed by the Mac OS X boomer installation package in /usr/bin. The program can be invoked by typing eOUT and then providing the OUT file name and other parameters. It can also be invoked as one command.
eOUT test.OUT key nskip nline tflagwhere test.out is the .OUT file name (test should also work), key is the string to seach (enclose with '' if it contains spaces), nskip is the number of lines to skip after key is found, nline is the number of lines to export and tflag is 0 or 1 (if equal to 1 converts run of spaces to one tab). It is useful to examine the .OUT file before running eOUT to determine the best values of the parameters, especially key, nskip and nline.
The auxilliary program eBEST can extract the best run from a multi-run .OUT file. The eBEST program is installed by the Mac OS X boomer installation package in /usr/bin. The program can be invoked by typing eBEST and then providing the OUT file name. It can also be invoked as one command.
eBEST test.OUTwhere test.out is the .OUT file name. The program then asks
Enter percent from Best WWWEnter a percent (x) and the program displays the number of runs with a WSS within x% of the best. The best run is placed in a file BESt.OUT which is opened by BBEdit (if installed) for review and saving.
-4) Multi run
This option allows either repeated runs of the BAT file or chaining to another BAT file. It will not be selected initially but only at the end of the first run. Repeatedly running a BAT file is useful when fitting data using the Simplex method as each run starts with a randomly generated simplex. A multi run .BAT file can be created from a 'ordinary' .BAT file by changing the -1 or -2 on the last line to -4 and entering the .BAT file name on the following line.
-3) Random Number Test
The Macintosh toolbox routine, random, is used to generate random numbers for the simplex method (to create the starting simplex), and the simulation with error values. This option gives the user the opportunity of testing these random numbers. The numbers are output to disk files. Random produces uniformly random numbers (URN) in the range -32,768 to 32,767 which are scaled to real numbers between the value of 0.0 and 1.0. Normally distributed random numbers (NDRN) with mean of zero and standard deviation of one are calculated using the equation [Abramowitz, 1964, page 953]:
-2) To close (or open) batch files
This option allows the user to break a batch file run down into smaller parts. Each batch file could then be combined in the editor. Alternately the first data set could be analyzed using the make batch file mode and closed with a -2 without exiting the program. The editor could then be used to edit this newly created batch file to analyze the remaining data sets. A -2 is put on the end of any batch file created. Any printer output can be sent to the printer at this point by choosing `Send output to Printer' from the `Boomer' or `MultiForte' menu. Batch files can be chained by manually editing the .BAT file with the editor. Replace the -1 or -2 at the end of the batch file with a -4. On the next line put the name of the batch file you want to run next. Enter the name without the .BAT extension. Make sure that there is a carriage return after the .BAT file name. All output will go to the .OUT file specified in the first .BAT files.
-1) To finish the program
All normally finishing runs will come back to this choice. If no further analyses are required a -1 is entered. With a batch file run this -1 is put on the last line of the batch file in the expected sequence. If editing a batch file remember to ensure that a carriage return must be present in the batch file after the -1 for it to be correctly interpreted. Actually all batch, data, or model files must end with a carriage return to be interpreted correctly.
Where do you want the output? 0) Terminal screen 1) Disk file Enter choice (0-1)
0) Terminal screen
All output is directed to the computer terminal screen. This is useful for a quick check of the input data or for rapid analysis which doesn't need to be kept.
1) Disk file
Again a disk file name is requested after selecting this option. The final results and requested plots can be directed to an .OUT file output file for later printing or enhancement. The output is saved as a text type file which could be read into a word processor for editing or other visual enhancement. This option can be particularly useful for receiving the results from a large series using batch input if there is any problem with a printer jamming while unattended. Of course the disk receiving a large output file will need to have enough space.
One limitation of the original implementation of MultiForte was that the user had to have and know how to use the Fortran compiler in order to study different pharmacokinetic models. This is still an option with the original version of MultiForte. A later section 'Demo Version' describes one such compiled program which will handle a limited number of different models. Later sections of this manual briefly describe the subroutines CP, MODEL, MODELOUT, and PART which need to be written by the user studying more involved models.
Boomer, a version of Multi-Forte described in this section, allows the user to describe a model in terms of time interrupts (e.g. dose times or lag times), dose, first order rate constants, zero order rate constants, constants, Michaelis-Menten kinetics elimination rate systems [Portmann, 1970; Michaelis, 1913], multipliers (i.e. inverse apparent volumes of distribution), sums of exponentials, Hill equation components [Rubinow, 1975; Hill, 1910], second order rate constants, or simple physiological models [Duddy, 1984]. This building of a model is similar to the approach taken with the mainframe simulation and modeling program, SAAM [Berman, 1968].
MODEL Definition and Parameter Entry * Allowed Parameter Types * -5) read model -3) display choices -2) display parameters 0) Time interrupt 1) Dose/initial amount 2) First order rate 3) Zero order 4-5) Vm and Km of Michaelis-Menten 6) Added constant 7) Kappa-Reciprocal volume 8-10) C = a * EXP(-b * (X-c)) 11-13) Emax (Hill) Eq with Ec(50%) & S term 14) Second order rate 15-17) Physiological Model Parameters (Q, V, and R) 18) Apparent volume of distribution 19) Dummy parameter for double dependence 20-22) C = a * SIN(2 * pi * (X - c)/b) Special Functions for First-order Rate Constants 23-24) k = a * X + b 25-27) k = a * EXP(-b * (X - c)) 28-30) k = a * SIN(2 * pi * (X - c)/b) 31,32-33) dAt/dt = - k * V * Cf (Saturable Protein Binding) 34-36) k * (1 - Imax * C/(IC(50%) + C)) Inhibition 0 or 1st order 37-39) k * (1 + Smax * C/(SC(50%) + C)) Stimulation 0 or 1st order 40) Uniform [-1 to 1] and 41) Normal [-3 to 3] Probability 42) Switch parameter 43) Clone component 44-47) Four parameter logistic model Enter type# for parameter 1 (-5 to 47)
With this approach the model is constructed step-by-step during the data entry process. Once a useful model is chosen it may be more convenient to use the batch mode, by first preparing a batch file, editing, and running the program using the edited batch file (see Data entry section 1 and 2). The types of parameters available are shown in the screen dump above.
-5) Read model from .mdl file
The model, described with adjustable or fixed constants, can be defined using a .mdl file. If a file has been previously defined it could be read at this point of the data entry.
-4) Write model from .mdl file (available after the model is defined)
The model, described with adjustable or fixed constants, can be defined using a .mdl file. Once the model has been defined a .mdl file can be written for later use.
0) Time interrupt
These are times during a simulation at which you may wish to make a change. For example, give a second or subsequent dose, start or stop an infusion, or incorporate a lag time into the calculation. The information requested for this parameter is simply a name and a value. It is necessary to specify any interrupt time before it can be used (e.g. to turn off an infusion). Also, when you request a list of time interrupt values they will be numbered excluding other parameter types. It is easier to specify all the parameters that don't change first then the time interrupt values, and finally the doses or other parameters that depend on a time. As with any of the parameter types these time interrupt values may be fixed, adjustable, or dependent.
1) Dose/initial amount
This would typically be the dose or it could be an initial concentration. Given as a dose it would generally be fixed, however, as an initial concentration it may be more appropriate to 'best-fit' this value. Name, value, and component to receive the dose should be specified.
2) First order rate constant
By specifying various rate constants it is possible to build a system of differential equations for simulation of the model required.
The model, shown above, would involve three first order rate constants, k10, k12, and k21. Each of these first order rate constants could be specified by name (k10), value, component to receive flux (0) and component to lose flux (1). In this model drug would flow from the central compartment (Component 1) to the outside world (the 'sink' designated as Component 0). The remaining rate constants can be specified in a similar manner. The data entry for k12 is shown below. Connection between the components of the model (specifically component 1) and the data to be fitted or simulated is made with the volume term. In this case the plasma data is equal to the amount in component 1 divided by the volume (of distribution).
Enter -3 to see choices, -1 or -4 (save model) to exit this section
Enter type# for parameter 2 (-5 to 43) 2
Enter parameter name k12
Enter k12 value 0.13
0) fixed, 1) adjustable, 2) single dependence
or 3) double dependence 1
Enter lower limit 0
Enter upper limit 1.0
Enter component to receive flux 2
Enter component to lose flux 1
Input summary for k12 (type 2)
Initial value 0.1300 float between 0.000 and 1.000
Transfer from 1 to 2
Enter 0 if happy with input, 1 if not, 2 to start over
3) Zero order rate constant
Zero order rate constants would normally be used to describe an infusion dose from one time interrupt value to another. A start time of 0 specifies that the infusion starts at the beginning of the simulation. A stop time of 0 represents an infusion which continues until the end of the simulation. As with the first order rate constant, components to receive and lose flux are specified. Typically, the infusion comes from a component which does not need to be modeled, thus a value of 0 (the outside world) can be used for the "from" component.
4-5) Michaelis-Menten Kinetic Parameters, Vm and Km
Rate processes can be specified as Michaelis-Menten processes, for example the rate of elimination may be as shown in the equation below.
where Vm is the maximum rate and Km is the Michaelis-Menten constant. X(1) is the amount (or concentration) in the from component. The Vm and Km are entered as two separate parameters, however, the to/from information remains the same for both. Be careful of the units. If doses (and thus the values in the components of the model) are given as mg or g, then, Vm units should be mg/hr or g/hr (mass/time). Dose/initial concentrations given as mg/L will mean Vm units of mg/L/hr (mass/volume/time).
The second half of the Michaelis-Menten entry process is the Km value. Only a value need by entered, unless the parameter is adjustable or dependent. Again, be careful of units. With a dose (in g or mg) entered the units of Km should also be g or mg. However if a Cp(zero) is used with mass/volume units, Km units should be the same.
6) Constant
This is the first of the summing or intensity terms. It is necessary to distinguish between the components (and their numbers) used to describe the system of differential equations and the observed data sets (and their numbers). For example in the pharmacokinetic model above there are two components not including the 'sink'. As shown, there is one data set, Plasma. We could also specify another, second data set called Tissue, and a third called Rate of Change. When a constant is specified only the name, value and data set line is required. We might use this to offset a fixed (or adjustable) background level? Or it could be used as a constant term in a Hill equation sequence or a linear fit sequence.
7) Kappa (Reciprocal volume of distribution)
Typically, if doses are entered in mg or g, then the kappa term becomes equivalent to the reciprocal of the apparent volume of distribution (see parameter type 18). That is:
OR
where X(1) is the amount in component 1 at the time specified. The units of kappa are reciprocal volume. If the dose/initial concentration is entered as an adjustable parameter the kappa term should probably be 1 and fixed.
We could also use this type of parameter with a constant term to fit a straight line. For example
X(observed) is specified by entering 0 for the component to be added.
Specifying a negative value for the 'from' component allows the simulation of rates of change. For example, urine data may be collected (and analysed) as rate of excretion. In this case Kappa would probably be fixed with a value of 1 but with the component from specified as -x, where x is the component representing the cumulative amount in urine.
8-10) Sum of Exponential Terms
Rather than specify a model in terms of components and rate constants (differential equations) it may be convenient (and generally faster) to use a sequence of exponential terms. The pre-exponential term (a) is the parameter before the exponent term. A two exponential equation could be:
The pre-exponential terms are 'A' and 'B'. Again, the Xobserved is specified by entering 0 for component from. The lambda or exponent terms (b) in the above equation are [[alpha]] and [[beta]]. A lag time can be specified for each exponential term. Usually these lag times would be the same. Second and subsequent lag times could be made dependent on the first if adjustable.
11-13) Emax (Hill Equation) Parameters
Effect versus concentration data can modeled using this sequence of parameters. The concentration data could be entered as Xobserved values, or they could be calculated from a system of differential equations. The full equation (Hill equation) is
where X(i) is the concentration in component 'i'. Specify component i = 0 if concentrations are entered as X-values (shown as time). Ec(50%) is the concentration which produces an effect of 50% of maximum. The third parameter is the power or slope term. You may want to try giving this a value of 1, fixed, before doing a full fit with all three parameters adjustable.
When plotting the Hill equation as Effect versus log(Concentration) it has been noticed that for the 20 - 80 % of maximum effect region a straight line may describe the Effect versus log(Concentration) well. Thus the simpler equation:
may be useful. Since version 2.7.9 Boomer will allow this equation by modifying the way parameter type 11 works. The third parameter, slope parameter, is used as the switch. If slope is greater than 0 the Hill equation is used in the calculation. If slope is -1 (remember to not let this parameter be adjustable) then a is entered instead of Emax, and b is entered instead of EC(50%).
For both versions of parameter type 11 a baseline 'effect' may be worth considering.
14) Second order rate constant
Second order reaction processes can be described using this parameter type. Input includes one or two components for reactants or products.
15-17) Physiological Flow Parameters
Physiologically based models can be defined by specifying blood flow rates between organs or vein/artery blood pools. A flow balance is calculated and displayed if any blood flow parameters are specified. Generally flows in and out of an organ will need to be specified. Volume and extraction ratio data which are subsequently entered describe characteristics of the source of flow.
Drug movement from this source is calculated using the differential equation segment:
where X(i) is the amount of drug in the source (component i), X(i)/V(i) is the concentration of drug in the source tissue or blood pool, and X(i)/(V(i)*R(i)) is the drug concentration in the blood leaving the source tissue. In the case of vein or artery blood pool the value of R(i) should be 1. In addition to a flow term (15 above) a volume term is specified for the source organ or blood pool. Thus for blood flow from kidney to vein, the volume term would be for the kidney tissue. This volume term is used in the differential equation above and also to convert the amount values calculated internally to calculated concentration values. If you don't need to calculate values for a particular tissue you can enter 0 for the data set number. The extraction ratio, R, for the source tissue is the final term specified. For vein or artery blood pools the value of R should be 1.
18) Apparent volume of distribution
Typically, if doses are entered in mg or g, the volume term can be used to relate the amount term, dose, with the concentration term measured. That is:
where X(1) is the amount in component 1 at the time specified. If the dose/initial concentration is entered as an adjustable parameter the volume term should probably be 1 and fixed.
19) Dummy parameter for double dependence
This parameter is not used in any calculations of y-value except for the calculation of double dependent parameters. One example might be the estimation of the bioavailability factor F. The initial amount administered may be expressed as F * DOSE where F is unknown and to be estimate while DOSE is fixed and known. Thus F and DOSE can both be entered as dummy parameters, one adjustable and the other fixed. The term F*DOSE (a DOSE/Initial amount, type 1 parameter) would be specified as double dependent with double dependence type 5.
20-22) SIN Function
These parameters (type 20-22) can be used to incorporate diurnal variation in a line function. Parameters include Amplitude, Period and Offset. The result of these three parameters is the addition of the value: Amplitude x SIN (2 x [[pi]] x (X - Offset)/Period) to the selected data set (line). This would result in a value varying between +/- Amplitude. Normally a constant type 6 parameter (possibly adjustable, etc) would be also be used for the same data line so that the value would vary between Constant +/- Amplitude. With both parameters defined the equation for line 1 becomes:
23-24) First Order Rate Constant as a Linear Function
On occasion it may be useful to add a rate constant to a model which varies with time (or the value in a component of the model) in a linear fashion. Thus a first order rate constant can be defined as: a * X + b. The `a' parameter is a slope term and `b' is an intercept. The X term can be the observed x, usually time or the value (or reciprocal - entered as a negative number) in any component. Thus zero or second order rate constant can be entered with this parameter type. For example, a zero order rate constant of 5.0 could be entered as a(SL) = 5.0 with drug from component 1 to 0, b = 0.0 and the F-dependence equal -1, the from component. The differential equation for this process is then:
k is now a zero order rate constant.
25-27) First Order Rate Constant as an Exponential Function
Under some circumstances the value for a rate constant may change progressively with time. For example an elimination half-life may change from 8 hours (kel = 0.0866 hr-1) to 6 hours (kel = 0.1155 hr-1) over a period of 60 hours (with a half-life of 12 hours; k = 0.05776 hr-1) because of enzyme induction. This could be simulated with two rate constants. The first would be a type 1 parameter with a value of 0.1155 hr-1. The second parameter would be type 25 with a = -0.0289 (=-( 0.1155 - 0.0866)), b = 0.05776, and c = 0). Thus the rate constant would change from 0.0866 to 0.1155 hr-1 over a period of 60 hours (5 half-lives). The differential equation would become:
or
F-dependence is also possible with this type of parameter. Specifying a positive (or negative) component for F-dependence multiplies (or divides) the above term by the amount in the specified component.
28-30) First Order Rate Constant as a SIN Function
On occasion it may be necessary to include diurnal variation in the rate constants of a model. Again by using two rate constants, one type 1 and the other type 28, it is possible to simulate or analyze this phenomenon. With just the type 28 parameter the differential equation would be:
F-dependence is also possible with this type of parameter. Specifying a positive (or negative) component for F-dependence multiplies (or divides) the above term by the amount in the specified component.
31-33) Saturable Protein Binding Elimination
Many drugs are bound to plasma protein. With large dose drugs with small apparent volumes of distribution this protein binding may be saturated. Further, for drugs which also have low extraction ratios, that is, clearance is not flow limited, elimination of total drug is directly related to free drug concentrations [Bourne, 1981; Bevill, 1982]. Because of the complexity of this model its specification is relatively rigid. For example, considering i.v. bolus administration with saturable protein binding controlled elimination a dose, amount, would be added to component 1. The first-order rate constant is entered with units of reciprocal time. The association constant, Ka, has units of reciprocal concentration. The total protein, Pt, is entered as a concentration. The amount of drug is converted into total concentration with a single apparent volume of distribution (type 18) linked to the protein binding parameter, type 31, through the `from' specification for both parameters. If required a second data set line equal to the free drug concentration can be specified after entering the Ka value.
By specifying parameter type 32, Ka and Pt parameters can entered and thus Cfree versus Ctotal data sets can be analyzed. This can be useful for the analysis of in vitro data sets.
34-36) First or zero order rate constants subject to inhibition
37-39) First or zero order rate constants subject to stimulation
40) Uniform [-1 to 1] Probability
This dummy parameter can be used for single or double dependence and provides the value 0 or 1. If the parameter value (in the range -1 to 1) is greater than a random number generated uniformally distributed between -1 and 1 it provides the value 1, otherwise it provides a value of 0. One use for this parameter might be the simulation of compliance. For example, if 20% of the doses are not taken you might use a parameter value of 0.6 (= 2 x (1 - 0.2) - 1 = 2 x (1-x) -1) where x is the 'fractional' probability for 0 or non-compliance.
41) Normal [-3 to 3] Probability
This dummy parameter can be used for single or double dependence and provides the value 0 or 1. If the parameter value (in the range -3 to 3) is greater than a random number, normally distributed with mean 0 and standard deviation of 1, it provides the value 1, otherwise it provides a value of 0. One use for this parameter might be the simulation of compliance. For example, if 68.3% of the doses are taken you might use a parameter value of 1 (= one standard deviation above the mean), if 95.4 use a value of 2).
42) Switch parameter
Similar to parameter type 40 or 41 except change depends on the value of a specified parameter. This parameter is set equal to 1 or 0 depending on the value of the control parameter.
43) Clone component
Clone a component using this parameter. The new component will have the same numerical value as the 'donor' component but will not influence the donor system. Uses ??
44-47) Four parameter logistic model
This selection allows the addition of the four parameter logistic model as describe by Kraupp et al., 1986. Normally this would be the only parameter selection for the analysis of radioimunoassay assay data. The equation is shown below.
The parameter A represents Bo, B is a shape factor, C is the midpoint of the curve (ED50 and D represents non specific binding (Nsb).
E. Specifying Parameter Values and Limits
0) Fixed Parameters
If the method selected is a simulation, then all the parameters will be fixed and you won't be asked to specify fixed, adjustable, dependent, or double dependent. Alternately if a normal or bayesian fitting is to be performed some parameters can be specified as fixed. Often dose and time interrupt values will be fixed. You may wish to hold some rate constants fixed for one run and then allow them to be adjustable in a subsequent run. Enter 0 (or press return) to specify fixed parameters.
1) Adjustable parameters
At least two adjustable parameters must be specified for a normal fitting problem. If 1 is entered to specify adjustable, a lower and upper limit are requested. The fitting process will then keep the adjustable parameter values between these limits. To reduce the influence of parameter limits you could specify a very wide range between the upper and lower limit.
2) Dependent parameters
Occasionally a parameter (value) may be used more than once in a complete model definition. If these values are adjustable it is appropriate to let the first value be adjustable but set subsequent parameters equal to the first parameter (or the negative of the first parameter). This entry is shown in the figure. Entering 0 will produce a listing of the parameters already specified.
3) Double Dependent Parameters
Parameters can be dependent on two other parameters. This is usually accomplished using at least one dummy parameter, see above for an example.
After entering all the information required for each parameter the program will ask if you wish to make any changes. Entering 0 will bring up the list of parameter types on the screen for the next entry.
1) Straight Line: ycalc = intercept + slope * xobs
2) Sum of exponentials 
3) Two compartment pharmacokinetic model, with infusion dose
The above examples indicate the entries required to define the models specified. If the run is for a fit to data, additional information as described in the 'Specifying Parameter Values and Limits' would be necessary.
-----
Typically the user is simply asked to select the required model. In the subroutine model parameters such as number of lines (ln), differential equations (ndeq), parameters (m), and constants (nc) are defined as well as names for the model (modst), parameters (ps), and constants (cst).
H. Numerical Integration Method
If the number of differential equations in the model is greater than zero, then the next question is choice of numerical integration method.
1. Classical 4th order Runge-Kutta
This is a relatively slow method but reasonably robust. The relative error is controlled by decreasing the step size and reevaluation.
2. Runge-Kutta-Gill
This is another fourth order method. Again it is slow, and the relative error is controlled by adjustment of step-size. This is the method used in MULTI(RUNGE) but with automatic control of step size added in Multi-Forte.
3. Fehlberg RKF45
This method involves a fourth order evaluation with a fifth order used for efficient selection of optimal step-size. It is quite efficient and requires relative and absolute error values.
4. Adam Predictor-Corrector
This method is incorporated into the DIFSUB subroutine and can be used for non-stiff systems of differential equations. The overhead involved in this subroutine makes it less efficient for small systems of differential equations. Possibly it would be better with bigger systems. An absolute error term is requested for error control.
5. Gear with S/R Pederv
The Gear method is a very good method for solving systems of stiff equations. Again an absolute error value is required. Also with this option a subroutine PART must be described by the user for use by the subroutine PEDERV. The partial derivative of each equation with respect to each dependent variable is entered in subroutine PART. Thus if there are n differential equations and therefore n dependent variables, there will be n*n terms defined in PART. Commonly many of these terms will be zero, but they must all be defined.
6. Gear without Pederv
This is also a good method for solving stiff equations. The subroutine PART is not required for this method as the partial derivatives are determined by numerical differencing. An absolute error term is also used in the step-size control process. This method may be slower than method 4 (with Pederv) but the PART subroutine is not required. Boomer automatically includes the PART subroutine so either method can be used.
I. Numerical Integration - Error Terms
Relative and/or absolute error values are requested if numerical integration is required. Since all calculations are performed in single precision (in this version) specifying a relative error of less than 0.00001 may cause problems. For the same sort of reason very small values for the absolute error are not recommended and may in fact cause a failure of the program. The default value for each parameter is 0.0001 which can be entered by pressing return. Larger values for the error term should also be avoided especially if the program is used to fit to a data set. With a tight convergence criteria (low PC value, see below) and a larger error term the iteration process may oscillate erratically. This can be readily induced using a non-stiff system of differential equations and Gears method of numerical integration for a fitting problem. A large absolute error term and low PC value also helps.
Unless the simulation method of analysis is chosen, the next request is for the specification of the fitting algorithm. These algorithms are as presented in the original MULTI programs, simply translated to FORTRAN.
1. Gauss-Newton
This is generally the fastest method, but it may diverge giving much slower times, or get lost giving incorrect answers. With good initial estimates for the parameters this can be a very useful method.
2. Damping Gauss-Newton
This is also a quick method with the tendency to diverge more or less under control. It may still get lost if given poor initial estimates. Rerunning the same analysis with different initial estimates can give more confidence in a consistently achieved minimum sum of squares. This method would generally be the method of choice.
3. Marquardt
In this method the Gauss-Newton approach has been altered to give better performance for oddly shaped sums of squares surfaces. For most problems it is slower than options 0 or 1.
4. Simplex
This method can be quite slow, although not too bad. It can be especially useful in a case with poorer initial estimates. It is supposed to be more robust. Better data and/or better initial estimates will also help. This method will occasionally fail to produce a best fit. As the initial simplex is produced randomly rerunning the analysis with the same initial estimates will produce a different analysis, with possibly the same end-point.
5. Simplex -> Damping G-N
With this analysis the simplex method is used at first because of its robustness. The final results from the simplex method are then used with the damping Gauss-Newton method for further refinement of the analysis.
One parameter method
When only one adjustable parameter is specified another fitting method is automatically selected. This method is essentially a damped Newton-Raphson method.
1. DT value
If a fitting method (i.e. not a simulation) other than the simplex method is chosen the program will request a value for DT. Entering 0.0 (or pressing return without an entry) will cause the default value of 0.001 to be used. The DT value is used in the calculation of the partial derivatives of the WSS with respect to the parameter values. This value is used as a fractional change in the parameter value during the calculation of the partial. (If the parameter value is zero, the DT value is used as an absolute change in value). Too large a number may lead to a distorted search for the minimum WSS, too small a value may over emphasize rounding errors.
2. PC value
If a fitting method is chosen (i.e. not a simulation) a value for PC will be requested. For the Gauss-Newton/Marquardt methods a default value of 0.00001 is suggested, whereas for the simplex method the default value is 0.0001. The PC value for the Gauss-Newton methods controls the convergence exit during the iteration process. If the relative change in the WSS is smaller than PC the iteration process is terminated normally. For the simplex method the convergence criteria is based on the relative change in the parameter values compared with the value of PC.
This is used to identify the data set. Up to 60 characters can be entered. The data set title is printed with the final output.
M. Individual Data Value Entry
Enter data from 0) Disk file 2) ...including weights 1) Keyboard 3) ...including weights
0) Disk File
If the data has been previously entered and stored in a data file it may be read by the program. This request is followed by a request for the .DAT filename. Again the .DAT part of the name should not be entered.
The data file format has changed from previous versions of MultiForte and Boomer. Data files now consist of two column data with a tab character between the `x' and the `y' value. A carriage return separates each data pair. This format is consistent with spreadsheet and graphing program formats. The current versions of Boomer and MultiForte will automatically read data in the older Boomer and MultiForte formats (as well as data stored in Cricket Graph(c) format).
1) Keyboard
With this option the individual time/concentration data are to be entered from the keyboard. Entering -1 for time will end this data entry.
2) and 3) are the same except the user enter estimates of the weight to be applied to each data point. This should usually be the reciprocal of the variance of the data point.
After data entry from the keyboard, or by choice after reading the data from a disk file, the data set may be viewed or edited. Editing is quite complete with options to change, delete, insert an individual data point or add an offset to the x-values of the data line. After data entry or editing the data set can be saved as a .DAT file by specifying the disk file name. Again reusing a file name will result in the loss of the older data.
If no weight is entered during the data entry process the program will ask for the user to select from various formulas. The relative error in each data point may or may not be known. The user may have some information about the relative error in each measurement. This may be independent of the value or proportional to some function of the value. The fitting process can be made aware of this information by use of a weighting scheme. Multi-Forte allows a number of different weighting schemes to be used.
Weighting function entry for [Drug] 0) Equal weights 1) Weight by 1/Cp(i) 2) Weight by 1/Cp(i)^2 3) Weight by 1/a*Cp(i)^b 4) Weight by 1/(a + b*Cp(i)^c) 5) Weight by 1/((a+b*Cp(i)^c)*d^(tn-ti))
0) Equal weight
An unweighted least squares fit can be specified by giving each data point an equal weight.
1) Weight by 1/Cp(Obs) or 1/Cp(Calc)
Each data point can be weighted by the reciprocal of the observed or calculated data value. This may be a reasonable alternative to the unweighted least squares.
2) Weight by 1/Cp(Obs)^2 or 1/Cp(Calc)^2
The square of the observed or calculated value is used in the weighting calculation.
3) Weight by 1/(a*Cp(Obs)^b) or 1/(a*Cp(Calc)^b)
This option allows non-integer powers of the observed or calculated concentration to be used in the weighting process. The values of a and b could be determined from the data by plotting log(average concentration) versus log(variance) from multiple subject data sets as proposed by Wagner [Wagner, 1975].
4) Weight by 1/(a + b*Cp(Obs)^c) or 1/(a + b*Cp(Calc)^c)
This is an even more flexible version of the above scheme.
5) Weight by 1/((a + b*Cp(Obs)^c)*d^(tlast - ti)) or 1/((a + b*Cp(Calc)^c)*d^(tlast - ti))
This allows incorporation of estimated assay sensitivity (a may
equal sensitivity
)
and coefficient of variation (b may equal CV)
as well as reducing emphasis on older samples (d, typically 1.01 to
1.05) as may be appropriate in a clinical setting. This type of
weighting scheme has been incorporated into some clinically useful
Bayesian estimation programs [Peck,
1984].
For simulations a weight of 1 is arbitrarily assigned to each data point. Observations with a value of zero are given a weight of zero. Selection of a 'best' weighting scheme is aided by consideration of the weighted residual plots which can be printed by Boomer or MultiForte.
If the MultiForte model chosen has associated constants, these are displayed by name and values requested.
P. Entry of Initial Parameter Estimates
1. Parameter values
With all analyses (except simulations with no parameters) initial values for each parameter are required. Generally the better the estimate the better the analysis.
2. Parameter limits
For all analyses (except for simulations) parameter limits, upper and lower, are requested. These limits can be wide or very close, it's up to the user. For those preferring no limits, very wide limits will suffice. These values are not used in anyway to change the parameter domain, but simply act as absolute limits for the parameter values. More like a box than a curved dish. When limits are specified some care may be necessary. For example, a lower limit for a volume of distribution value entered as zero may cause problems if the limit is approached during the early part of an iterative run. A value for volume of zero will give infinite calculated plasma concentrations and can cause an early exit to the fitting. A lower limit of zero for rate constants is usually better tolerated.
3. Population values
If the Bayesian estimation procedure is chosen then the user must supply estimates of the population mean parameter value. Also required is an estimate of the uncertainty in the population value. This is entered as the population standard deviation. If a larger uncertainty is to be used a larger value for the population standard deviation would be entered.
Q. Simulation with Error Input
Simulation with random ERROR 0) no error for line ( 1) 1) error = uniform error 2) error = normal error 3) error = uniform error x true value 4) error = normal error x true value 5) error = log normal error about true value 6) error = EXP(normal error) Enter choice (0-6) 6 Enter error intensity factor
If a simulation is chosen, parameter limits etc are not required. When simulation with error is specified information about the type and magnitude of the error is requested. Uniform, normally distributed, log normal and exponentially distributed error can be specified. This error may be simply added to the calculated value (choices 1 and 2) or multiplied by and added to the calculated value (choices 3 and 4). The intensity factor refers to the magnitude of the error. With choices 1 and 2 this would be the maximum error or the standard deviation of the error added. For choices 3, 4, and 6 this is a fractional values. For example, a value of 0.1 with choice 3 would mean that the added error would have a standard deviation of 10% about the calculated value.
1) Valn = Valo + Uni x Int
2) Valn = Valo + Nor x Int
3) Valn = Valo x (1.0 + Uni x Int)
4) Valn = Valo x (1.0 + Nor x Int)
5) Valn = 10log(Valo) + Nor x Int
6) Valn = Valo x e(Nor x Int)
where
Valn = new value
Valo = original value
Int = intensity factor (standard deviation)
Uni = uniformly distributed random number between -1 and + 1
Nor = normally distributed random number with mean zero and standard deviation equal to the intensity value
Log = log normally distributed random number
Exp = exponentially distributed random number
The program will automatically calculate the area under the curve (AUC), the area under the first moment curve (AUMC), and the mean residence time (MRT) of any data set. For accurate numbers a zero (dummy or otherwise) data point is needed. The program uses the trapezoidal rule to calculate the areas. The program uses the last two calculated data points to calculate a terminal slope for extrapolation to zero. Calculated zero time values are used is an observed zero time value is zero. If no zero time value is given a warning is issued. Zero time values with zero value are given zero weight in any fitting analysis. Enter 0 to exit this section.
S. Additional Output including Saving Data, Graphs, Extra Files and Sensitivity Ouput
Additional Output Enter -> 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Save Data x x x x x x x x Graphs x x x x x x x x Extra files x x x x x x x x Sensitivity x x x x x x x x Save calculated data files, Produce linear, semi-log and weighted residual plots, save extra data files or perform sensitivity or optimal sampling analysis Enter choice (0-15)
Enter 0 - Continue without additional output
This leads to the next analysis without any plotting.
Enter 1, 3, 5, 7, 9, 11, 13, or 15 - Save data on disk for later graphing
Select FORMAT for x-value
0) Don't save data 1) G14.1 2) G14.2 3) G14.3 4) G14.4 5) G14.5 6) G14.6 7) F14.0 8) F14.1 9) F14.2 10) F14.3 11) F14.4 12) F14.5 13) F14.6
Enter choice (0-13)
This option allows writing the calculated values and/or the standardized residuals to disk files (.DAT files) for later analysis. Suitable file name(s) are requested as before. These data files are output as x-value, tab, y-value, return. Via the clipboard, this format should be compatible with a variety of Macintosh programs. Many programs will read the .DAT files directly. The data format can be specified in the dialog box. Choices include a `G' general format with 1 through 6 digits or a `F' floating point format with 0 through 6 decimal places. [Overflow with the `F' format will produce 14 *'s]. The file name (without extension) should also be provided.
Enter 2, 3, 6, 7, 10, 11, 14, or 15 - Produce a printer plot now
Printer-type plots (not high resolution) can be prepared for quick evaluation of the fit to the data. Four types of plots are produced. Linear and semi-log observed versus time, standardized residual versus time and standardized residual versus observed y value. These plots can be useful in the evaluation of the model or the selected weighting scheme [Draper, 1966].
When there is more than one line of data the program also prints a combined plot for standardized residuals versus time or observed y value.
Example Plot Output
Plots of observed (*) and calculated values (+)
versus time for Cp . Superimposed points (X)
24.74 Linear 24.74 Semi-log
|* |*
|+ |+
| | +
| | *
| |
| + | X
| * |
| |
| | X*
| | +
| |
| X |
| | X
| |
| X | X
| * |
| |
| + | X
| |
| |
| X |
| |
| X | X
| |
| X |
| |
| | X
| X |
| |
| X |
| X | X
|_____________________________________ |_____________________________________
2.730 2.730
0 <--> 15. 0 <--> 15.
Example Standardized Weighted Residual Plot
Plot of Std Wtd Residuals (X) Plot of Std Wtd Residuals (X)
versus time for Heart_Rate versus log(calc Cp(i)) for Heart_Rate
1.742 1.742
|X | X
| |
| |
| X | X
| X X X |X X
|X X | X X
| |
0======================X============== 0X====================================
| XX | XX
| |
|X X X |XX X
|X | X
| |
| |
| |
| |
| X | X
| |
-2.586 -2.586
0.0 <--> 48. 67. <--> 94.
Enter 4, 5, 6, 7, 12, 13, 14, or 15 - Extra output files
Choosing this option produces extra output files. ***.htm, ***.csv and ***.mdl files. The htm file is a html formatted file for web display or input as external data by OpenOffice. The csv file can be read by a spreadsheet program such as Excel or OpenOffice. The mdl file contains details of the model with parameter values. The mdl file may be used as is for input or edited.
Enter 8, 9, 10, 11, 12, 13, 14, or 15 - Sensitivity or Optimal Sampling Output
With this option the program break the x scale (time) into 100 segments and calculate the sensitivity value according to the formula.
For kel
Max for line (A) Cp is 494.3 at 48.00
Max for line (B) Heart_Rate is 55.58 at 13.44
Max 494.3
| AA
| AAA
| AA
| AA
| AAA
| AA
| AA
| AA
| AAA
| AA
| AA
| AAA
| A
| AAA
| AAA
| AA
| AAA
| AA
| AA
| AA
| AA
| AA
| AA
| AAA
| AA
| AAA
| AAA
| A
| AAA
| AA
| AA
| AAA
| AA
| AA
| AA BBBBBBB
| AA BBBB BBBB
| A BB BBBBB
| A BBBB BBBBB
| BBBBBBBB BBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB
|B------------------------------------------------------------------------------
Min 0.000
0 48.00
The weight scheme for these data were:
Weighting for Cp by 1/a*Cp(Obs )^b
With a = 0.2500E-02 and b = 2.000
Weighting for Heart_Rate by 1/a*Cp(Obs )^b
With a = 4.000 and b = 0.000
Three Dimensional WSS Output Table
If the grid search method is chosen the program may output the WSS at each point on the specified grid with the parameter values to two disk files per pair of parameters. One disk file is in SYLK format for entry into the DeltaGraph(c) graphing program. The other file is in tab-return delimited format for entry into spreadsheets etc. With this output WSS surface plots can be generated.
1. Batch Files
Creating a batch file on the run could be useful if you have a large number of data sets requiring similar analyses. Alternately, a single data set may be analyzed by different weighting schemes. You could enter the first set of information from the keyboard, creating a batch file as you go. The batch file could then be edited and the section containing data repeatedly duplicated using the copy and paste commands. An entry of -1 is used to terminate the data entry. This will need to be changed to 0 at the end of each set except the last. (This final -1 must be followed by a carriage return or it will not be properly recognized.) When editing the batch file, be sure to change data file names and initial estimates for parameter values. Do not change the integer or real number description starting in column 40.