.. _ddmore0061_conv_example:


DDMoRe0061 Conversion Example
####################################


Overview
==========

A real world model, based on a gentamicin |pk| study [Harling2015]_. Combined |bolus| + |infusion| |pk| model. Uses ADVAN3 TRANS4 in |nonmem|. Some covariate effects on |pk| parameters. 210 individuals. 1949 total data rows. See :ddmore_link:`0061`

Data Conversion
================

We describe the conversion of this data set in |nonmem| format:-

.. code-block:: console

    c:\PoPy\validation\ddmore0061\Simulated_gentamicin_pk.csv
    
To this data set in |popy| format:-

.. code-block:: console

    c:\PoPy\validation\ddmore0061\popy\popy_data.csv

Using the :ref:`n2pdat_script` located here:-

.. code-block:: console

    c:\PoPy\validation\ddmore0061\popy\n2p_script.pyml
    
The original |nonmem| data set looks something like :numref:`table_nm_fields_ddmore0061`.

.. _table_nm_fields_ddmore0061:

.. csv-table:: Main data fields for |nonmem| (first six rows)
    :file: ddmore0061_nm_fields_first_6_rows.csv
    
The critical entries for the :ref:`n2pdat_script` are as follows:-

.. code-block:: pyml
    
   INPUT_NONMEM_FIELDS:
        obs_cmt_numbers: [1]
        dose_cmt_numbers: [1]
    
These two entries specify that the observations and doses all occur in compartment one. Technically the |nonmem| data set above has no CMT field. However CMT defaults to one in |nonmem| if it is NOT provided.

The dataset output by :ref:`n2pdat_script` is the same as the input, but with the 3 extra columns shown in :numref:`table_popy_fields_ddmore0061`.

.. _table_popy_fields_ddmore0061:

.. csv-table:: Extra data fields for |popy| (first six rows)
    :file: ddmore0061_popy_fields_first_6_rows.csv

Here the new 'DRUG_CONC' field will be the target value for the new |popy| script. The 'DRUG_CONC_FLAG' switches the observations on/off appropriately. Note that 'DRUG_CONC_FLAG' is the inverse of the |nonmem| 'MDV' flag. The new 'TYPE' field is formatted to process either |bolus| or |infusion| doses in |popy|.

Script Conversion
====================

We describe the manual conversion of the |nonmem| script file:-

.. code-block:: console

    c:\PoPy\validation\ddmore0061\Executable_gentamicin_pk.mod
    
To create the |popy| :ref:`fit_script` here:-

.. code-block:: console

    c:\PoPy\validation\ddmore0061\popy\fit_script.pyml

The |nonmem| script loads the data as follows:-

.. code-block:: none

    $DATA Simulated_gentamicin_pk.csv IGNORE=@
       
Whereas the |popy| script uses this syntax:-

.. code-block:: pyml

    FILE_PATHS: {input_data_file: popy_data.csv}
    
In the |nonmem| script, the theta/omega/sigma |fe| definitions are:-

.. code-block:: none

    $THETA  
        (0,4.1) ; CL
        (1,8.63) ; V1
        (0,1.21) ; Q
        (0,8.12) ; V2
        
        (-0.01,0.00982,0.016) ; BCLC
        (-0.409) ; ALB
        (0,0.00683) ; DCLC
        
    $OMEGA  
        0.0465 ; CL
        0.376 ; Q
        
    $SIGMA  
        0.0297 ; E1
        0.0486 ; E2
        

In |popy| these |fe| are defined in |effects|:-

.. code-block:: pyml

    EFFECTS:
        POP: |
            f[CL] ~ P 4.1
            f[V1] ~ unif(1, +inf) 8.63
            f[Q] ~ P 1.21
            f[V2] ~ P 8.12
            f[BCLC_eff] ~ unif(-0.01, 0.016) 0.00982
            f[ALB_eff] = -0.409
            f[DCLC_eff] ~ P 0.00683
            
            f[CL_isv, Q_isv] ~ diag_matrix() [ [ 0.0465, 0.376 ] ]
            
            f[PNOISE_var] ~ P 0.0297
            f[ANOISE_var] ~ P 0.0486
                
Note the ranges and initial values of the |fx| variables are defined using different syntax from the $THETA variables. Note also that the |fx| have explicit names, whereas the |nonmem| list of thetas has to be clumsily annotated with comments to make it human readable.

In |popy| you can specify a :ref:`diagonal matrix<pos_diag_mat_dist>`, which can be explicitly used as a covariance matrix input to a :ref:`mnorm_dist` distributed |rx| vector, as follows:-

.. code-block:: pyml

    EFFECTS:
        ID: |
            r[CL_isv, Q_isv] ~ mnorm([0,0], f[CL_isv, Q_isv])
            
There :pyml:`r[CL_isv]` is the equivalent of 'ETA(1)' and :pyml:`r[Q_isv]` is the equivalent of 'ETA(2)'. In |nonmem| you just have to remember the index values, whereas |popy| uses actual real life variable names.

In |popy| the sigmas are |fes| and also defined in |effects|.
                 
The |nonmem| PK section is as follows:-
       
.. code-block:: none
  
    $PK       
        BSA = 71.84*(WT**0.425)*(HT**0.725)
        BSA = BSA/10000
        HT1 = HT/100
        BMI = WT/HT1**2
        
        IF(NEWIND.NE.2) BCLC = CLC
        DCLC = CLC-BCLC
        IF(NEWIND.NE.2) BBSA = BSA
        
        TVCL = THETA(1) * (1 + THETA(5)*(BCLC-81) + THETA(7)*DCLC)
        TVV1 = THETA(2) * BBSA*(ALB/34)**THETA(6)
        TVQ  = THETA(3)
        TVV2 = THETA(4)
        
        CL = TVCL * EXP(ETA(1))
        V1 = TVV1
        Q  = TVQ  * EXP(ETA(2))
        V2 = TVV2
        
        S1 = V1
    
This maps almost exactly to the |popy| |model_params| section:-

.. code-block:: pyml

    MODEL_PARAMS: |
        BSA = 71.84 * (c[WT]**0.425) * (c[HT]**0.725)
        BSA = BSA / 10000
        HT1 = c[HT] / 100
        BMI = c[WT] / HT1**2
   
        BCLC = c[CLC]
        DCLC = c[CLC] - BCLC
        BBSA = BSA
        
        TVCL = f[CL] * (1 + f[BCLC_eff]*(BCLC-81) + f[DCLC_eff]*DCLC)
        TVV1 = f[V1] * BBSA * (c[ALB]/34)**f[ALB_eff]
        TVQ  = f[Q]
        TVV2 = f[V2]
        
        m[CL] = TVCL * EXP(r[CL_isv])
        m[V1] = TVV1
        m[Q]  = TVQ  * EXP(r[Q_isv])
        m[V2] = TVV2
        
        m[PNOISE_var] = f[PNOISE_var]
        m[ANOISE_var] = f[ANOISE_var]
    
Note in the above |nonmem| specifies 'S1 = V1', to scale the amounts in compartment one. This is not necessary in |popy|, as the compartment amounts are converted to concentrations in the |predictions| section explicitly.

Also note that |popy| requires the lines:-

.. code-block:: pyml

    m[PNOISE_var] = f[PNOISE_var]
    m[ANOISE_var] = f[ANOISE_var]
    
To propagate the |fx| noise |fes| which are similar to the |nonmem| builtin sigma variables. 
   
Note in this example the |nonmem| 'NEWIND' keyword is superfluous, as the :pyml:`c[CLC]` covariate and BSA variable are constant within each individual.

The |nonmem| control file prescribes the compartment model on this line:-

.. code-block:: none

    $SUBROUTINE ADVAN3 TRANS4
    
That uses the inbuilt 'ADVAN3' model with 'TRANS4' parametrisations. This uses the analytic solution for a two compartment model that expects the CL, V1, Q, V2 parameters to be defined in the |nonmem| PK section.
    
The |popy| control file here specifies the compartment model equations explicitly in the |derivatives| section:-

.. code-block:: pyml

    DERIVATIVES: |
        dose[DOSE_bolus] = @bolus{amt:c[AMT], lag:0.0}
        dose[DOSE_infrate] = @inf_rate{amt:c[AMT], lag:0.0, rate:c[RATE]}
        d[CENTRAL] = dose[DOSE_bolus] + dose[DOSE_infrate] - s[CENTRAL]*m[Q]/m[V1] + s[PERI]*m[Q]/m[V2] - s[CENTRAL]*m[CL]/m[V1]
        d[PERI]    =                                         s[CENTRAL]*m[Q]/m[V1] - s[PERI]*m[Q]/m[V2]

        
Note that |popy| has an equivalent builtin version of 'ADVAN3 TRANS4' as follows:-

.. code-block:: pyml

    DERIVATIVES: |
        s[CENTRAL, PERI] = @iv_two_cmp_cl{ dose: @bolus{amt:c[AMT]}, CL: m[CL], V1: m[V1], Q: m[Q], V2: m[V2] }
       
However this function is unable to process two different types of dose (bolus + infusion), unlike the long hand version above. Note that |popy| will be using an |ode| solver here (instead of a built in analytic solution), so we also need to specify this |popy| field:-

.. code-block:: pyml

    ODE_SOLVER: {SCIPY_ODEINT: {}}
    
Which tells |popy| to use the scipy 'odeint' solver, which is based on LSODA (like |nonmem| Advan13). 
       
The |nonmem| ERROR section is as follows:-
       
.. code-block:: none

    $ERROR
        Y = F*EXP(ERR(1)) + ERR(2)
        IPRED=F
        IRES=DV-IPRED

The |popy| equivalent is shown below in the |predictions| section:-
        
.. code-block:: pyml

    PREDICTIONS: |
        p[DV_CENTRAL] = s[CENTRAL]/m[V1]
        var = m[PNOISE_var] * p[DV_CENTRAL]**2 + m[ANOISE_var]
        c[DRUG_CONC] ~ norm(p[DV_CENTRAL], var)
        
Note here the line 'Y = F*EXP(ERR(1)) + ERR(2)' defines an proportional exponential normal noise model with an additional additive normal noise component. It is not clear (to the |popy| developers) what the distribution of a log normal + a normal distribution is. |popy| requires a known distribution to calculate a likelihood. If you approximate the exponential (assuming ERR(1) is small) and get 'Y = F*(1+ ERR(1)) + ERR(2)'. Then this is the standard proportional + additive noise model. In this case, the |popy| :ref:`fit_script` return a very similar objective function to the |ddmore| |nonmem| objective function. Indicating that |nonmem| also makes the same approximation. 
              
The |nonmem| control file specifies using FOCE fitting (METHOD=1) below:-
        
.. code-block:: none
        
   $ESTIMATION MAXEVALS=0 SIG=3 PRINT=10 METHOD=1 INTER
   
The |popy| equivalent, specifying the JOE fitting method (roughly equivalent to FOCE) is:-

.. code-block:: pyml

    FIT_METHODS: [JOE: {max_n_main_iterations: 0}]
      
Setting 'max_n_main_iterations' to zero, means that |popy| will optimise the |rx| parameters and leave the |fx| unchanged.