ed_laplace_ofv: Evaluate the expectation of determinant the Fisher Information Matrix (FIM) using the Laplace approximation.

Description

Compute the expectation of the det(FIM) using the Laplace approximation to the expectation. Computations are made based on the model, parameters, distributions of parameter uncertainty, design and methods defined in the PopED database or as arguments to the funciton.

Usage

ed_laplace_ofv(model_switch, groupsize, ni, xtopto, xopto, aopto, bpopdescr,
  ddescr, covd, sigma, docc, poped.db, method = 1, return_gradient = FALSE,
  optxt = poped.db$optsw[2], opta = poped.db$optsw[4], x = c())

Arguments

The design parameters to compute the gradient on.

xtopto

the sampling times

xopto

the discrete design variables

optxt

If sampling times are optimized

opta

If continuous design variables are optimized

aopto

the continuous design variables

method

If 0 then use an optimization routine translated from poped code written in MATLAB to optimize the parameters in the Laplace approximation. If 1 then use optim to compute both k and the hessian of k (see Dod

return_gradient

Should the gradient be returned.

model_switch

A matrix that is the same size as xt, specifying which model each sample belongs to.

groupsize

A vector of the numer of individuals in each group.

A vector of the number of samples in each group.

poped.db

A PopED database.

covd

Matrix defining the covariances of the IIV variances. Set to zero if not defined.

sigma

Matrix defining the variances can covariances of the residual variability terms of the model. can also just supply the diagnonal parameter values (variances) as a c().

docc

Matrix defining the IOV, the IOV variances and the IOV distribution

bpopdescr

Matrix defining the fixed effects, per row (row number = parameter_number) we should have:

column 1 the type of the distribution for E-family designs (0 = Fixed, 1 = Normal, 2 = Uniform, 3 = User Defined Distribution, 4 = lognormal and 5 = trunca

ddescr

Matrix defining the diagnonals of the IIV (same logic as for the bpopdescr).

Value

The FIM and the hessian of the FIM.

Details

This computation follows the method outlined in Dodds et al, "Robust Population Pharmacokinetic Experiment Design" JPP, 2005, equation 16.

Typically this function will not be run by the user. Instead use evaluate.e.ofv.fim.

Examples

Run this code

## Warfarin example from software comparison in:
## Nyberg et al., "Methods and software tools for design evaluation 
##   for population pharmacokinetics-pharmacodynamics studies", 
##   Br. J. Clin. Pharm., 2014. 

## Optimization using an additive + proportional reidual error to 
##   avoid sample times at very low concentrations (time 0 or very late samoples).
library(PopED)

## find the parameters that are needed to define from the structural model
ff.PK.1.comp.oral.sd.CL

## -- parameter definition function 
## -- names match parameters in function ff
sfg <- function(x,a,bpop,b,bocc){
  parameters=c(CL=bpop[1]*exp(b[1]),
               V=bpop[2]*exp(b[2]),
               KA=bpop[3]*exp(b[3]),
               Favail=bpop[4],
               DOSE=a[1])
  return(parameters) 
}

######################
# Normal distribution
######################
bpop_vals <- c(CL=0.15, V=8, KA=1.0, Favail=1)
bpop_vals_ed_n <- cbind(ones(length(bpop_vals),1)*1, # normal distribution
                        bpop_vals,
                        ones(length(bpop_vals),1)*(bpop_vals*0.1)^2) # 10\% of bpop value
bpop_vals_ed_n["Favail",]  <- c(0,1,0)
bpop_vals_ed_n

## -- Define initial design  and design space
poped.db.n <- create.poped.database(ff_file="ff.PK.1.comp.oral.sd.CL",
                                    fg_file="sfg",
                                    fError_file="feps.add.prop",
                                    bpop=bpop_vals_ed_n, 
                                    notfixed_bpop=c(1,1,1,0),
                                    d=c(CL=0.07, V=0.02, KA=0.6), 
                                    sigma=c(0.01,0.25),
                                    groupsize=32,
                                    xt=c( 0.5,1,2,6,24,36,72,120),
                                    minxt=0,
                                    maxxt=120,
                                    a=70,
                                    mina=0,
                                    maxa=100)


## ED evaluate using LaPlace approximation 
tic()
output <- evaluate.e.ofv.fim(poped.db.n,use_laplace=TRUE)
toc()
output$E_ofv

## expected value (roughly)
  tic()
  e_ofv_mc_n <- evaluate.e.ofv.fim(poped.db.n,ED_samp_size=500)
  toc()
  e_ofv_mc_n$E_ofv
  
  
  ## Using ed_laplce_ofv directly
  ed_laplace_ofv(model_switch=poped.db.n$global_model_switch,
                 groupsize=poped.db.n$groupsize,
                 ni=poped.db.n$gni,
                 xtopto=poped.db.n$gxt,
                 xopto=poped.db.n$gx,
                 aopto=poped.db.n$ga,
                 bpopdescr=poped.db.n$gbpop,
                 ddescr=poped.db.n$gd,
                 covd=poped.db.n$covd,
                 sigma=poped.db.n$sigma,
                 docc=poped.db.n$docc, 
                 poped.db.n)
  
  
  ######################
  # Log-normal distribution
  ######################
  
  # Adding 10% log-normal Uncertainty to fixed effects (not Favail)
  bpop_vals <- c(CL=0.15, V=8, KA=1.0, Favail=1)
  bpop_vals_ed_ln <- cbind(ones(length(bpop_vals),1)*4, # log-normal distribution
                           bpop_vals,
                           ones(length(bpop_vals),1)*(bpop_vals*0.1)^2) # 10% of bpop value
  bpop_vals_ed_ln["Favail",]  <- c(0,1,0)
  bpop_vals_ed_ln
  
  ## -- Define initial design  and design space
  poped.db.ln <- create.poped.database(ff_file="ff.PK.1.comp.oral.sd.CL",
                                       fg_file="sfg",
                                       fError_file="feps.add.prop",
                                       bpop=bpop_vals_ed_ln, 
                                       notfixed_bpop=c(1,1,1,0),
                                       d=c(CL=0.07, V=0.02, KA=0.6), 
                                       sigma=c(0.01,0.25),
                                       groupsize=32,
                                       xt=c( 0.5,1,2,6,24,36,72,120),
                                       minxt=0,
                                       maxxt=120,
                                       a=70,
                                       mina=0,
                                       maxa=100)
  
  
  
  ## ED evaluate using LaPlace approximation 
  tic()
  output <- evaluate.e.ofv.fim(poped.db.ln,use_laplace=TRUE)
  toc()
  output$E_ofv
  
  ## expected value (roughly)
  tic()
  e_ofv_mc_ln <- evaluate.e.ofv.fim(poped.db.ln,ED_samp_size=500)[["E_ofv"]]
  toc()
  e_ofv_mc_ln
  
  ## Using ed_laplce_ofv directly
  ed_laplace_ofv(model_switch=poped.db.ln$global_model_switch,
                 groupsize=poped.db.ln$groupsize,
                 ni=poped.db.ln$gni,
                 xtopto=poped.db.ln$gxt,
                 xopto=poped.db.ln$gx,
                 aopto=poped.db.ln$ga,
                 bpopdescr=poped.db.ln$gbpop,
                 ddescr=poped.db.ln$gd,
                 covd=poped.db.ln$covd,
                 sigma=poped.db.ln$sigma,
                 docc=poped.db.ln$docc, 
                 poped.db.ln)

Run the code above in your browser using DataLab

Get 50% off unlimited learning