NIBR data analysis

vcat(pk_log_plots, avg_pk_plots) would be more efficient usually.

1 Like

@Michael I tried your suggestion. see below:

all_plots = vcat(pk_log_plots, avg_pk_plots) #also tried embedding this command in report().
report(all_plots, 
      title = "Individual Semi-Log Exploratory Plots", 
      output = "./practice/nibr/data/plots",
      header = "Single Ascending Dose",
      footer = "Confidential")
7-element Vector{Any}:
 Figure()
 Figure()
 Figure()
 Figure()
 Figure()
 Figure()
 FigureAxisPlot()

[ Info: Calculating required `inspect` results.
ERROR: LoadError: Unknown argument type `Figure`.
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:33
  [2] (::PumasReports.var"#categorise#4"{Vector{Any}, Vector{Any}, Vector{Any}, Vector{Any}})(f::Figure)
    @ PumasReports /builds/PumasAI/PumasSystemImages-jl/.julia/packages/PumasReports/1jlEU/src/PumasReports.jl:293
  [3] foreach(f::PumasReports.var"#categorise#4"{Vector{Any}, Vector{Any}, Vector{Any}, Vector{Any}}, itr::Tuple{Figure})
    @ Base ./abstractarray.jl:2141
  [4] PumasReports.ModelBundle(args::Tuple{Figure}; inspect::Bool, infer::Bool)
    @ PumasReports /builds/PumasAI/PumasSystemImages-jl/.julia/packages/PumasReports/1jlEU/src/PumasReports.jl:294
  [5] PumasReports.ModelBundle(args::Figure; inspect::Bool, infer::Bool)
    @ PumasReports /builds/PumasAI/PumasSystemImages-jl/.julia/packages/PumasReports/1jlEU/src/PumasReports.jl:333
  [6] (::PumasReports.var"#14#15"{Base.Iterators.Pairs{Symbol, Bool, Tuple{Symbol, Symbol}, NamedTuple{(:inspect, :infer), Tuple{Bool, Bool}}}})(::Tuple{Int64, Figure})
    @ PumasReports ./none:0
  [7] iterate
    @ ./generator.jl:47 [inlined]
  [8] collect(itr::Base.Generator{Base.Iterators.Enumerate{Vector{Any}}, PumasReports.var"#14#15"{Base.Iterators.Pairs{Symbol, Bool, Tuple{Symbol, Symbol}, NamedTuple{(:inspect, :infer), Tuple{Bool, Bool}}}}})
    @ Base ./array.jl:678
  [9] #_named_objects#13
    @ /builds/PumasAI/PumasSystemImages-jl/.julia/packages/PumasReports/1jlEU/src/PumasReports.jl:343 [inlined]
 [10] report(fitted_models::Vector{Any}; output::String, title::String, version::VersionNumber, date::Dates.DateTime, author::String, categorical::Vector{Any}, force::Bool, inspect::Bool, infer::Bool, clean::Bool, header::String, footer::String, plot_fontsize::Int64, plot_resolution::Tuple{Int64, Int64})
    @ PumasReports /builds/PumasAI/PumasSystemImages-jl/.julia/packages/PumasReports/1jlEU/src/PumasReports.jl:177
 [11] top-level scope
    @ ~/data/code/practice/nibr/data/dataprep

The last element in that all_plots vector will be what’s confusing report. It’s a FigureAxisPlot object rather than the expected Figure that report accepts. Was your avg_pk_plots created with paginate = true? You can either set paginate = true or “split apart” the output of summary_observations_vs_time like so:

avg_pk_plots, _, _ = summary_observations_vs_time(...)

which will just keep the Figure part of FigureAxisPlot and drop the Axis and Plot.

This is something that could probably be automated by report internally. I’ll have a look at adding support for that in a future release of Pumas.

@Michael paginate=false for the avg plots. See below.

#
pk_log_plots = observations_vs_time(pop_nca, 
                                axis = (xlabel = "Time (hours)", 
                                        ylabel = "Concentration (μg/mL)",
                                        yscale=log10, ytickformat=x -> string.(round.(x; digits=1)), 
                                        yticks = [0.1, 1, 10],
                                        ygridwidth = 3, 
                                        yminorticksvisible = true,
                                        yminorgridvisible = true,
                                        yminorticks = IntervalsBetween(10),
                                        xminorticksvisible = true,
                                        xminorgridvisible = true,
                                        xminorticks = IntervalsBetween(5),
                                        limits = (nothing, nothing, nothing, 20),
                                        spinewidth = 2),
                                        columns = 4, rows = 3,
                                        paginate = true,
                                        markersize = 8,
                                        markercolor = :blue,
                                        linewidth = 2,
                                        color = :green,
                                        facet = ( combinelabels = true,))
#
pk_log_plots[1]
pk_log_plots[6]

avg_pk_plots = summary_observations_vs_time(pop_nca, 
                                axis = (xlabel    = "Time, hr", 
                                        ylabel    = "Concentration, ug/L",
                                        xlim      = (0,72),
                                        xticks    = 0:12:72),
                                markersize = 1,
                                paginate  = false,
                                columns   = 3, rows = 2,
                                facet     = (combinelabels = true,))
#
all_plots = vcat(pk_log_plots, avg_pk_plots)
report(all_plots, 
      title = "Individual Semi-Log Exploratory Plots", 
      output = "./practice/nibr/data/plots",
      header = "Single Ascending Dose",
      footer = "Confidential")

Setting paginate = true should be what’s needed for that summary_observations_vs_time plot to work correctly in report.

@Michael That worked! Thank you. Is there a way to make the size of the plots constant. You saw my code above.
See the last plot.

@julius I tried your trick, In general it seems to work. Changing the range of the yaxis for semi-log plot seems to be a challenge. Tried few approaches to make the ylims between 0.01, 100; yticks [0.1,1,10,100]. The code runs without any error; but plot is unaffected. The top part of the plot is chopped off (see below).

pk_log_plots = observations_vs_time(pop_nca, 
                                axis = (xlabel = "Time (hours)", 
                                        ylabel = "Concentration (μg/mL)",
                                        yscale=log10, 
                                        ytickformat=x -> string.(round.(x; digits=1)), 
                                        xlims = (0,72),
                                        xticks = (0:12:72),
                                        ylims = (0.01,100),
                                        yticks = [0.1, 1, 10, 100],
                                        #yticks = (0.01:10:100),
                                        #yticks = LogTicks(IntegerTicks()),
                                        ygridwidth = 3, 
                                        yminorticksvisible = true,
                                        yminorgridvisible = true,
                                        yminorticks = IntervalsBetween(10),
                                        xminorticksvisible = true,
                                        xminorgridvisible = true,
                                        xminorticks = IntervalsBetween(5),
                                        limits = (nothing, nothing, nothing, 20),
                                        spinewidth = 2),
                                        columns = 2, rows = 2,
                                        paginate = true,
                                        markersize = 8,
                                        markercolor = :blue,
                                        linewidth = 2,
                                        color = :green,
                                        facet = ( combinelabels = true,))
#

as in my previous post, use limits = (nothing, nothing, nothing , 100) instead, where each argument is the xlow, xhigh, ylow, yhigh

Worked like a charm. Thank you.

1 Like

@vijay @Michael any advice on the squished graph on the last page?

it has to do with the combination of paginate, limit and facet. Can you try running this code that I posted earlier

using Pumas
using PumasUtilities
using Dates
using CairoMakie
using AlgebraOfGraphics
using DataFramesMeta
df = CSV.read("Single_Ascending_Dose_Dataset2.csv", DataFrame, missingstrings = ["NA"])
df = @rsubset df :TIME >= 0
@rtransform! df route = "ev"
##  map NCA data from CSV
ncapop = read_nca(df,
                id  = :ID, time = :NOMTIME, amt = :AMT, route = :route,
                observations = :LIDV, 
                group = [:DOSE])

#plot means               
avg_pk_plots = summary_observations_vs_time(ncapop, 
                                  color = :black, linewidth = 2, whiskerwidth = 8,
                                  paginate = true,
                                  limit = 1,
                                  axis = (xlabel = "Time (hours)", 
                                            ylabel = "Concentration (μg/mL)",
                                            yscale=log10, ytickformat=x -> string.(round.(x; digits=1)), 
                                            yticks = [0.1, 1, 10],
                                            ygridwidth = 3, 
                                            yminorgridcolor = :darkgrey,
                                            yminorticksvisible = true,
                                            yminorgridvisible = true,
                                            yminorticks = IntervalsBetween(10),
                                            xminorticksvisible = true,
                                            xminorgridvisible = true,
                                            xminorticks = IntervalsBetween(5),
                                            limits = (nothing, nothing, nothing, 30),
                                            spinewidth = 2),
                                    #facet = ( combinelabels = true,),
                                    figure = (resolution = (800,800),
                                                        fontsize = 22),)
                                  
pk_log_plots = observations_vs_time(ncapop, paginate = true, color = :black,
                                axis = (xlabel = "Time (hours)", 
                                        ylabel = "Concentration (μg/mL)",
                                        yscale=log10, ytickformat=x -> string.(round.(x; digits=1)), 
                                        yticks = [0.1, 1, 10],
                                        ygridwidth = 3, 
                                        yminorticksvisible = true,
                                        yminorgridvisible = true,
                                        yminorticks = IntervalsBetween(10),
                                        xminorticksvisible = true,
                                        xminorgridvisible = true,
                                        xminorticks = IntervalsBetween(5),
                                        limits = (nothing, nothing, nothing, 20),
                                        spinewidth = 2),
                                        columns = 4, rows = 4,
                                        facet = ( combinelabels = true,))

report([pk_log_plots..., avg_pk_plots...],
      title = "Summary and Individual Semi-Log Exploratory Plots", 
      output = "./practice/nibr/data/plots",
      header = "Single Ascending Dose",
      footer = "Confidential")

@vijay I ran your code. The last page graphs being different persists. My code also rendered similar discrepancy.

For the avg plots, you opted for one plot per page. That works sure. But I was trying 2x2 format to be concise.
Bob

@vijay I tried to consolidate your suggested data formatting code, like so. It worked, it took me a good 30min to figure out that > requires a dot and the other commands do not. I do see other commands such as ‘transform’ also are to be performed per row, but they do not require a dot. Anyway it works for me. except for the graph on the last page.

df = @chain df begin
        @subset(:TIME .>= 0.0)
        @transform(:ROUTE = "ev")
        @transform(:AMT_UG = :AMT * 1000)
        @select(:ID,:TIME,:NOMTIME,:AMT_UG,:CONC = :LIDV,:CMT,:EVID,:WEIGHTB,:SEX,:DOSE,:ROUTE)
end

@vijay the above code for avg PK profiles produces a plot for each dose. I am used to facet_wrap() type command. I did not know which argument to use to switch-off a facet plot. That is, to have the avg plots across all doses on one plot.
Bob

@vijay @storopoli There appears to be one extra minortick while plotting logscale. Is there a way to fix this?

it would be good to see this in a visual. There are many plots in this thread. Can you please markup and show where exactly?