04-Inputs.Rmd

# Pmetrics Input Files


## Data.csv Files

Pmetrics accepts input as a spreadsheet "matrix" format. It is designed
for input of multiple records in a concise way. **Please keep the number
of characters in the file name ≤ 8.**

In <span class="r6">R6</span> use `PM_data$new("filename")` to create a `PM_data` object by reading the file.

```{r echo=T, eval=FALSE}
#ensure that data.csv is in the working directory
data1 <- PM_data$new("data.csv")
```

You can also build an appropriate data frame in R and provide that as an argument to `PM_data$new()`.

```{r echo=T, eval=FALSE}
#assume df is data frame with at least these columns:
#id, time, dose, out
data1 <- PM_data$new(df)
```

Once you have created the `PM_data` object, you never need to create it again during your R session. You also don't have to bother copying the data file to the Runs folder each time you run the model.

In <span class="legacy">Legacy</span> you must always have the the data file in the current working directory. You can manually copy it there from a previous run or some other folder or use the shortcut of providing a prior run number as an argument to `NPrun` or `ITrun`.

```{r echo=T, eval=FALSE}
#Run 1 - ensure that data.csv is in the working directory
NPrun("data.csv", "model.txt")

#run 2 - use the data from run 1 in this run
#note that the file model.txt still has to be copied
# into the working directory in this example
NPrun(data = 1, "model.txt")
```

### Data file format

If you create your data as a file, which is the most
common scenario, the file should be in comma-separated-values (.csv) format. It is possible
to use other separators, like the semicolon, by setting the appropriate argument with `setPMoptions(sep = ";")`. 

Examples of programs
that can save .csv files are any text editor (e.g. TextEdit on Mac,
Notepad on Windows) or spreadsheet program (e.g. Excel). 

    

The file or data frame format is much more flexible in R6. The only required columns are id, time, dose, and out. You may also specify time as clock time if you include a date column. The default format of the date column is YYYY-MM-DD and HH:MM for the time column by default, but other formats can be specified. See `?PM_data` for more details. There is no header required, the column order can be anything you wish, but the names should be the same as in the legacy format below. Ultimately, `PM_data$new()` converts all valid data into the format used in legacy Pmetrics.

<span class="legacy">Legacy</span>

All .csv files must be formatted as below. As for R6, other separators are possible by using `setPMoptions()`. 

**IMPORTANT**: The order, capitalization and names of the header and the
first 12 columns are fixed. All entries must be numeric, with the
exception of ID and "." for non-required placeholder entries.

***POPDATA DEC\_11***
```{r echo=F, results='asis'}
library(knitr)
tab <- read.csv("Data/mdata.csv",na.strings=".")
names(tab)[1] <- "#ID"
tab$OUT <- as.character(tab$OUT)
kable(tab)
```

* ***POPDATA DEC\_11*** This is the fixed header for the file and must be
in the first line. It identifies the version. It is not the date of your
data file.

* ***\#ID*** This field must be preceded by the "\#" symbol to confirm
that this is the header row. It can be numeric or character and
identifies each individual. All rows must contain an ID, and all
records from one individual must be contiguous. Any subsequent row
that begins with "\#" will be ignored, which is helpful if you want to
exclude data from the analysis, but preserve the integrity of the
original dataset, or to add comment lines. IDs should be 11 characters
or less but may be any alphanumeric combination. **There can be at
most 800 subjects per run.**

* ***EVID*** This is the event ID field. It can be 0, 1, or 4. Every row
must have an entry.

    + 0 = observation

    + 1 = input (e.g. dose)

    + 2, 3 are currently unused

    + 4 = reset, where all compartment values are set to 0 and the time
counter is reset to 0. This is useful when an individual has multiple
sampling episodes that are widely spaced in time with no new
information gathered. This is a dose event, so dose information needs
to be complete.

* ***TIME*** This is the elapsed time in decimal hours since the first
event. It is not clock time (e.g. 09:30), although the `PMmatrixRelTime`
function can convert dates and clock times to decimal hours. 
Every row must have an entry, and within a given ID, rows
must be sorted chronologically, earliest to latest.

* ***DUR*** This is the duration of an infusion in hours. If EVID=1,
there must be an entry, otherwise it is ignored. For a bolus (i.e. an
oral dose), set the value equal to 0.

* ***DOSE*** This is the dose amount. If EVID=1, there must be an entry,
otherwise it is ignored.

* ***ADDL*** This specifies the number of additional doses to give at
interval II. It may be missing for dose events (EVID=1 or 4), in which
case it is assumed to be 0. It is ignored for observation (EVID=0)
events. Be sure to adjust the time entry for the subsequent row, if
necessary, to account for the extra doses. If set to -1, the dose is
assumed to be given under steady-state conditions. ADDL=-1 can only be
used for the first dose event for a given subject, or an EVID=4 event,
as you cannot suddenly be at steady state in the middle of dosing
record, unless all compartments/times are reset to 0 (as for an EVID=4
event). To clarify further, when ADDL=-1, all compartments in the
model will contain the predicted amounts of drug at the end of the
100th II interval.

* ***II*** This is the interdose interval and is only relevant if ADDL
is not equal to 0, in which case it cannot be missing. If ADDL=0 or is
missing, II is ignored.

* ***INPUT*** This defines which input (i.e. drug) the DOSE corresponds
to. Inputs are defined in the model file.

* ***OUT*** This is the observation, or output value. If EVID=0, there
must be an entry; if missing, this must be coded as -99. It will be
ignored for any other EVID and therefore can be ".". There can be at
most 150 observations for a given subject.

* ***OUTEQ*** This is the output equation number that corresponds to the
OUT value. Output equations are defined in the model file.

* ***C0, C1, C2, C3*** These are the coefficients for the assay error
polynomial for that observation. Each subject may have up to one set
of coefficients per output equation. If more than one set is detected
for a given subject and output equation, the last set will be used. If
there are no available coefficients, these cells may be left blank or
filled with "." as a placeholder.

* ***COV***\... Any column after the assay error coefficients is assumed
to be a covariate, one column per covariate. The first row for any subject
must have a value for all covariates, since the first row is always a dose.
Covariate values are applied at the time of doses.

### Manipulation of CSV files




* As we have seen, in <span class="r6">R6</span> `PM_data$new("filename")` will 
read an appropriate data file in the current working directory to create a 
new `PM_data` object. In <span class="legacy">Legacy</span>  
`PMreadMatrix("filename", ...)` reads *filename* and
creates a PMmatrix object in R which can be plotted (see
`?plot.PMmatrix`) or otherwise analyzed. However, unlike R6, it cannot be used
to run a model. For that, you need to copy the file into the working directory
each time, either yourself or by using the `NPrun(data = 1, ...)` shortcut, for 
example.

* <span class="r6">R6</span> `PM_data$write("filename")` will write the `PM_data` object to a file called
"filename". This can be useful if you have loaded or created a data file and then
changed it in R. <span class="legacy">Legacy</span> 
`PMwriteMatrix(data.frame, "filename", ...)` writes an
appropriate data frame as a new .csv file. It will first check the
data.frame for errors via the `PMcheck()` function below, and writing will
fail if errors are detected. This can be overridden with `override=T`.

* <span class="r6">R6</span> `PM_data$new()` automatically calls 
`PMcheck(filename / PMmatrix, model,...)`, which is also available as an
indepedent function in either <span class="r6">R6</span> or 
<span class="legacy">Legacy</span>. This function will check a .csv
file named *filename*, a `PM_data` object in <span class="r6">R6</span> or a
`PMmatrix` data frame in <span class="legacy">Legacy</span>  containing a previously
loaded .csv file (the output of `PMreadMatrix`) for errors which would
cause the analysis to fail. If a model file is provided, and the data
file has no errors, it will also check the model file for errors. If it
finds errors, it will generate a new errors.xlsx file with all errors
highlighted and commented so that you can find and correct them easily.
See `?PMcheck` for details in R.

* <span class="r6">R6</span> accepts data files with calendar dates and clock
times by calling the <span class="legacy">Legacy</span> function `PMmatrixRelTime()` 
any time `PM_data$new()` is invoked. This function converts dates and clock times of
specified formats into relative times for use in the NPAG, IT2B and
Simulator engines. This means Pmetrics in <span class="r6">R6</span> does not 
require any action by the user to handle calendar dates and clock times. In
contrast, for <span class="legacy">Legacy</span> `PMmatrixRelTime()` must
be called by the user and the output used to create a data frame with relative
times that can be saves as a new .csv file with `PMwriteMatrix()`, which in turn
serves as input to a run. See `?PMmatrixRelTime` for details.

The following functions are the same in either <span class="r6">R6</span> or 
<span class="legacy">Legacy</span>.

* `PMwrk2csv()` This function will convert old-style, single-drug USC\*PACK
.wrk formatted files into Pmetrics data .csv files. Details are
available with `?PMwrk2csv` in R.

* `PMmb2csv()` This function will convert USC\*PACK .mb files into Pmetrics
data .csv files. Details are available with `?PMmb2csv` in R.

* `NM2PM()` Although the structure of Pmetrics data files are similar to
NONMEM, there are some differences. This function attempts to
automatically convert to Pmetrics format. It has been tested on several
examples, but there are probably NONMEM files which will cause it to
crash. Running `PMcheck()` afterwards is a good idea. Details can be found
with `?NM2PM` in R.

### Specifying Data Objects in R6

<span class="r6">R6</span>

R6 introduces a new concept, the data object. The idea of this object is to represent a dataset that is going to be modeled/simulated.
All its behaviour is represented by the class `PM_data`. This class allows datasets to be checked, plotted, written to disk and more.

To create a PM_data object, users should have their data in the .csv format explained [here](#csv_format), and use the function `PM_data$new()` i.e:

```{r echo=T, eval=F}
drug_data <- PM_data$new("my_data.csv")
```

Where `my_data.csv` is a valid .csv data file compatible with Pmetrics.

## Specifying Models in R6

<span class="r6">R6</span>

In R6 Pmetrics, use the `PM_model` function to create models directly
in R. See `?PM_model` for help on this object class. Blocks in the legacy model.txt files which were delimited with the "#" character become listsin R6.

The R6 model components are:

* [PRImary](#priR6)
* [COVariate](#covR6)
* [SECondary](#secR6)
* [BOLus](#bolR6)
* [INItial conditions](#iniR6)
* [F (bioavailability)](#FaR6)
* [LAG time](#lagR6)
* [DIFferential equations](#diffR6)
* [OUTputs](#outR6)


### PRImary variables {#priR6}
Primary variables are the model parameters that are to be estimated by
Pmetrics or are designated as fixed parameters with user specified
values. It should be a list of variable names, one name to a line.
Variable names should be 11 characters or fewer. Some variable names are 
[reserved](#reserved) for use by Pmetrics and cannot be used as
primary variable names. **The number of primary variables must be
between 2 and 32, with at most 30 random or 20 fixed.**

Each variable can be specified by `range` or `msd`. The first defines the absolute search space for that parameter for NPAG/NPOD. For IT2B/RPEM, it defines the mean of the piror as the midpoint of the range, and the range covers 6 standard deviations, e.g. ±3 SD above and below the mean, or 99.7% of the piror distribution. `msd` is the companion function that specifies a mean and SD in IT2B and RPEM. For NPAG/NPOD, it will be converted in to a range in the reverse fashion as described for `range`. For both specifying functions, `gtz` is an argument to force the parameter value to be positive, i.e. `gtz=T`, which is the default. To allow negative parameters, set `gtz=F`.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = list(
    Ke = range(0,5),
    V = msd(100,20),
    eff = range(-2,2,gtz=F)
  )
))
```


### COVariates {#covR6}
Covariates are subject specific data, such as body weight, contained in
the data .csv file. The covariate names, which are the column names in
the data file, must be declared, even if not used in the model object. Once declared, they can be used in secondary variable
and differential equations. The order and names should be the same as in the data file.

Covariates are applied at each dose event. The first dose event for each
subject must have a value for every covariate in the data file. 

<span class="update">Update</span> By default, missing covariate values for subsequent dose events are linearly interpolated between existing values, or carried forward if the first value is the only non-missing entry.

To specify a new covariate value at a time other than a dose, enter a dose event in the data file with 0 dose amount and the new covariate value.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = list(...),
  cov = list("wt","age")
))
```

### SECondary variables {#secR6}
Secondary variables are those that are defined by equations that are
combinations of primary, covariates, and other secondary variables. If
using other secondary variables, define them first within this block.
Equation syntax must be Fortran. Specify each variable equation as a character vector. It is permissible to have conditional
statements, but because expressions in this block are translated into
variable declarations in Fortran, expressions other than of the form \"X = function(Y)\" must be on a new line, prefixed by \"&\" and contain only variables which have been previously defined in the Primary, Covariate, or Secondary blocks. 

In the example below, V0 is the primary parameter which will be estimated, but internally, the model uses V as V0*wt, unless age is >18, in which case weight is capped at 75 kg. It's the same for CL0. Note that the conditional statement is not named.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = list(
    CL0 = range(0,5),
    V0 = msd(10,3),
    eff = range(-2,2,gtz=F)
  ),
  cov = list("wt","age"),
  sec = list(
    V = "V0*wt",
    "&IF(age >18) V = V0 * 75",
    CL = "CL0 * wt"
  )
))
```

### BOLus inputs {#bolR6}

By default, inputs with DUR (duration) of 0 in the data .csv file are
\"delivered\" instantaneously to the model compartment equal to the
input number, i.e. input 1 goes to compartment 1, input 2 goes to
compartment 2, etc. This can be overridden with NBOLUS(input number) =
compartment number.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  bol = list("NBCOMP(1) = 2")
))
```
### INItial conditions {#iniR6}

By default, all model compartments have zero amounts at time 0. This can be changed by specifying the compartment amount as X(.) = expression, where \".\" is the compartment number. Primary and secondary variables and covariates may be used in the expression, as can conditional statements in Fortran code. An \"&\" continuation prefix is not necessary in this block for any statement, although if present, will be ignored.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = pri = list(
    Ke = range(0,5),
    V = msd(100,30),
    IC3 = range(0,1000)
  ),
  cov = list("wt","age","IC2"),
  ini = list(
    X2 = "IC2*V",
    X3 = "IC3"
  )
))
```

In the example above, IC is a covariate with the measured trough
concentration prior to an observed dose and IC3 is a fitted primary parameter specifying an initial amount in unobserved compartment 3.

In the first case, the initial condition for compartment 2 becomes the
value of the IC covariate (defined in `cov` list) multiplied by
the current estimate of V during each iteration. This is useful when a
subject has been taking a drug as an outpatient, and comes in to the lab for PK sampling, with measurement of a concentration immediately prior to a witnessed dose, which is in turn followed by more sampling. In this case, IC or any other covariate can be set to the initial measured concentration, and if V is the volume of compartment 2, the  initial condition (amount) in compartment 2 will now be set to the measured concentration of drug multiplied by the estimated volume for each iteration until convergence.

In the second case, the initial condition for compartment 3 becomes
another variable, IC3 defined in the `pri` list, to fit in the
model, given the observed data.

### FA (bioavailability) {#FaR6}

Specify the bioavailability term, if present. Use the form FA(.) =
expression, where \".\" is the input number. Primary and secondary
variables and covariates may be used in the expression, as can
conditional statements in Fortran code. An \"&\" continuation prefix is not necessary in this block for any statement, although if present, will be ignored.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = pri = list(
    Ke = range(0,5),
    V = msd(100,30),
    FA1 = range(0,1)
  ),
  fa = list(
    fa1 = "FA1"
  )
)
```

### LAG time {#lagR6}

Specify the lag term, if present, which is the delay after an absorbed
dose before observed concentrations. Use the form TLAG(.) = expression, where \".\" is the input number. Primary and secondary variables and covariates may be used in the expression, as can conditional statements in Fortran code. An \"&\" continuation prefix is not necessary in this block for any statement, although if present, will be ignored.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = pri = list(
    Ke = range(0,5),
    V = msd(100,30),
    lag1 = range(0,4)
  ),
  lag = list(tlag1 = "lag1")
  )
)
```

#### Differential equations

Specify a model in terms of ordinary differential equations, in Fortran format. XP(.) is the notation for dX(.)/dt, where \".\" is the
compartment number. X(.) is the amount in the compartment. **There can
be a maximum of 20 such equations.**

Specify equations as elements in a list, with XP(1) replaced by XP1, for example, to name the list, and the list value a character vector in Fortran.


```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = pri = list(
    Ka = range(0,5),
    Ke = range(0,5),
    V = msd(100,30),
    Kcp = range(0,5),
    Kpc = range(0,5)
  ),
  diff = list(
    xp1 = "-Ka * X(1)",
    xp2 = "RATEIV(1) + Ka * X(1) - (Ke + Kcp) * X(2) + Kpc * X(3)",
    xp3 = "Kcp * X(2) - Kpc * X(3)"
  )
))
```

RATEIV(1) is the notation to indicate an infusion of input 1 (typically drug 1). The duration of the infusion and total dose is defined in the data.csv file. **Up to 7 inputs are currently
allowed.** These can be used in the model file as RATEIV(1), RATEIV(2), etc. The compartments for receiving the inputs of oral (bolus) doses are defined in the `bol` list, but can be accessed by using the B(1), B(2), etc notation in equations.


### OUTputs {#outR6}

Output equations are in Fortran format. Outputs are of the form Y(.) =
expression, where \".\" is the output equation number. Primary and
secondary variables and covariates may be used in the expression, as can conditional statements in Fortran code. An \"&\" continuation prefix is not necessary in this block for any statement, although if present, will be ignored. **There can be a maximum of 6 outputs.** 

They are referred
to as Y(1), Y(2), etc. These equations may also define a model
explicitly as a function of primary and secondary variables and
covariates.

The `out` list is a series of nested lists.  The outer list defines all the outputs. The next level defines each output equation. Within the equation list is the equation and the error model. Within the error model for that output, is the last list comprising model and assay error specifications. 

* To name output equations, Y(1) is replaced by Y1
```{r echo=T, eval=F}
out = list(
  Y1 = list(...)
)
```
* The output equation is a character vector followed by an error list. 

```{r echo=T, eval=F}
out = list(
  Y1 = list(
    "X(1)/V",
    err = list(...)
  )
)
```
* The error model for an output equation has two elements. The first is the `model` error, which can be one of three functions: `proportional`, `additive`, or `combination`. The arguments to these functions are a number and optionally `fixed`, which defaults to `FALSE`. If `fixed` is `FALSE`, the number serves as the starting estimate for the model error. If `fixed` is `TRUE`, the number serves as the model error, with no estimation. Note that you can only fix $\lambda$ currently to zero.

The second element is the `assay` error model. It is a vector of 4 numbers that define a polynomial equation to permit calculation of the standard deviation of an observation, based on the noise of the asaay. The four terms estimate SD according to the folowing equation: $C0 + C1 * [obs] + C2 * [obs]^2 + C3 * [obs]^3$ and $[obs]$ is the observation. The values for the coefficients should ideally come from the analytic lab in the form of inter-run standard deviations or coefficients of variation at standard concentrations. You can use the Pmetrics function `makeErrorPoly` to choose the best set of coefficients that fit the data from the
laboratory. Alternatively, if you have no information about the assay,
you can use the Pmetrics ERRrun engine as an argument to the `run` function for `PM_fit` objects  (i.e., `$run(engine="err")`) to estimate the coefficients
from the data. Finally, you can use a generic set of coefficients. We
recommend that as a start, $C0$ be set to half of the lowest
concentration in the dataset and $C1$ be set to 0.1. $C2$ and $C3$ can
be 0.

The proportional model weights each observation by $1/(\gamma * SD)^2$, where $\gamma$ is either fixed or estimated. The additive model weights each observation by $1/(\lambda + SD)^2$, where $\lambda$ is either fixed or estimated. The combination model uses $1/((\gamma * SD)^2 + (\lambda + SD)^2)$.

In the proportional model, $\gamma$ is a scalar on assay SD. In general,
well-designed and executed studies and models with low mis-specification will have data with $\gamma$ values approaching 1. Values <1 suggest over inflated assay noise. Poor quality, noisy data will result in $\gamma$ of 5 or more. A good starting value for $\gamma$ is usually 5, and sometimes 10 if data are particularly complex or noisy.

In the additive model, $\lambda$ is additive to assay SD. In general,
well-designed and executed studies and models with low mis-specification will have data with $\lambda$ values approaching 0. Values of 0 may suggest over inflated assay noise. Poor quality, noisy data will result in $\lambda$ of $5*C0$ or more. A good starting value for $\lambda$ is usually $3*C0$. Note, that $C0$
should generally not be 0, as it represents machine noise (e.g. HPLC or
mass spectrometer) that is always present.


```{r echo=T, eval=F}
out = list(
  Y1 = list(
    "X(1)/V",
    err = list(
      model = proportional(1),
      assay = c(0.15, 0.1, 0, 0)
    )
  )
)
```


```{r echo=T, eval=F}
out = list(
  Y1 = list(
    "X(1)/V",
    err = list(
      model = additive(1, fixed = TRUE)
      assay = c(0.05, 0.1, 0, 0)
    )
)
```

More complete examples.

```{r echo=T, eval=F}
mod <- PM_model$new(list(
  pri = pri = list(
    Ke = range(0,5),
    V = msd(100,30),
  ),
  out = list(
    y1 = list(
      "X(1)/V",
      err = list(
        model = proportional(5),
        assay = c(0.05, 0.1, 0, 0)
      )
    )
  )
))

mod2 <- PM_model$new(list(
  pri = pri = list(
    kin = range(0,5),
    kout = range(0,5),
    tpd = range(0,5),
    V = msd(100,30),
  ),
  sec = list("RES = B(1) * KIN/(KIN-KOUT) * (EXP(-KOUT*TPD)-EXP(-KIN*TPD))"),
  out = list(
    y1 = list(
      "RES/V",
      err = list(
        model = combination(0.4,3) #additive, proportional
        assay = c(0.3, 0.15, 0, 0)
      )
    )
  )
))
```

This last example is known as the Bateman equation for a model with
linear absorption (KIN) into and elimination (KOUT) from a central
compartment, and a time post-dose (TPD) or lag time. Here B(1) is the
oral bolus dosing vector for drug 1, and V is the volume of the central
compartment.


## Specifying Models in Legacy

<span class="legacy">Legacy</span>

In legacy Pmetrics, models are text files that are ultimately translated into Fortran text files with a header version of TSMULT\... 

However, after Pmetrics version 0.30, we  adopted a
very simple user format that Pmetrics will use to generate the Fortran
code automatically for you. Version 0.4 additionally eliminates the
previously separate instruction file. A model library is available on
our website at
[http://www.lapk.org/pmetrics.php](http://www.lapk.org/pmetrics.php).

**Naming your model files.** The default model file name is "model.txt,"
but you can call them whatever you wish. However, **please keep the
number of characters in the model file name ≤ 8**. When you use a model
file in NPrun(), ITrun(), ERRrun(), or SIMrun(), Pmetrics will make a
Fortran model file of the same name, temporarily renaming your file. At
the end of the run, your original model file will be in the /inputs
subfolder of the run folder, and the generated Fortran model file will
be called "model.for" and moved to the /etc subfolder of the run folder.
If your model is called "mymodel.txt", then the Fortran file will be
"mymodel.for".

You can still use appropriate Fortran model files directly, but we
suggest you keep the .for extension for all Fortran files to avoid
confusion with the new format. If you use a .for file as your model, you
will have to specify its name explicitly in the NPrun(), ITrun,
ERRrun(), or SIMrun() command, since the default model name again is
"model.txt." If you use a .for file directly, it will be in the /inputs
subfolder of the run folder, not in /etc, since you did not use the
simpler template as your model file.

**Structure of model files.** The new model file is a text file with 11
blocks, each marked by \"\#\" followed by a header tag. Only details which are different than the R6 documentation are included below.

* [\#PRImary variables](#pri)
* [\#COVariates](#cov)
* [\#SECcondary variables](#sec)
* [\#BOLus inputs](#bol)
* [\#INItial conditions](#ini)
* [\#Fa (bioavailability)](#fa)
* [\#LAG time](#lag)
* [\#DIFferential equations](#dif)
* [\#OUTputs](#out)
* [\#ERRor](#err)
* [\#EXTra](#extra)

For each header, only the capital letters are required for recognition
by Pmetrics. The blocks can be in any order, and header names are
case-insensitive (i.e. the capitalization here is just to show which
letters are required). Fortran is also case-insensitive, so in variable
names and expressions case is ignored. Details of each block are next,
followed by a [complete example](#completeEx).

Important: Sometimes it is important to preserve spacing and formatting in Fortran code that you might insert into blocks, particularly the \#EXTRA block. If you wish to do this, insert [format] and [\/format] before and after any code that you wish to reproduce verbatim with spacing in the fortran model file.

Comments: You can insert comments into your model text file by starting a line with a capital "C" followed by a space. These lines will be removed/ignored in the final fortran code.

### Primary variables {#pri}

Primary variables are the model parameters that are to be estimated by
Pmetrics or are designated as fixed parameters with user specified
values. It should be a list of variable names, one name to a line.
Variable names should be 11 characters or fewer. Some variable names are [reserved](#reserved) for use by Pmetrics and cannot be used as
primary variable names. **The number of primary variables must be
between 2 and 32, with at most 30 random or 20 fixed.**

On each row, following the variable name, include the range for the
parameter that defines the search space. These ranges behave slightly
differently for NPAG, IT2B, and the simulator.

- For all engines, the format of the limits is *min, max*. A single
    value will indicate that the parameter is to be fixed but unknown in
    the population, i.e. the value is taken as the starting point for
    the optimization, but the final value will depend on the model and
    data and will be the same across the population. A single value
    followed by an "!" will indicate that this value is to be held
    constant (i.e. fixed and known) across the population, and not to be
    estimated.

- For **NPAG**, when *min/max* limits are specified, they are
    absolute, i.e. the algorithm will not search outside this range.

- For **IT2B**, the range defines the Bayesian prior distribution of
    the parameter values for cycle 1. For each parameter, the mean of
    the Bayesian prior distribution is taken as the middle of the range,
    and the standard deviation is *xsig*\*range (see [[IT2B
    runs]{.underline}](\l)). Adding a plus sign (+) to a line will
    prevent that parameter from being assigned negative values. NPAG and
    the simulator will ignore the pluses as the ranges are absolute for
    these engines. Note that prior to version 1.5.0, this used to be an
    exclamation point (!) but to be consistent throughout the model
    file, the exclamation point is now used when fixed values are
    desired.

- The **simulator** will ignore the ranges with the default value of
    NULL for the *limits* argument. If the simulator *limits* argument
    is set to NA, which will mean that these ranges will be used as the
    limits to truncate the simulation (see [[Simulator
    Runs]{.underline}](\l)).

Example:
<div class="script">
\#Pri\
KE, +0, 5\
V, 0.01, 100\
KA, 0, 5\
KCP, 5\
KPC, 0 , 5\
Tlag1, 0, 2\
IC3, 0, 10000\
FA1, 0.5!
</div>

### Covariates {#cov}

By default, missing covariate values for subsequent dose events are
linearly interpolated between existing values, or carried forward if the
first value is the only non-missing entry. To suppress interpolation and
carry forward the previous value in a piece-wise constant fashion,
include an exclamation point (!) in any declaration line.

**Note** that any covariate relationship to any parameter may be
described as the user wishes by mathematical equations and Fortran code,
allowing for exploration of complex, non-linear, time-dependent, and/or
conditional relationships. This is accomplished in the [\#Sec](#sec) block.

Example:

<div class="script">
\#Cov\
wt\
cyp\
IC!
</div>

where IC will be piece-wise constant and the other two will be linearly
interpolated for missing values.

### Secondary variables {#sec}

Equation syntax must be Fortran. It is permissible to have conditional
statements, but because expressions in this block are translated into
variable declarations in Fortran, expressions other than of the form \"X
= function(Y)\" must be prefixed by a \"&\" and contain only variables
which have been previously defined in the Primary, Covariate, or
Secondary blocks. Note that prior to version 1.5.0, the continuation
symbol was "+" before each line, but to avoid confusion with the use of
"+" in the Primary block for IT2B models, and to be consistent with
Fortran notation, the "&" is used henceforth.

Example:

<div class="script">
\#Sec\
CL = Ke \* V \* wt\*\*0.75\
& IF(cyp .GT. 1) CL = CL \* cyp
</div>

### Bolus inputs {#bol}

Example:

<div class="script">
\#Bol\
NBCOMP(1) = 2
</div>

### Initial conditions {#ini}

Example:

<div class="script">
\#Ini\
X(2) = IC\*V\
X(3) = IC3\
</div>

In the first case, the initial condition for compartment 2 becomes the
value of the IC covariate (defined in \#Covariate block) multiplied by
the current estimate of V during each iteration. This is useful when a
subject has been taking a drug as an outpatient, and comes in to the lab
for PK sampling, with measurement of a concentration immediately prior
to a witnessed dose, which is in turn followed by more sampling. In this
case, IC or any other covariate can be set to the initial measured
concentration, and if V is the volume of compartment 2, the initial
condition (amount) in compartment 2 will now be set to the measured
concentration of drug multiplied by the estimated volume for each
iteration until convergence.

In the second case, the initial condition for compartment 3 becomes
another variable, IC3 defined in the \#Primary block, to fit in the
model, given the observed data.

### Fa (bioavailability) {#Fa}

Specify the bioavailability term, if present. Use the form FA(.) =
expression, where \".\" is the input number. Primary and secondary
variables and covariates may be used in the expression, as can
conditional statements in Fortran code. An \"&\" continuation prefix is
not necessary in this block for any statement, although if present, will
be ignored.

Example:

<div class="script">
\#Fa\
FA(1) = FA1
</div>

### Lag time {#lag}

Example:

<div class="script">
\#Lag\
TLAG(1) = Tlag1
</div>

### Differential equations {#dif}

Example:

<div class="script">
\#Dif\
XP(1) = -KA\*X(1)\
XP(2) = RATEIV(1) + KA\*X(1) - (KE+KCP)\*X(2) + KPC\*X(3)\
XP(3) = KCP\*X(2) - KPC\*X(3)
</div>

### Outputs {#out}

Output equations, in Fortran format. Outputs are of the form Y(.) =
expression, where \".\" is the output equation number. Primary and
secondary variables and covariates may be used in the expression, as can
conditional statements in Fortran code. An \"&\" continuation prefix is
not necessary in this block for any statement, although if present, will
be ignored. **There can be a maximum of 6 outputs.** They are referred
to as Y(1), Y(2), etc. These equations may also define a model
explicitly as a function of primary and secondary variables and
covariates.

Examples:

<div class="script">
\#Out\
Y(1) = X(2)/V

\#OUT

RES = B(1) \* KIN/(KIN-KOUT) \* (EXP(-KOUT\*TPD)-EXP(-KIN\*TPD))

Y(1) = RES/VD
</div>

### Error {#err}

Unlike the R6, the error block is separate from the output block. In Legacy Pmetrics, this block contains all the information Pmetrics requires for the
structure of the error model. 

To specify the model in this block, the first line needs to be either
L=number or G=number for a $\lambda$ or $\gamma$ error model. The
number is the starting value for $\lambda$ or $\gamma$. If you include an exclamation point (!) in the declaration,
then $\lambda$ or $\gamma$ will be fixed and not estimated. Recall that you can
only fix $\lambda$ currently to zero.

The next line(s) contain the values for $C0$, $C1$, $C2$, and $C3$,
separated by commas. There should be one line of coefficients for each
output equation. By default Pmetrics will use values for these
coefficients found in the data file. If none are present or if the model
declaration line contains an exclamation point (!) the values here will
be used.

Example 1: estimated $\lambda$, starting at 0.4, one output, use data file
coefficients but if missing, use 0.1,0.1,0,0.

<div class="script">
\#Err\
L=0.4\
0.1,0.1,0,0\
</div>

Example 2: fixed $\gamma$ of 2, two outputs, use data file coefficients but
if missing, use 0.1,0.1,0,0 for the first output, but use 0.3, 0.1, 0, 0
for output 2 regardless of what is in the data file.

<div class="script">
\#Err\
G=2!\
0.1,0.1,0,0\
0.3,0.1,0,0!\
</div>

### Extra

This block is for advanced Fortran programmers only. 
<span class="update>Update</span>It is not yet implemented in R6.
Occasionally, for very complex models, additional Fortran subroutines are required. They can be placed here. The code must specify complete Fortran subroutines
which can be called from other blocks with appropriate call functions.
As stated earlier, sometimes it is important to preserve spacing and
formatting in Fortran code that you might insert into blocks,
particularly the \#EXTRA block. If you wish to do this, insert [format] and [\/format] in the fortran model file around the affected code.

### Reserved Names {#reserved}

The following cannot be used as primary, covariate, or secondary
variable names. They can be used in equations, however.

```{r echo=F}
res <- read.csv("Data/reserved.csv")
knitr::kable(res, col.names = gsub("[.]", " ", names(res)))
```


### Complete Example {#completeEx}

Here is a complete example of a model file, as of Pmetrics version 0.40
and higher:

<div class="script">
\#Pri\
KE, 0, 5\
V0, 0.1, 100\
KA, 0, 5\
Tlag1, 0, 3\

\#Cov\
wt\
C this weight is in kg\

\#Sec\
V = V0*wt\

\#Lag\
TLAG(1) = Tlag1\

\#Out\
Y(1) = X(2)/V\

\#Err\
L=0.4\
0.1,0.1,0,0\
</div>

*Notes*:

By omitting a \#Diffeq block with ODEs, Pmetrics understands that you
are specifying the model to be solved algebraically. In this case, at
least KE and V must be in the Primary or Secondary variables. KA, KCP,
and KPC are optional and specify absorption, and transfer to and from
the central to a peripheral compartment, respectively.

The comment line "C this weight is in kg" will be ignored.

### Brief Fortran Tutorial

Much more detailed help is available from
http://www.cs.mtu.edu/\~shene/COURSES/cs201/NOTES/fortran.html.

```{r echo=F}
for1 <- read.csv("Data/fortran1.csv")
knitr::kable(for1, col.names = gsub("[.]", " ", names(for1)))

```



------------------------- ----------------- -----------------------
**Relational Operator**   **Alternative**   **Meaning**
\<                        .LT.              less than
\<=                       .LE.              less than or equal
\>                        .GT.              greater than
\>=                       .GE.              greater than or equal
==                        .EQ.              equal
/=                        .NE.              not equal
------------------------- ----------------- -----------------------

+---------------------------------------+------------------------+
| **Selective Execution**               | **Example**            |
+---------------------------------------+------------------------+
| IF (logical-expression) one-statement | IF (T \>= 100) CL = 10 |
+---------------------------------------+------------------------+
| IF (logical-expression) THEN          | IF (T \>= 100) THEN    |
|                                       |                        |
| statements                            | CL = 10                |
|                                       |                        |
| END IF                                | V = 10                 |
|                                       |                        |
|                                       | END IF                 |
+---------------------------------------+------------------------+
| IF (logical-expression) THEN          | IF (T \>= 100) THEN    |
|                                       |                        |
| statements-1                          | CL = 10                |
|                                       |                        |
| ELSE                                  | ELSE                   |
|                                       |                        |
| statements-2                          | CL = CL                |
|                                       |                        |
| END IF                                | END IF                 |
+---------------------------------------+------------------------+


## Using legacy models in R6

<span class="r6">R6</span>
<span class="to-be-reviewed">TBR</span>

Pmetrics allows users to re-use old text-file based models withing the R6 framework by using the following syntax:

```{r echo=T, eval=F}
mod <- PM_model$new("model.txt")
```

Where model.txt is a model file present in your current working directory. The object obtained will a PM_model and have access to all methods of that class. see `?PM_model` for more information.