# Calculating the power of covariates in population non-linear mixed effects models: the Monte Carlo Mapped Power approach

This post is based on the work of, among others, Camille Vong and the hands-on course about the MCMP given by Rob ter Heine and Elin Svensson. Read and cite the following publications when using an MCMP analysis in your next project:

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4791488/

### Power

You have probably heard people talk about the power of a study. In short, what we mean with statistical power in drug research is: what is the probability of finding a significant (e.g. p < 0.05) drug effect when that drug effect is really present.

We can derive this power of a study by knowing the effect size of the study, the variability and the sample size. For example, if the effect size is large, if there is a low degree of variability in the effect and if we have unlimited patients in our study, you can imagine that we will have a very high power. Hence, we are able to get a p < 0.05 from the resulting analysis of that study in almost all trials that we execute. In reality, a study design resulting in a power of 80% is already considered acceptable in most cases.

FYI, there are many resources available online on how to do a ‘standard’ power analysis, also in R:

https://www.statmethods.net/stats/power.html

### Power in population non-linear mixed effects models

But what if we are not analyzing the data from a study with a t-test or a linear mixed effects model? What happens when we start using our population non-linear mixed effects models? What if we are interested in differences in the PK parameters between populations?

For example, does this new drug formulation really change the bioavailability? Does this sub-population of CYP metabolizers really have a different clearance of the drug?

What is the probability of getting a significant drop in the objective function value (OFV) in NONMEM when we analyse the results of this new study and what sample size do we need? That is where the Monte Carlo Mapped Power (MCMP) method comes in!

### Objective

In this post, I will show you how to calculate the power of identifying a 20% inhibition in clearance between 2 populations using the MCMP method.

This method is implemented in Perl-speaks-NONMEM and a detailed user guide is available online:

https://uupharmacometrics.github.io/PsN/index.html

### NONMEM models

In the MCMP analysis, we need to code two different NONMEM models. An input (full) model, which includes the covariate effect and the population parameters of your (literature) PK model, and an alternative (reduced) model, without any covariate effect in the model.

Full model – run1.mod

```\$PROBLEM Full model with covariate

\$INPUT ID TIME CMT AMT DV MDV GROUP
\$DATA Dataset_sim_MCMP.csv

\$MODEL
COMP(DEPOT, DEFDOSE)		; Depot compartment
COMP(CENTRAL, DEFOBS)		; Central compartment

\$PK
TVCL = THETA(1)*(1-GROUP*THETA(4))
CL = TVCL*EXP(ETA(1))

TVV = THETA(2)
V  = TVV*EXP(ETA(2))

KA = THETA(3)*EXP(ETA(3))

K20= CL/V

; Scaling factor compartment 2
S2= V

\$DES
; Depot compartment
; Central compartment

\$ERROR
IPRED=F
Y=IPRED*(1+EPS(1))

\$THETA
(0,4.5)    ; 1, CL (L/h)
(0,31.1)	; 2, V (L)
(0,6.07)	; 3, KA (/h)
(0,0.2)	; 4, CL effect

\$OMEGA
0.06
0.06
0 FIX

\$SIGMA
0.01

\$EST PRINT=5 MAX=9999 METHOD=1 INTERACTION POSTHOC NOABORT
\$COV COMP PRINT=E

\$TABLE ID PRED IPRED CWRES DV TIME EVID NOAPPEND
\$TABLE ID CL V ETA1 ETA2 NOAPPEND NOPRINT ONEHEADER FILE=patab1
```

Reduced model – run1red.mod

```\$PROBLEM Reduced model without covariate

\$INPUT ID TIME CMT AMT DV MDV GROUP
\$DATA Dataset_sim_MCMP.csv

\$MODEL
COMP(DEPOT, DEFDOSE)		; Depot compartment
COMP(CENTRAL, DEFOBS)		; Central compartment

\$PK
TVCL = THETA(1)*(1-GROUP*THETA(4))
CL = TVCL*EXP(ETA(1))

TVV  = THETA(2)
V  = TVV*EXP(ETA(2))

KA = THETA(3)*EXP(ETA(3))

K20= CL/V

; Scaling factor compartment 2
S2	= V

\$DES
; Depot compartment
; Central compartment

\$ERROR
IPRED=F
Y=IPRED*(1+EPS(1))

\$THETA
(0,4)    ; 1, CL (L/h)
(0,30)	; 2, V (L)
(0,7)	; 3, KA (/h)
0 FIX 	; 4, CL effect

\$OMEGA
0.01
0.01
0 FIX

\$SIGMA
0.01

\$EST PRINT=5 MAX=9999 METHOD=1 INTERACTION POSTHOC NOABORT
\$COV COMP PRINT=E

\$TABLE ID PRED IPRED CWRES DV TIME EVID NOAPPEND
\$TABLE ID CL V ETA1 ETA2 NOAPPEND NOPRINT ONEHEADER FILE=patab1
```

As you can see in the full model, we assume that there is a 20% inhibitory effect on clearance in GROUP 1 of the population which we would like to identify in our new study. There is inter-individual variability on clearance and the central volume of distribution.

## Monte Carlo Mapped Power – MCMP

So what is the MCMP method? Lets first quote the introduction of the user guide:

The mcmp tool (Monte-Carlo Mapped Power method) is a tool for power computations. The method is based on the use of individual Objective Function Values (iOFV) and aims to provide a fast and accurate prediction of the power and sample size relationship without any need for adjustment of the significance criterion. The principle of the iOFV method is as follows:

1. A large dataset (e.g. 1000 individuals) is simulated with a full model and subsequently the full and reduced models are re-estimated with this data set.
2. iOFV’s are extracted from the NONMEM .phi-files, and for each subject the difference in iOFV between the full and reduced model is computed (delta iOFV )
3. Delta iOFV’s are sampled according to the design for which power is to be calculated and a starting sample size (N)
4. The delta iOFV sum for each sample is calculated
5. Steps 3 and 4 are repeated many times
6. The percentage of sums of the OFV greater than the significance criterion (e.g. 3.84 for one degree of freedom and p = 0.05) is taken as the power for sample size N,
7. Steps 3-6 are repeated with increasing N to provide the power at all sample sizes of interest.

As you can read at point 1, we only perform one simulation of a very large dataset and estimate the full and the reduced model on this simulated dataset.

### Dataset creation

First, we need to create a NM dataset that has the same setup as we want to investigate in our new study. The sampling times can be picked based on clinical feasibility or, preferably, be optimized with an optimal design approach.

In this scenario, we want to use a parallel study design, a grouping variable (GROUP) in a 1:1 ratio and a 24h sampling window. We create 1000 subjects per group for the MCMP.

The drug is dosed (AMT=100) orally in compartment 1.

```######################################################
##########
########## Create dataset for MCMP
########## M.J. van Esdonk
##########
library(dplyr)

## Number of subjects per group
n = 1000

sampling_times <- c(0.25,0.5,1,2,4,6,12,24)

dose <- 100

################################# Group 1
## Create observation dataframe
df_obs <- expand.grid(ID=1:n,TIME=sampling_times) %>%
mutate(CMT=2,AMT=0,DV=NA,MDV=0,GROUP=0)

## Create observation dataframe
df_dose <- data.frame(ID=1:n,TIME=0,AMT=dose,CMT=1,DV=NA,MDV=1,GROUP=0)

## Bind datasets together
df_1 <- df_obs %>%
bind_rows(df_dose)

################################# Group 2
## Create observation dataframe
df_obs <- expand.grid(ID=(n+1):(n*2),TIME=sampling_times) %>%
mutate(CMT=2,AMT=0,DV=NA,MDV=0,GROUP=1)

## Create observation dataframe
df_dose <- data.frame(ID=(n+1):(n*2),TIME=0,AMT=dose,CMT=1,DV=NA,MDV=1,GROUP=1)

## Bind datasets together
df_2 <- df_obs %>%
bind_rows(df_dose)

#################################################
## Bind group 1 and group 2 datasets together
df <- df_1 %>%
bind_rows(df_2) %>%
arrange(ID,TIME)

write.table(df,file='Dataset_sim_MCMP.csv',na='.',row.names = F,sep=',',quote = F)

```

### Execution and results

Lets execute the MCMP function. We include the -stratify_on argument in this command to pick the same number of individuals from each group. More information can be found in the MCMP user guide, also on how to include unbalanced designs in your investigation.

`mcmp -full_model=run1.mod -reduced_model=run1red.mod -dir=mcmp_run1 -threads=5 -rplots=2 -stratify_on=GROUP `

Since only 1 simulation is performed and 2 models are re-estimated, the MCMP method is a fast method for the calculation of the power! The run time for this analysis on 2 threads was approximately 13 minutes.

After completion of the MCMP analysis, a mcmp_results.csv file is created showing the power at a range of total sample sizes for different alpha’s. The R script that is generated by PsN, PsN_mcmp_plots.R, visualizes these results and can anwer our questions on the power of the study design. If no pdf file is created, manually execute the R script since some errors may occur during execution which terminates the automatic execution of the script after analysis.

Based on the study design we used, the power of our parallel design study with 8 subjects in each group would only be 47.5% at an alpha of 0.05. We would need to increase the study size to a total of 36 subjects (18 per cohort) to reach a power of 80%.

### Conclusion

In conclusion, the MCMP is a fast method to estimate the power of non-linear mixed effects models. This method is implemented in PsN and therefore free and easy to use.

In this post, the analysis was performed for only 1 study design. The simulation dataset and the reduced model can be changed to explore different scenario’s. One could also add new alternative models with certain timepoints ignored to explore that in one MCMP run.