Render Models in LaTeX

Author

Vijay Ivaturi

1 Introduction

When building pharmacometric or systems pharmacology models, it is often useful to share the mathematical form of the model, either for publication, collaboration, or documentation. However, transcribing equations by hand is error-prone and time-consuming.

The Pumas model macros capture enough structured information so that Latexify can automatically generate LaTeX (or MathJax) representations of a model’s parameters and equations. The equation can then be copied, adapted, or directly rendered in any environment.

This document walks through the Warfarin PKPD model below, showing how to:

  • Extract individual blocks (like the @param or @dynamics blocks).
  • Configure the look-and-feel of the generated LaTeX.
  • Stitch multiple blocks together into a single model representation.

2 Defining the Warfarin PKPD Model

First, the Warfarin model is restated. Ensure that Pumas is loaded, and then define the model:

using Pumas

warfarin_pkpd_model = @model begin
    @param begin
        # PK parameters
        """
        Clearance (L/h/70kg)
        """
        pop_CL  RealDomain(lower = 0.0, init = 0.134)

        """
        Central Volume L/70kg
        """
        pop_V  RealDomain(lower = 0.0, init = 8.11)

        """
        Absorption time (h)
        """
        pop_tabs  RealDomain(lower = 0.0, init = 0.523)

        """
        Lag time (h)
        """
        pop_lag  RealDomain(lower = 0.0, init = 0.1)

        # PD parameters
        """
        Baseline
        """
        pop_e0  RealDomain(lower = 0.0, init = 100.0)

        """
        Emax
        """
        pop_emax  RealDomain(init = -1.0)

        """
        EC50
        """
        pop_c50  RealDomain(lower = 0.0, init = 1.0)

        """
        Turnover
        """
        pop_tover  RealDomain(lower = 0.0, init = 14.0)

        # Inter-individual variability
        """
          - ΩCL
          - ΩVc
          - ΩTabs
        """
        pk_Ω  PDiagDomain([0.01, 0.01, 0.01])

        """
        Ωlag
        """
        lag_ω  RealDomain(lower = 0.0, init = 0.1)

        """
          - Ωe0
          - Ωemax
          - Ωec50
          - Ωturn
        """
        pd_Ω  PDiagDomain([0.01, 0.01, 0.01, 0.01])

        # Residual variability
        """
        Proportional residual error for drug concentration
        """
        σ_prop  RealDomain(lower = 0.0, init = 0.00752)

        """
        Additive residual error for drug concentration (mg/L)
        """
        σ_add  RealDomain(lower = 0.0, init = 0.0661)

        """
        Additive error for PCA
        """
        σ_fx  RealDomain(lower = 0.0, init = 0.01)
    end

    @random begin
        pk_η ~ MvNormal(pk_Ω)
        lag_η ~ Normal(0.0, lag_ω)
        pd_η ~ MvNormal(pd_Ω)
    end

    @covariates begin
        FSZV
        FSZCL
    end

    @pre begin
        # PK
        CL = FSZCL * pop_CL * exp(pk_η[1])
        Vc = FSZV * pop_V * exp(pk_η[2])
        tabs = pop_tabs * exp(pk_η[3])
        Ka = log(2) / tabs

        # PD
        e0 = pop_e0 * exp(pd_η[1])
        emax = pop_emax * exp(pd_η[2])
        c50 = pop_c50 * exp(pd_η[3])
        tover = pop_tover * exp(pd_η[4])
        kout = log(2) / tover
        rin = e0 * kout
    end

    @dosecontrol begin
        lags = (Depot = pop_lag * exp(lag_η),)
    end

    @init begin
        Turnover = e0
    end

    @vars begin
        cp := Central / Vc
        ratein := Ka * Depot
        pd := 1 + emax * cp / (c50 + cp)
    end

    @dynamics begin
        Depot' = -ratein
        Central' = ratein - CL * cp
        Turnover' = rin * pd - kout * Turnover
    end

    @derived begin
        """
        Warfarin Concentration (mg/L)
        """
        conc ~ @. Normal(cp, sqrt((σ_prop * cp)^2 + σ_add^2))

        """
        PCA
        """
        pca ~ @. Normal(Turnover, σ_fx)
    end
end
PumasModel
  Parameters: pop_CL, pop_V, pop_tabs, pop_lag, pop_e0, pop_emax, pop_c50, pop_tover, pk_Ω, lag_ω, pd_Ω, σ_prop, σ_add, σ_fx
  Random effects: pk_η, lag_η, pd_η
  Covariates: FSZV, FSZCL
  Dynamical system variables: Depot, Central, Turnover
  Dynamical system type: Nonlinear ODE
  Derived: conc, pca
  Observed: conc, pca

3 Generating LaTeX with Latexify

3.1 Loading Latexify

Latexify comes with Pumas so it is not necessary to install it separately. It can be loaded with:

using Pumas.Latexify

3.2 Latexifying Individual Blocks

Each Pumas model is composed of blocks that capture different stages of the model definition (@param, @random, @pre, @init, @dynamics, @derived, etc.). LaTeX can be generated for each block independently.

For instance, LaTeX for the parameter block can be generated as follows:

latex_param = latexify(warfarin_pkpd_model, :param)

popCLRealDomain(init=0.134;lower=0.0)popVRealDomain(init=8.11;lower=0.0)poptabsRealDomain(init=0.523;lower=0.0)poplagRealDomain(init=0.1;lower=0.0)pope0RealDomain(init=100.0;lower=0.0)popemaxRealDomain(;init=1.0)popc50RealDomain(init=1.0;lower=0.0)poptoverRealDomain(init=14.0;lower=0.0)pkΩPDiagDomain([0.010.010.01])lagωRealDomain(init=0.1;lower=0.0)pdΩPDiagDomain([0.010.010.010.01])σpropRealDomain(init=0.00752;lower=0.0)σaddRealDomain(init=0.0661;lower=0.0)σfxRealDomain(init=0.01;lower=0.0)

The above output is generated as the param block is latexified as seen below:

\begin{align}
pop_{CL} &\in \mathrm{RealDomain}\left( init = 0.134; lower = 0.0 \right) \\
pop_{V} &\in \mathrm{RealDomain}\left( init = 8.11; lower = 0.0 \right) \\
pop_{tabs} &\in \mathrm{RealDomain}\left( init = 0.523; lower = 0.0 \right) \\
pop_{lag} &\in \mathrm{RealDomain}\left( init = 0.1; lower = 0.0 \right) \\
pop_{e0} &\in \mathrm{RealDomain}\left( init = 100.0; lower = 0.0 \right) \\
pop_{emax} &\in \mathrm{RealDomain}\left( ; init = -1.0 \right) \\
pop_{c50} &\in \mathrm{RealDomain}\left( init = 1.0; lower = 0.0 \right) \\
pop_{tover} &\in \mathrm{RealDomain}\left( init = 14.0; lower = 0.0 \right) \\
pk_{\Omega} &\in \mathrm{PDiagDomain}\left( \left[
\begin{array}{c}
0.01 \\
0.01 \\
0.01 \\
\end{array}
\right] \right) \\
lag_{\omega} &\in \mathrm{RealDomain}\left( init = 0.1; lower = 0.0 \right) \\
pd_{\Omega} &\in \mathrm{PDiagDomain}\left( \left[
\begin{array}{c}
0.01 \\
0.01 \\
0.01 \\
0.01 \\
\end{array}
\right] \right) \\
\sigma_{prop} &\in \mathrm{RealDomain}\left( init = 0.00752; lower = 0.0 \right) \\
\sigma_{add} &\in \mathrm{RealDomain}\left( init = 0.0661; lower = 0.0 \right) \\
\sigma_{fx} &\in \mathrm{RealDomain}\left( init = 0.01; lower = 0.0 \right)
\end{align}

Similarly, the random effects can be examined:

latex_random = latexify(warfarin_pkpd_model, :random)

pkηMvNormal(pkΩ)lagηNormal(0.0,lagω)pdηMvNormal(pdΩ)

The above output is generated as the random block is latexified as seen below:

\begin{align}
pk_{\eta} &\sim \mathrm{MvNormal}\left( pk_{\Omega} \right) \\
lag_{\eta} &\sim \mathrm{Normal}\left( 0.0, lag_{\omega} \right) \\
pd_{\eta} &\sim \mathrm{MvNormal}\left( pd_{\Omega} \right)
\end{align}

And the same can be done for :pre, :init, :vars, :dynamics, and :derived. The name of each block corresponds exactly to the model macro blocks.

3.3 Rendering and Viewing the Equations

3.3.0.1 In a Julia REPL or VSCode

  • In the REPL, calling latexify(...) directly will show the LaTeX string.
  • In VSCode, the built-in plot pane or an extension can be used to render LaTeX. Alternatively, println(latexify(...)) can be called to see the code, then copy-paste into a LaTeX document or a Markdown cell with MathJax.

Example in a notebook cell:

latexify(warfarin_pkpd_model, :dynamics)

dDepot(t)dt=KaDepot(t)dCentral(t)dt=CLCentral(t)Vc+KaDepot(t)dTurnover(t)dt=koutTurnover(t)+rin(1+emaxCentral(t)Vc(c50+Central(t)Vc))

This will produce a nicely rendered system of ODEs corresponding to the @dynamics block.

Tip

The rendered output in the VSCode plot or in the rendered notebook will be either a CHTML image or an SVG image, depending on the MathJax preferences set by the user. By right-clicking on the image, various options are available to copy the underlying LaTeX code in different formats. The MathML format, once copied, can be directly pasted into a Microsoft Word document, where it will render as an Equation Editor object.

If the output is rendered as an SVG image, it can be saved similarly to a plot by using the “Julia: save plot” icon located at the top of the plot pane. This feature is particularly useful for creating presentations in PowerPoint, as it allows for easy integration of high-quality images.

In cases where the “Julia: save plot” icon does not function, it indicates that the output is being rendered in CHTML format. To change this setting, right-click on the rendered image, navigate to Math Settings, then Math Renderer, and select SVG. This setting is persistent, so it does not need to be adjusted frequently.

3.4 Customizing the Output

Latexify provides additional keyword arguments to tweak how the equations look.

An example customization is shown below:

latexify(warfarin_pkpd_model, :pre; index = :bracket, snakecase = true)

CL=FSZCLpop_CLepk_η[1]Vc=FSZVpop_Vepk_η[2]tabs=pop_tabsepk_η[3]Ka=log(2)tabse0=pop_e0epd_η[1]emax=pop_emaxepd_η[2]c50=pop_c50epd_η[3]tover=pop_toverepd_η[4]kout=log(2)toverrin=e0kout

3.5 Stitching Multiple Blocks Together

A common desire is to produce a single LaTeX document that shows all the relevant equations. Because latexify(model, blockname) returns a LaTeXString, they can be easily combined. For example:

# Collect the blocks needed:
blocks = [:param, :random, :pre, :dynamics, :derived]

# Map each block to its latexify result, then join
model_representation = Latexify.LaTeXString(
    join(
        [
            latexify(warfarin_pkpd_model, block, index = :bracket, snakecase = true) for
            block in blocks
        ],
        "\n\n", # separate by new lines
    ),
);

This will produce a single LaTeX string that can be copied and pasted into a LaTeX document or a Markdown cell with MathJax.

\begin{align}
pop_{CL} &\in RealDomain\left( init = 0.134; lower = 0.0 \right) \\
pop_{V} &\in RealDomain\left( init = 8.11; lower = 0.0 \right) \\
pop_{tabs} &\in RealDomain\left( init = 0.523; lower = 0.0 \right) \\
pop_{lag} &\in RealDomain\left( init = 0.1; lower = 0.0 \right) \\
pop_{e0} &\in RealDomain\left( init = 100.0; lower = 0.0 \right) \\
pop_{emax} &\in RealDomain\left( ; init = -1.0 \right) \\
pop_{c50} &\in RealDomain\left( init = 1.0; lower = 0.0 \right) \\
pop_{tover} &\in RealDomain\left( init = 14.0; lower = 0.0 \right) \\
pk_{\Omega} &\in PDiagDomain\left( \left[
\begin{array}{c}
0.01 \\
0.01 \\
0.01 \\
\end{array}
\right] \right) \\
lag_{\omega} &\in RealDomain\left( init = 0.1; lower = 0.0 \right) \\
pd_{\Omega} &\in PDiagDomain\left( \left[
\begin{array}{c}
0.01 \\
0.01 \\
0.01 \\
0.01 \\
\end{array}
\right] \right) \\
\sigma_{prop} &\in RealDomain\left( init = 0.00752; lower = 0.0 \right) \\
\sigma_{add} &\in RealDomain\left( init = 0.0661; lower = 0.0 \right) \\
\sigma_{fx} &\in RealDomain\left( init = 0.01; lower = 0.0 \right)
\end{align}


\begin{align}
pk_{\eta} &\sim MvNormal\left( pk\_\Omega \right) \\
lag_{\eta} &\sim Normal\left( 0.0, lag\_\omega \right) \\
pd_{\eta} &\sim MvNormal\left( pd\_\Omega \right)
\end{align}


\begin{align}
CL &= FSZCL \cdot pop\_CL \cdot e^{pk\_\eta\left[1\right]} \\
Vc &= FSZV \cdot pop\_V \cdot e^{pk\_\eta\left[2\right]} \\
tabs &= pop\_tabs \cdot e^{pk\_\eta\left[3\right]} \\
Ka &= \frac{\log\left( 2 \right)}{tabs} \\
e0 &= pop\_e0 \cdot e^{pd\_\eta\left[1\right]} \\
emax &= pop\_emax \cdot e^{pd\_\eta\left[2\right]} \\
c50 &= pop\_c50 \cdot e^{pd\_\eta\left[3\right]} \\
tover &= pop\_tover \cdot e^{pd\_\eta\left[4\right]} \\
kout &= \frac{\log\left( 2 \right)}{tover} \\
rin &= e0 \cdot kout
\end{align}


\begin{align}
\frac{\mathrm{d} \cdot Depot(t)}{\mathrm{d}t} &=  - Ka \cdot Depot(t) \\
\frac{\mathrm{d} \cdot Central(t)}{\mathrm{d}t} &= \frac{ - CL \cdot Central(t)}{Vc} + Ka \cdot Depot(t) \\
\frac{\mathrm{d} \cdot Turnover(t)}{\mathrm{d}t} &=  - kout \cdot Turnover(t) + rin \cdot \left( 1 + \frac{emax \cdot Central(t)}{Vc \cdot \left( c50 + \frac{Central(t)}{Vc} \right)} \right)
\end{align}


\begin{align}
conc &\sim Normal\left( cp, \sqrt{\left( \sigma\_prop \cdot cp \right)^{2} + \sigma\_add^{2}} \right) \\
pca &\sim Normal\left( Turnover, \sigma\_fx \right)
\end{align}

To render it, use:

model_representation

popCLRealDomain(init=0.134;lower=0.0)popVRealDomain(init=8.11;lower=0.0)poptabsRealDomain(init=0.523;lower=0.0)poplagRealDomain(init=0.1;lower=0.0)pope0RealDomain(init=100.0;lower=0.0)popemaxRealDomain(;init=1.0)popc50RealDomain(init=1.0;lower=0.0)poptoverRealDomain(init=14.0;lower=0.0)pkΩPDiagDomain([0.010.010.01])lagωRealDomain(init=0.1;lower=0.0)pdΩPDiagDomain([0.010.010.010.01])σpropRealDomain(init=0.00752;lower=0.0)σaddRealDomain(init=0.0661;lower=0.0)σfxRealDomain(init=0.01;lower=0.0)

pkηMvNormal(pk_Ω)lagηNormal(0.0,lag_ω)pdηMvNormal(pd_Ω)

CL=FSZCLpop_CLepk_η[1]Vc=FSZVpop_Vepk_η[2]tabs=pop_tabsepk_η[3]Ka=log(2)tabse0=pop_e0epd_η[1]emax=pop_emaxepd_η[2]c50=pop_c50epd_η[3]tover=pop_toverepd_η[4]kout=log(2)toverrin=e0kout

dDepot(t)dt=KaDepot(t)dCentral(t)dt=CLCentral(t)Vc+KaDepot(t)dTurnover(t)dt=koutTurnover(t)+rin(1+emaxCentral(t)Vc(c50+Central(t)Vc))

concNormal(cp,(σ_propcp)2+σ_add2)pcaNormal(Turnover,σ_fx)

The model_representation can be wrapped in sections or environment tags:

tex_param = latexify(warfarin_pkpd_model, :param; index = :bracket, snakecase = true)
tex_random = latexify(warfarin_pkpd_model, :random; index = :bracket, snakecase = true)
tex_pre = latexify(warfarin_pkpd_model, :pre; index = :bracket, snakecase = true)
tex_dynamics = latexify(warfarin_pkpd_model, :dynamics; index = :bracket, snakecase = true)
tex_derived = latexify(warfarin_pkpd_model, :derived; index = :bracket, snakecase = true)

println("""
        ### Warfarin Model: Parameter Block

        $tex_param

        ### Random Effects

        $tex_random

        ### Pre Block

        $tex_pre

        ### Dynamics Block

        $tex_dynamics

        ### Derived Block

        $tex_derived
        """)

3.5.1 Warfarin Model: Parameter Block

popCLRealDomain(init=0.134;lower=0.0)popVRealDomain(init=8.11;lower=0.0)poptabsRealDomain(init=0.523;lower=0.0)poplagRealDomain(init=0.1;lower=0.0)pope0RealDomain(init=100.0;lower=0.0)popemaxRealDomain(;init=1.0)popc50RealDomain(init=1.0;lower=0.0)poptoverRealDomain(init=14.0;lower=0.0)pkΩPDiagDomain([0.010.010.01])lagωRealDomain(init=0.1;lower=0.0)pdΩPDiagDomain([0.010.010.010.01])σpropRealDomain(init=0.00752;lower=0.0)σaddRealDomain(init=0.0661;lower=0.0)σfxRealDomain(init=0.01;lower=0.0)

3.5.2 Random Effects

pkηMvNormal(pk_Ω)lagηNormal(0.0,lag_ω)pdηMvNormal(pd_Ω)

3.5.3 Pre Block

CL=FSZCLpop_CLepk_η[1]Vc=FSZVpop_Vepk_η[2]tabs=pop_tabsepk_η[3]Ka=log(2)tabse0=pop_e0epd_η[1]emax=pop_emaxepd_η[2]c50=pop_c50epd_η[3]tover=pop_toverepd_η[4]kout=log(2)toverrin=e0kout

3.5.4 Dynamics Block

dDepot(t)dt=KaDepot(t)dCentral(t)dt=CLCentral(t)Vc+KaDepot(t)dTurnover(t)dt=koutTurnover(t)+rin(1+emaxCentral(t)Vc(c50+Central(t)Vc))

3.5.5 Derived Block

concNormal(cp,(σ_propcp)2+σ_add2)pcaNormal(Turnover,σ_fx)

When generating a true .tex file, it may be necessary to escape some backslashes appropriately or place them in a valid LaTeX environment (e.g., \begin{document} ... \end{document}).

3.6 A Example in Action

Here’s a compact snippet that:

  1. Defines the Warfarin model (using the code above).
  2. latexifys the :dynamics block.
  3. Renders it in a Quarto notebook:
dynamics_equations = latexify(warfarin_pkpd_model, :dynamics)

# In a notebook, just do:
dynamics_equations

dDepot(t)dt=KaDepot(t)dCentral(t)dt=CLCentral(t)Vc+KaDepot(t)dTurnover(t)dt=koutTurnover(t)+rin(1+emaxCentral(t)Vc(c50+Central(t)Vc))

This will display nicely typeset ODEs for Depot, Central, and Turnover, including the definitions of the rates.

3.7 Summary and Next Steps

  • Latexify offers a powerful way to automatically generate mathematical expressions from Pumas models.
  • Focus can remain on model building while letting Latexify handle the tedious transcription.
  • Different blocks of the model (parameters, random effects, pre, init, dynamics, derived, etc.) can be latexified individually or combined into a single string.
  • Customization options (keyword arguments) allow control over styles.

For more advanced usage, refer to the Latexify Documentation for formatting and environment options.