POPT - Installation and user guide.

 

Table of contents

Introduction to Design

Preamble about POPT

Installation

The POPT Suite of Programs

To run POPT and evaluate sampling windows

Input files for running POPT and evaluate sampling windows

Interpreting Output

P2NM

Optimization Algorithms

Copyright Statement

Bibliography

 

 

Introduction                                                                                     BACK TO TOP

(from WinPOPT User Guide)

 

A brief description of optimal design is provided here.  For those readers who are not familiar with optimal design an introduction is provided in: Monteleone JPR and Duffull SB.  Choice of best design. In: Kimko HC and Duffull SB.  Simulation for Designing Clinical Trials.  A pharmacokinetic-pharmacodynamic modeling perspective.  Marcel Dekker, NY. 2003, pp173-200. 

 

The purpose of optimal experimental design is to obtain data that provides the most information about a model and its parameter values from a new experiment.  This contrasts with modelling which proposes to provide the best, unbiased parameter estimates, given a model and a design from which the current data arose.

 

Optimal design is a phrase used to describe a general set of design methods that enable an optimum design to be constructed.  These methods are often described as information theoretic since they do not rely on the availability of real or simulated data.  However, it should be noted that these methods are widely used in other industries and are gaining popularity in drug development with numerous prospective applications of optimally designed studies.

 

The general approach involves computation of the Fisher information matrix (see bibliography for selected relevant readings).  For nonlinear models, the Fisher information matrix can be calculated analytically given any model, set of parameter values and a candidate design.  The inverse of the Fisher information matrix is an estimate of the covariance matrix.  The square roots of the diagonal elements of the covariance matrix are the expected standard errors of the parameter estimates (which is approximately equivalent to the output of NONMEM when asking for $COV).  The optimal design is therefore a design that optimises a summary measure of the Fisher information matrix.  Generally, in pharmacokinetics (PK) and pharmacokinetics-pharmacodynamics (PKPD) it is usual to use the determinant of the matrix as the summary measure.  Statistically it can be shown that maximizing the determinant of the Fisher information matrix is the same as minimizing the determinant of the covariance matrix.  Therefore the optimum design that maximises the determinant of the Fisher information matrix will provide the smallest standard errors.

 

Since PK and PKPD models in populations are hierarchical and nonlinear then the standard errors are approximate and asymptotic.  The design is also dependent on the values of the unknown parameters.  It is common therefore in optimal design to consider that the parameter values and perhaps also that the model have some level of uncertainty.

 

Choosing the optimal design

The concept optimal design engenders the notion that there is one single design which is better than all the rest and that this design is in some manner the most desirable design.  While it is true that the optimal design is the best single design it is not necessarily the most desirable design.  Indeed, there is, as yet, no proof that it is even possible to locate the optimal design for a nonlinear mixed effects model. 

 

Similarly searching for the optimal design is a non-trivial process.  The determinant of the Fisher information matrix over the design space is often a very convoluted surface and many standard (local) search strategies, e.g. gradient methods, can become easily trapped in local minima.  It is usual therefore to choose search methods which are considered to be global, or near global procedures.  Global means the algorithm given enough time would look almost everywhere in the design space.  Local means that the algorithm follows a specific path with little if any deviation to explore surrounding areas.  Two main advantages of global search algorithms are 1) that it doesn’t matter where you start (not dependent on initial estimates of the sampling times) and 2) they generally do not have convergence problems.  The principle disadvantage is their very slow convergence.  WinPOPT has two (near) global search algorithms.  These are simulated annealing (an asymptotically global search algorithm) and an exchange algorithm.  Simulated annealing can be very slow but provides good convergence properties.  The exchange algorithm is much quicker but has less rigorous convergence properties.  The exchange algorithm can be used for quick exploration of the design space, but the final design should be based on a run from simulated annealing.

 

It should be noted that simulated annealing and the exchange algorithm always converge.  However, if their search is terminated, then the sampling times provided at this juncture are not subject to bias.  This contrasts with modelling using maximum likelihood estimation methods (that most users will be familiar with) where a lack of convergence implies that the parameters are not unbiased best least squares estimators and therefore probably contain some element of bias.  So to recap, sampling times gained during a non-converged optimisation run using WinPOPT are just that, they are sampling times that provide the current estimates of the standard errors of the parameters.

 

 

Interpreting the design output

When exploring the design space using the determinant of the Fisher information matrix there are various metrics that are provided in the output.  Not all have practical significance.

 

 

Metric

What is it?

Practical use?

Determinant

The determinant is calculated mathematically from the Fisher information matrix and is a summary measure of the volume of a matrix.

 

The optimal design will be the one that provides the highest value of the determinant.  This metric has no practical utility and cannot be interpreted by itself.

Criterion

The determinant to the power of 1 over the number of parameters in the model.

The ratio of 2 criterion for any given model can be expressed as a fraction (or %) and is called efficiency; e.g. design A is 30% more efficient than design B.  A criterion by itself has no meaning.

 

Eigenvalues

These are the eigenvalues of the Fisher information matrix.

The ratio of the highest to lowest eigenvalue can be expressed as a condition number.  The condition number is a measure of the amenability of the matrix to mathematical computation.  (Smaller numbers are better).

 

Standard errors

These are the asymptotic estimates of the estimation standard deviation for a parameter.  They are the square roots of the diagonal elements of the inverse of the Fisher information matrix.

Values less than some predefined limit (e.g. < 20% for fixed effects parameters) can be used to judge the utility of a design

 

From a practical perspective, the goodness of any particular design should be judged on its ability to solve the problem at hand; i.e. provide precise estimates of the parameters of interest and is able to be performed clinically.  For comparison between designs then criterion values can be used as their ratio to provide an efficiency rating.

 

 

Preamble                                                                                          BACK TO TOP

 

POPT is a set of programs that have been written for MATLAB to help in the design of clinical studies in which PK or PKPD data is to be collected for the intention of developing models.  POPT has mostly been written and developed in MATLAB ver 6.x, but is backwardly and upwardly compatible with MATLAB ver 5.x and ver 7.x.  POPT does not require any specialised MATLAB toolboxes to be available.  The provision of toolboxes will not improve the performance of POPT code as it is currently written, but it is likely that some amendments could be made to increase performance. 

 

The programs will run on any operating system platform (e.g. Windows, Linux ...) under which MATLAB has been successfully installed.  We have tried both Windows and Linux without problems.

 

POPT is available free of charge and may be downloaded at

http://www.uq.edu.au/pharmacy/sduffull/popt.htm

POPT is subject to copyright.  Please refer to section xx for details of the copyright license.

 

 

Installation                                                                                       BACK TO TOP

 

There are no special instructions for installation of POPT.  POPT is made available in a zipped format which contains a directory of files "POPT Distribution copy" and a subdirectory "Example files".  These directories can be unzipped to any location on your computer.

 

It is recommended that you then copy the contents of the directory "POPT Distribution copy" to a new directory in your current work space for MATLAB in order to ensure that you do not inadvertently change any of the download files. From your Current Directory in MATLAB you can highlight and run the appropriate programs.

 


 

The POPT suite of programs                                                      BACK TO TOP

 

There are 21 objects in the main directory of POPT.

 

Programs for optimal design - DO NOT EDIT THESE

 

disp_2_scrn_sa.m

graphical interim results of simulated annealing algorithm

display_to_screen.m

graphical interim results of exchange algorithm

exch_setp.m

initialization of the exchange algorithm

exchange.m

exchange algorithm

exchange_with_v6_patch.m

this is the same as the exchange algorithm but is designed to work with version 6 of MATLAB (or any version that does not have access to the subroutine unidrnd()). If the usual version of the exchange algorithm does not work then rename exchange.m to exchange_v7.m and rename exchange_with_v6_patch.m to exchange.m

p_error.m

checks for some user errors in POPT_INI.m

p_expect.m

computes expectation of det(information matrix) over upper and lower bound

p_nars.m

searches for the max(det(information matrix)) by random sampling

p_output0.m

produces output to screen and file for single model optimizations

p_output1.m

produces output to screen and file for multiple model optimizations

p_sa_a.m

simulated annealing algorithm for approximate design

p_sa_dn.m

simulated annealing algorithm for optimizing dose number for samples

p_sa_m8.m

general simulated annealing algorithm

pop_run.m

sets up matrices for computation of the information matrix

pop_runs.m

evaluates information matrix and computes its determinant

popt.m

RUNS POPT (RUN THIS FILE)

 

Programs for setting up optimal design - EDIT THESE AS REQUIRED

 

P_MODEL.m

enter your structural model here

POPT_INI.m

enter your design specification and model details here

 

Programs for evaluating sampling windows - DO NOT EDIT THESE

 

p_robutst_1.m

RUNS sampling window estimation (RUN THIS FILE)

p_robust_popt.m

sets up sampling window estimation

repeat_part_3.m

resamples over sampling windows

 

Programs for evaluating sampling windows - EDIT AS REQUIRED

 

POPT_INI_robust.m

enter your sampling window and design specification and model details here

 

Miscellaneous

variable_names.txt

a list of variable names

P2NM

Creates a NONMEM data file from the current design in POPT_INI.m

 

 

To run POPT                                                                                    BACK TO TOP

 

To run POPT change MATLAB to the current directory that the POPT programs are located and either type POPT at the command window prompt or open POPT and click on the run button.

 

POPT will then run and produce output to a file which is echoed to the screen.  Outputs to files are appended - so no previous runs will be overwritten.  This means that all runs will be appended one after the other and the output file can get very large. The output file is generated automatically and its naming is based on the optimization method, e.g. method0.res is the file name when no optimization method is selected (i.e. method=0).

 

To evaluate sampling windows run P_ROBUST_1.

 

Input files for running POPT                                                       BACK TO TOP

 

            Example of POPT_INI file

            Example of POPT_INI_ROBUST file

 

POPT does not have an advanced or particularly user friendly method for entering your models and designs.  The structural model must be written in the file P_MODEL.m and the design component must be written in the file POPT_INI.m.  Because of this inefficiency it is recommended that you create a new directory (copying all POPT files) for each new design project.  The complete design is printed out with the output files and therefore you will always have a record of the design that you tried.

 

Examples of P_MODEL are provided for

  • 2 compartment model with first order input
  • 1 compartment model with zero order input
  • 1 compartment model with first order input defined by ODES
  • multiple models where a 1 compartment model and 2 compartment model are written
  • multiple response (PK and PD)

 

Each P_MODEL.m file contains a detailed description of how to construct the structural model for each of these occasions and further description is not provided here.

 

Examples of POPT_INI are provided for

  • multiple models where a 1 compartment model and 2 compartment model are written
  • multiple response (PK and PD)

 

POPT_INI is a relatively detailed (and therefore complex) input file.  A description of each item is provided.

 

 

MULTIPLE MODELS:  This situation arises when you want to optimize a design for more than 1 model at the same time.  This may be because

  1. you have 2 competing models and are not sure which is true,
  2. you have 2 drugs which you want to optimize the design for at the same time,
  3. you have 1 drug but 2 different sets of parameter values (e.g. for a paediatric and adult population)

 

MULTIPLE RESPONSES:  This situation arises when you have more than 1 response from a model.  Consider the framework that A à B à C, then C is dependent on B which is dependent on A.  The information matrix for C will therefore also have information about parameters that describe the formation of B from A as well as those that describe B à C.  Situations where multiple response models occur can include:

  1. parent, metabolite PK models
  2. plasma, urine PK models
  3. PKPD models

 

Example POPT_INI.m file.                                                   BACK TO TOP

 

There are 5 parts to this file.

Part 1 deals with model features

Part 2 design features

Part 3 optimization features

Part 4 optimization methods

Part 5 prior information

 

The description of the POPT_INI file will be divided into these sections.

 

MODEL FEATURES                                                            Back to this example

 

weighting=[];

Use this function only if you have either multiple models or multiple responses.  In this case you need to enter a weighting value for each candidate model.  The weightings need not add to 1, but this is convenient.  e.g. use: 

            weighting = [0.5 0.5];

 

mr = 0;

POPT cannot determine the difference between a multiple response model and multiple competing models.  Use this flag variable to indicate to POPT which it is.  If mr = 0 then this is a multiple model and if mr = 1 then it is a multiple response.

 

BETAI=[];

These are your prior values of the model parameters.  If you have more than one set of parameters (multiple models or multiple response models) then enter a row for each model.  You cannot have more than 2 rows for a multiple response model.  You can have any number of rows for a multiple model.  BETAI must be symmetric.  If there are an uneven number of parameters between model 1 and model 2 then you need to insert place holders, e.g. notation: 

            BETAI = [1 4 20 1 1

          1 2 20 4 5];

in this example the first row has two 1's as placeholders.  Do not use 0's as placeholders.

 

POP_PAR(1,:)=();

Enter the parameter names.  This variable is free text.  Example use

            POP_PAR(1,:)=('         Ka  CL   V');

If you have 2 models then you will need two rows of names.  The total length of the string vector must be the same for each, example use

            POP_PAR(1,:)=('         Ka  CL   V              ');

            POP_PAR(2,:)=('         Ka  CL   V   Emax   EC50');

here the spaces in POP_PAR(1,:) are placeholders to make the POP_PAR matrix square.

 

FXFI=[];

Use this vector or matrix to fix none, one or more fixed effects parameters (i.e. to remove them from consideration for optimization).  Enter the position number of the parameter in BETAI that you wish to fix.  For a single model an example use is

            FXFI=[2];

this will fix parameter 2.  If BETAI has more than 1 row then FXFI should also have more than 1 row. You can use 0's as placeholders, example use

            FXFI=[2

           0];

in this syntax position 2 of the first row of BETAI is fixed but no positions are fixed in the second row.

 

OI=[];

OI are the diagonal elements of the variance-covariance matrix of random effects.  Only the diagonal elements can be entered, example use

            OI=[0.25 0.1 0.1];

There needs to be an element in OI for every element in BETAI.  The value 0.00001 should be used as placeholders for values that are not intended to be estimated - and also FXI should indicate a fixed position.

 

FXI=[];

Use this vector or matrix to fix none, one or more variance of the between subject random effects parameters (i.e. to remove them from consideration for optimization).  Enter the position number of the parameter in OI that you wish to fix.  For a single model an example use is

            FXI=[1];

This fixes position 1 in the OI matrix.

 

BSVMOD=[];

This vector (or matrix if it is a multiple model) indicates the structural form of the between subject variance-covariance matrix.  A 0 indicates an additive (normal) model and a 1 an exponential (lognormal) model.  Example use

            BSVMODI=[1 1 1];

The dimensions of BSVMOD must equal the dimensions of OI which must equal the dimensions of BETAI.  You can use either a 0 or 1 for a placeholder for BSVMOD - they are simply ignored if not needed (which is determined by FXI).

 

sd_prior=[];

sd_addi=[];

These two terms are used to provide the prior values for the standard deviation of the residual variance.  Either or both additive and proportional error terms can be added.  It is always recommended to include an additive term, even without prior evidence, as this ensures that the optimization routine does not equate smaller values of the response variable with greater precision

in its estimation, which can lead to the erroneous situation where an infinitely small concentration is without error.  An example of use is

            sd_propi=[0.15];

     sd_addi=[0.1];

If you have multiple models then an example of use is

            sd_propi=[0.15

               0.2];

            sd_addi=[0.1

              0.1];

These variables need to be the same number of rows as BETAI.         

 

FXRI=[];

Use this vector or matrix to fix none, one or more variance of the residual random effects parameters (i.e. to remove them from consideration for optimization).  Enter the position number of the parameter that you wish to fix.  For a single model an example use is

            FXRI=[2];

This fixes sd_addi.  Whereas:

            FXRI=[1];

fixes sd_propi

For a multiple model the syntax is not as obvious but to fix the additive component from model 1 and the proportional component from model 2 you would write:

            FXRI=[2

           1];

 


 

DESIGN FEATURES                                                           Back to this example

 

ns =

maxs =

These variables are used with (method=6) and allow for an approximate design to be constructed.  ns is the total number of subjects in the study and maxs is the total number of samples allowed.  POPT will then determine the best arrangement of dividing the samples up amongst the patients.  Example syntax:

ns = 100;

maxs = 500;

The arrangement will be based on the allocation of samples per group provided by NT. This method optimizes the number of patients per group and the sampling times.

 

NUM_SUB=[];

This variable is used with all other methods (i.e. method <> 6) and replaces ns and maxs.  NUM_SUB is used to define the number of subjects in each group.  Example use:

            NUM_SUB=[100 100];

Where there are 100 patients allocated to group 1 and 100 patients allocated to group 2.  Each group has a different design.

 

NT=[]

NT defines the number of sampling times per group.  Example use:

            NT=[3 3];

Indicates that group 1 has 3 sampling times as does group 2.

 

D=[];

This signifies the dose to be used for each group. Example use:

D=[100 100];

If you have multiple models (or a multiple response model) then dose needs to be entered for each model.  This is especially important for multiple models where the models may represent different drugs.  Example use:

D=[100 100

   100 100];

 

dose_zero_fix=;

This is a special feature that allows the user to set the information matrix to a matrix of zeros when no dose is to be given. This is of importance when model responses are present in the absence of dosing (e.g. a baseline value).  If the value is 1 then a dose of zero yields an information matrix of zeros.  If the value is 0 then a dose of zero yields an information matrix according to whatever components of the model can be estimated for this dose level.  Example use

Dose_zero_fix=0;

 

 

DI=[];

DI is the dose interval.  There must be a dose interval for every element in the dose matrix (D).  Example use

            DI=[24 24];

Example use for multiple models:

            DI=[24 24

         24 24];

 

DPM=[];

DPM are the initial estimates of the design points.  In this case design points are sampling times.  This should be written as a row vector, which can be expressed as a matrix (although still read as a column vector) by the use of the  ”, see example, the syntax

DPM=[0.1 1 4 ...

     2 6 12];

Is equivalent to

DPM=[0.1 1 4 2 4 12];

DPM must contain the same number of elements as the sum of NT.  All elements of DPM must fall between the lower (LBA) and upper bounds (UBA) provided below.

 

DPE=[];

DPE are discrete potential values of the design variables that may be chosen for optimization. These are only used with method=9 and is optional.  Example use

            DPE=[0.01 0.5 1 2 3 4 5 6 7 8 9 10 11 12 18 24 …

0.01 0.5 1 2 3 4 5 6 7 8 9 10 11 12];

In this setting the exchange algorithm would select possible best sampling times from only those provided in the list above.  The row vector may be any length. Do not use 0 (zero) as it is a reserved variable in the exchange algorithm.  It is usual not to enter these values, but rather let the exchange algorithm select discrete values from within the sampling space itself. Note in this example that the second sampling time cannot exceed 12.

 

DPE_index=[];

This row vector provides the break points for DPE so that the exchange algorithm can determine which set of discrete sampling times should be used for each group.  Example use (based on DPE example above):

            DPE_index=[16 14];

This indicates that group 1 has a total of 16 times from which to choose but group 2 has only 14. The number of elements in DPE_index must equal the number of elements in NT.

 

FX_DPM=[];

FX_DPM provides the index for any sampling times in DPM that are fixed and are not permitted to be optimized.  The default is not to fix any sampling times, see example

FX_DPM=[];

Leaving the row vector blank is read by POPT as entering a row vector of 0’s the same length as DPM. To fix a single sampling time (e.g. the 12 hour sampling time from group 2 only) then a 1 is used in the same position in the row vector for FX_DPM as the position of 12 is in DPM.  The other non-fixed elements in DPM are given a 0.  Example use

FX_DPM=[0 0 0

   0 0 1];

 

DN=[];

The dose number for which the samples are to be taken is provided by DN.  For instance a DN value of 5 would indicate that the sample should be taken on the 5th dose, i.e. on the dose following 4 previous doses given at the specified dose interval.  DN is a row vector that should be of equal length as DPM.  The default value of DN is a row vector of 1’s, indicating that all samples are taken from the first dose.  Example use

            DN=[];

If samples are to be taken from anything other than the first dose then the dose number for each sample must be provided.  The length of DN must be the same as the length of DPM.  Example use

            DN=[4 4 4

    4 4 1];

In this example all samples are taken from the 4th dose, except the 3rd sampling time from group 2 which is taken on the 1st dose.

MX_DN=[];

The dose number that any sample is taken can be optimized using method=7.  In this setting MX_DN provides the upper bound for the dose number that doses can be taken from.  If it is allowed that samples are taken from any dose up to steady state then it is recommended to set the MX_DN to be the dose where approximately 5 half-lives have elapsed, otherwise the optimization algorithm will have difficulty distinguishing the difference between a dose after 10 half-lives vs a dose after 20 half-lives.  Example use

            MX_DN= [5 5 5

   1 1 1];

In this example, the first group of subjects can have their samples taken on any dose up to and including the 5th dose, while group 2 subjects will have all of their samples taken on the first dose.

 

UBA=[];

This is the upper bound on the optimized sampling times.  The default is the upper bound equals the dose interval.  Example use

            UBA=[];

Any elements in