Model Representation

Once your model is defined and your analyses run, you might want to show it off. Transcribing a computational model to equations can be a tedious affair. Furthermore, it is easy to get something wrong. We, therefore, sought to automate that process.

Models created by the Pumas model macros contain all the information we need to automatically generate LaTeX equations from them using Latexify.

Let us first define a model so that we can demonstrate this:

using Pumas

model = @model begin
    @param begin
        tvcl ∈ RealDomain(; lower = 0)
        tvvc ∈ RealDomain(; lower = 0)
        tvq ∈ RealDomain(; lower = 0)
        tvvp ∈ RealDomain(; lower = 0)
        Ω_pk ∈ PDiagDomain(4)
        σ_prop_pk ∈ RealDomain(; lower = 0)
        # PD parameters
        tvturn ∈ RealDomain(; lower = 0)
        tvebase ∈ RealDomain(; lower = 0)
        tvec_50 ∈ RealDomain(; lower = 0)
        Ω_pd ∈ PDiagDomain(1)
        σ_add_pd ∈ RealDomain(; lower = 0)
    end

    @random begin
        ηpk ~ MvNormal(Ω_pk)
        ηpd ~ MvNormal(Ω_pd)
    end

    @pre begin
        CL = tvcl * exp(ηpk[1])
        Vc = tvvc * exp(ηpk[2])
        Q = tvq * exp(ηpk[3])
        Vp = tvvp * exp(ηpk[4])

        e_base = tvebase * exp(ηpd[1])
        ec_50 = tvec_50
        e_max = 1
        turn = tvturn
        k_out = 1 / turn
        k_in0 = e_base * kout
    end

    @init begin
        Resp = e_base
    end

    @vars begin
        conc := Central / Vc
        e_drug := e_max * conc / (ec_50 + conc)
        k_in := k_in0 * (1 - e_drug)
    end

    @dynamics begin
        Central' = -(CL / Vc) * Central + (Q / Vp) * Peripheral - (Q / Vc) * Central
        Peripheral' = (Q / Vc) * Central - (Q / Vp) * Peripheral
        Resp' = k_in - k_out * Resp
    end

    @derived begin
        dv ~ @. Normal(conc, sqrt(conc^2 * σ_prop_pk))
        resp ~ @. Normal(Resp, sqrt(σ_add_pd))
    end
end

PumasModel
  Parameters: tvcl, tvvc, tvq, tvvp, Ω_pk, σ_prop_pk, tvturn, tvebase, tvec_50, Ω_pd, σ_add_pd
  Random effects: ηpk, ηpd
  Covariates: 
  Dynamical system variables: Central, Peripheral, Resp
  Dynamical system type: Nonlinear ODE
  Derived: dv, resp
  Observed: dv, resp

If we want to extract LaTeX-formatted equations from the model then we first need to load Latexify:

using Latexify

We can now latexify the different model blocks using calls like latexify(model, blockname) where blockname is a Symbol that indicates which model block to latexify. For example:

latexify(model, :random)

\begin{align*} {\eta}pk &\sim \mathrm{MvNormal}\left( \Omega_{pk} \right) \\ {\eta}pd &\sim \mathrm{MvNormal}\left( \Omega_{pd} \right) \end{align*}

latexify(model, :pre)

\begin{align*} CL &= tvcl \cdot e^{{\eta}pk_{1}} \\ Vc &= tvvc \cdot e^{{\eta}pk_{2}} \\ Q &= tvq \cdot e^{{\eta}pk_{3}} \\ Vp &= tvvp \cdot e^{{\eta}pk_{4}} \\ e_{base} &= tvebase \cdot e^{{\eta}pd_{1}} \\ ec_{50} &= tvec_{50} \\ e_{max} &= 1 \\ turn &= tvturn \\ k_{out} &= \frac{1}{turn} \\ k_{in0} &= e_{base} \cdot kout \end{align*}

latexify(model, :dynamics)

\begin{align*} \frac{dCentral(t)}{dt} =& \frac{Q \cdot Peripheral(t)}{Vp} - \frac{CL \cdot Central(t)}{Vc} - \frac{Q \cdot Central(t)}{Vc} \\ \frac{dPeripheral(t)}{dt} =& \frac{Q \cdot Central(t)}{Vc} - \frac{Q \cdot Peripheral(t)}{Vp} \\ \frac{dResp(t)}{dt} =& k_{in0} \cdot \left( 1 - \frac{e_{max} \cdot Central(t)}{\left( ec_{50} + \frac{Central(t)}{Vc} \right) \cdot Vc} \right) - k_{out} \cdot Resp(t) \end{align*}

latexify uses the information you input via the model macro, so some models contain information that others do not. Which blocknames are valid thus differ depending on the model that you use. In this example, the blocknames :param, :random, :pre, :init, :vars, :dynamics and :derived would be valid since these are the model blocks that we explicitly defined.

Using the latexify output

The latexify function returns a LaTeXString, which is pretty much a standard String with some overloads for making it display nicely. You can print this string to see the generated LaTeX code:

println(latexify(model, :dynamics))

\begin{align}
\frac{\mathrm{d} \cdot Central(t)}{\mathrm{d}t} &= \frac{Q \cdot Peripheral(t)}{Vp} + \frac{ - CL \cdot Central(t)}{Vc} + \frac{ - Q \cdot Central(t)}{Vc} \\
\frac{\mathrm{d} \cdot Peripheral(t)}{\mathrm{d}t} &= \frac{ - Q \cdot Peripheral(t)}{Vp} + \frac{Q \cdot Central(t)}{Vc} \\
\frac{\mathrm{d} \cdot Resp(t)}{\mathrm{d}t} &= k_{in0} \cdot \left( 1 + \frac{ - e_{max} \cdot Central(t)}{Vc \cdot \left( ec_{50} + \frac{Central(t)}{Vc} \right)} \right) - k_{out} \cdot Resp(t)
\end{align}

Different editors/environments have different support for more complex rendering. latexify automatically displays as nicely rendered equations in some environments, notably including notebooks. In other environments, like VSCode, you can call render to render the equations:

render(latexify(model, :dynamics))

When we created Latexify, we first thought that rendering would be somewhat frivolous, but quick access to readable equations has proven helpful during model development. Try for yourself.

Configurations

You can tune the latexify output to different keyword arguments. Many are provided directly from Latexify package, and some are specific to Pumas models.

For example:

latexify(model, :dynamics; show_t = false, italicize = false, index = :bracket)

\begin{align*} \frac{\mathrm{dCentral}}{dt} =& \frac{Q \cdot \mathrm{Peripheral}}{Vp} - \frac{CL \cdot \mathrm{Central}}{Vc} - \frac{Q \cdot \mathrm{Central}}{Vc} \\ \frac{\mathrm{dPeripheral}}{dt} =& \frac{Q \cdot \mathrm{Central}}{Vc} - \frac{Q \cdot \mathrm{Peripheral}}{Vp} \\ \frac{\mathrm{dResp}}{dt} =& k_{in0} \cdot \left( 1 - \frac{e_{max} \cdot \mathrm{Central}}{\left( ec_{50} + \frac{\mathrm{Central}}{Vc} \right) \cdot Vc} \right) - k_{out} \cdot \mathrm{Resp} \end{align*}

For more options and features, have a look at Latexify's documentation.

Stitching things together

Since the latexify output works well with other strings, you can quickly stitch different model parts together

model_representation = """
  \\section{Model dynamics}
  Here we might want to describe the model dynamics, as given by
  $(latexify(model, :dynamics)).

  \\section{Random effects}
  The random effects are given by
  $(latexify(model, :random)).

  \\section{Conclusion}
  You have a lot of power to automate your model representation.
  Don't forget to escape your backslashes.
  """

println(model_representation)

\section{Model dynamics}
Here we might want to describe the model dynamics, as given by
\begin{align}
\frac{\mathrm{d} \cdot Central(t)}{\mathrm{d}t} &= \frac{Q \cdot Peripheral(t)}{Vp} + \frac{ - CL \cdot Central(t)}{Vc} + \frac{ - Q \cdot Central(t)}{Vc} \\
\frac{\mathrm{d} \cdot Peripheral(t)}{\mathrm{d}t} &= \frac{ - Q \cdot Peripheral(t)}{Vp} + \frac{Q \cdot Central(t)}{Vc} \\
\frac{\mathrm{d} \cdot Resp(t)}{\mathrm{d}t} &= k_{in0} \cdot \left( 1 + \frac{ - e_{max} \cdot Central(t)}{Vc \cdot \left( ec_{50} + \frac{Central(t)}{Vc} \right)} \right) - k_{out} \cdot Resp(t)
\end{align}
.

\section{Random effects}
The random effects are given by
\begin{align}
{\eta}pk &\sim \mathrm{MvNormal}\left( \Omega_{pk} \right) \\
{\eta}pd &\sim \mathrm{MvNormal}\left( \Omega_{pd} \right)
\end{align}
.

\section{Conclusion}
You have a lot of power to automate your model representation.
Don't forget to escape your backslashes.

A note about future releases

The overall functionality described here is considered part of our public API and should thus not break during anything but a major release (Pumas follow semantic versioning). However, only the intention behind a latexify command should be considered stable across non-breaking releases, not the exact formatting of the output string.