How to manipulate diagnostic plots

Diagnostics
1

Hello,

I was wondering if there is a way to manipulate diagnostic plots (for example, observations_vs_ipredictions) in Pumas?

I want to change x and y axis labels, ticks and range, marker shape etc.

2

Hello @pawangpt700

All plots generated using Pumas have extensive support for customization since they use Makie.jl as their backend.

Here’s an example of how you could customize your results from observations_vs_ipredictions. I’ll use the nlme_sample dataset from PharmaDatasets, and the following code to get the dataset, define a one-compartment model, and then fit it and inspect it:

using Pumas, PumasUtilities
using PharmaDatasets # Get dataset
using CairoMakie # Customize plots

# Get dataset
pkdata = dataset("nlme_sample")

# Create popoulation with read_pumas
population = read_pumas(
    pkdata;
    id = :ID,
    time = :TIME,
    amt = :AMT,
    covariates = [:WT, :AGE, :SEX, :CRCL, :GROUP],
    observations = [:DV],
    cmt = :CMT,
    evid = :EVID,
    rate = :RATE,
)

# Create the model
model = @model begin
    @metadata begin
        desc = "base model: 1comp"
        timeu = u"hr"
    end

    @param begin
        "Clearance (L/hr)"
        tvcl ∈ RealDomain(; lower = 0.0001)
        "Volume (L)"
        tvvc ∈ RealDomain(; lower = 0.0001)
        """
        - ΩCL
        - ΩVc
        """
        Ω ∈ PDiagDomain(2) #covariance matrix 
        "Additive RUV"
        σ²_add ∈ RealDomain(; lower = 0.0001) # variance 
        "Proportional RUV"
        σ²_prop ∈ RealDomain(; lower = 0.0001) # variance  
    end

    @random begin
        η ~ MvNormal(Ω)
    end

    @covariates WT AGE SEX CRCL GROUP

    @pre begin
        CL = tvcl * exp(η[1])
        Vc = tvvc * exp(η[2])
    end

    @dynamics Central1

    @derived begin
        CONC := @. Central / Vc
        """
        DrugY Concentration (ng/mL)
        """
        DV ~ @. Normal(CONC, sqrt(CONC^2 * σ²_prop + σ²_add)) # using variance
    end
end

params = (;
    tvvc = 5, 
    tvcl = 0.2, 
    Ω = Diagonal([0.01, 0.01]), 
    σ²_add = 0.01, 
    σ²_prop = 0.01
)

# Fit the model
fitted_model = fit(model, population, params, FOCE())

# Inspect the model
inspection = inspect(fitted_model)

Now we can use inspection to generate our observations_vs_ipredictions with the following line:

observations_vs_ipredictions(inspection)

However, if you want to customize the appearance of the resulting plot, you have to use the supported style keyword arguments for the function. For instance, we can use markercolor and markersize to change the look of the markers:

observations_vs_ipredictions(
    inspection;
    markercolor=:grey,
    markersize=5
)

We can also customize the look of the identity line, OLS fit, and loess fit. For example, we could hide the OLS line, display the identity line as a dashed line and change the color of the loess fit line:

observations_vs_ipredictions(
    inspection;
    markercolor=:grey,
    markersize=5,
    ols=false,
    identity_color=:black,
    identity_linestyle=:dash,
    loess_color=:blue
)

Lastly, a lot of customization can be achieved through the axis and figure keyword arguments. You can pass a named tuple containing any argument supported by Makie.jl’s Axis or Figure objects. As an example, we will change the axis’ labels, title and ticks, and the figure’s resolution and font size:

observations_vs_ipredictions(
    inspection;
    markercolor=:grey,
    markersize=5,
    ols=false,
    identity_color=:black,
    identity_linestyle=:dash,
    loess_color=:blue,
    axis=(;
        title="Diagnostic plot",
        xlabel="Individual predictions (ng/mL)",
        ylabel="Observed Conc. (ng/mL)",
        xticks=50:50:550
    ),
    figure=(;
        resolution=(1200, 800),
        fontsize=25
    ),
)

There is a nice tutorial on Customization of AlgebraOfGraphics.jl Plots where you can learn much more about the axis and figure keyword arguments.

One last thing you might want to do is to get a legend for your plot. You can achieve this by using the figurelegend function:

fig, ax, plt = observations_vs_ipredictions(
    inspection;
    markercolor=:grey,
    markersize=5,
    ols=false,
    identity_color=:black,
    identity_linestyle=:dash,
    loess_color=:blue,
    axis=(;
        title="Diagnostic plot",
        xlabel="Individual predictions (ng/mL)",
        ylabel="Observed Conc. (ng/mL)",
        xticks=50:50:550
    ),
    figure=(;
        resolution=(1200, 800),
        fontsize=25
    ),
)
figurelegend(fig)
fig

Notice that I had to store the Figure object returned by observations_vs_ipredictions (fig variable). You can do this with all Pumas and Makie.jl functions since they return a tuple containing the figure, axis, and plot. This should produce the following plot (using the nlme_sample dataset and the one-compartment model):

plot_24
plot_241200×800 67.5 KB

1 Like