Initial feedback

Overall good first experience, basic functionality works, looks very promising. Below some initial feedback / questions.

Fitting

  • Overall performance is fast, fits feels broadly similar in speed to NONMEM (tested only on very simple problems so far). Will do more benchmarking later.
  • Loading the Pumas package does take some time, about 35secs on a 2017 i7 MBP, but not a problem in most use-cases since you only have to recompile the package once in each session.
  • Is there any way I can get information on the iterations / warnings / etc that is going on when doing a fit? In the JuliaPro interface I’m only seeing the “cogs” icon when something is running.
  • Although time-varying covariates seem to be parsed fine using read-pumas(), when using this parse dataset in the fit() it complains with MethodError: no method matching ^(::Array{Float64,1}, ::Float64) which seems to indicate vectors of covariates are not allowed yet? Am I doing something wrong, or is this a feature that is coming later?

VPC

I realize this functionality is alpha, but here are some comments anyhow:

  • why are there two light-red areas in the VPC plot? I guess one is for the prediction interval for the median (of the simulated data) and the other the uncertainty of the observed data? If so, I wouldn’t do that, at least not as default, as it makes the plot too cluttered and hard to interpret.
  • I suppose the red line is the median of the observed, but what are the pink dots then? Most pharmacometricians will intuitively think these are the actual observations, but that doesn’t seem to be the case?
  • binning is an essential feature of the VPC, I don’t see that as an option yet? Should be easy to add at least simple k-means clustering and custom binning options.
  • personal preference: i would keep the colors of the prediction intervals (e.g. 5th and 95 percentile) the same (now blue and green), as they provide similar types of information.

Misc suggestions:

  • Would you be willing to share the development roadmap, or at least a broad version of it? This would help us decide when / whether this is likely to become a possible alternative in production use for us.
  • in general, as someone who is not too very familiar with Julia yet, it requires quite a bit of tinkering to get started beyond the few provided examples (e.g. using NONMEM datasets). In the R world, it is common that R package developers guide the e.g. by adding useful error / warning message that catch common mistakes when calling functions. Not sure if this is common in the Julia world, but would make transition easier for novices.
  • what exactly is the difference is between read_pumas() and process_nmtran()? After some tinkering I managed to read a NONMEM file with read_pumas, it seems the latter is not available yet? I think it would be useful for Julia novices to show with a simple example how to get a NONMEM dataset (csv) into Julia and fit it.
  • BayesMCMC: can you add some more detailed info in the manual on how to run using BayesMCMC, just swapping Pumas.FOCEI for BayesMCMC doesn’t work. Questions e.g. :
    • how do I specify the priors?
    • what MCMC samplers are available, which one is used by default?
    • does this use a separate library in the backend (MCMC?), does it need to be installed?

Dear Ron,

Thank you for the detailed feedback. While I will take a shot at responding to most of the things, I am sure our team will give more details.

Fitting

  • Overall performance is fast, fits feels broadly similar in speed to NONMEM (tested only on very simple problems so far). Will do more benchmarking later.

We have an entire suite of benchmarks that we will be putting out soon. We are tinkering on some internals and once that is done, we will be able run the comparisons easily.

  • Loading the Pumas package does take some time, about 35secs on a 2017 i7 MBP, but not a problem in most use-cases since you only have to recompile the package once in each session.

Yes, we realize this compile-time issue. @pkofod and @ChrisRackauckas are currently working on cutting them down. It is WIP

  • Is there any way I can get information on the iterations / warnings / etc that is going on when doing a fit? In the JuliaPro interface I’m only seeing the “cogs” icon when something is running.

In our initial releases we turned off the output printing by default. @pkofod is working on getting that back after our discussions at JuliaCon. But, if you would like to see the traces after fitting, you can access the fields of the fit object. E.g. if you have myfit = fit(..), then myfit.model, myfit.optim.trace etc… We have documentation coming out on that soon.

  • Although time-varying covariates seem to be parsed fine using read-pumas() , when using this parse dataset in the fit() it complains with MethodError: no method matching ^(::Array{Float64,1}, ::Float64) which seems to indicate vectors of covariates are not allowed yet? Am I doing something wrong, or is this a feature that is coming later?

Time-varying covariates work with ODE’s as of now. There is a special syntax in the @model macro to define time-varying covaraites. Again, a tutorial is ready and will put this out soon. Time-varying covariates on analytical solutions is WIP and we anticipate that to be ready by ACoP

VPC

I realize this functionality is alpha, but here are some comments anyhow:

There is much to be done here, and we just got started. Your experience and feedback would be really helpful here. We will answer these questions separately.

  • why are there two light-red areas in the VPC plot? I guess one is for the prediction interval for the median (of the simulated data) and the other the uncertainty of the observed data? If so, I wouldn’t do that, at least not as a default, as it makes the plot too cluttered and hard to interpret.
  • I suppose the red line is the median of the observed, but what are the pink dots then? Most pharmacometricians will intuitively think these are the actual observations, but that doesn’t seem to be the case?
  • binning is an essential feature of the VPC, I don’t see that as an option yet? Should be easy to add at least simple k-means clustering and custom binning options.
  • personal preference: i would keep the colors of the prediction intervals (e.g. 5th and 95 percentile) the same (now blue and green), as they provide similar types of information.

Misc suggestions:

  • Would you be willing to share the development roadmap, or at least a broad version of it? This would help us decide when / whether this is likely to become a possible alternative in production use for us.

Yes, we plan to do that before ACoP.

  • in general, as someone who is not too very familiar with Julia yet, it requires quite a bit of tinkering to get started beyond the few provided examples (e.g. using NONMEM datasets). In the R world, it is common that R package developers guide the e.g. by adding useful error / warning message that catch common mistakes when calling functions. Not sure if this is common in the Julia world, but would make transition easier for novices.

This has been on our to-do list. Providing appropriate error messages is hard and comes with user feedback and testing. We are working on getting more and more simple errors out, but this may take a little bit of time. Maybe @ChrisRackauckas or @andreasnoack can address this question better.

  • what exactly is the difference is between read_pumas() and process_nmtran() ? After some tinkering I managed to read a NONMEM file with read_pumas , it seems the latter is not available yet? I think it would be useful for Julia novices to show with a simple example how to get a NONMEM dataset (csv) into Julia and fit it.

While we started with process_nmtran which does everything that nmtran datasets require, we realized that a lot of the requirements can be overriden by easier Julia based features. read_pumas will be compatible with nmtran style datasets, but will eventually morph into its own style of data set preparations that will utilize the power of Julia. @ChrisRackauckas can expand more here.

  • BayesMCMC: can you add some more detailed info in the manual on how to run using BayesMCMC, just swapping Pumas.FOCEI for BayesMCMC doesn’t work. Questions e.g. :
    • how do I specify the priors?
    • what MCMC samplers are available, which one is used by default?
    • does this use a separate library in the backend ( MCMC ?), does it need to be installed?

Documentation on this is sparse and there is some work going on currently. To be honest, the only change here is that in the @param block, you direclty specify the distribution of the prior instead of the bounds. And then call the Pumas.BayesMCMC() in the fit function. We are working on post-processing the Bayes fits and other details. We plan to have a ready to use Bayes in month or two.

Please look out for the roadmap for other details.

Indeed, the benchmarks will be automatically runnable and updated to each version. I hope we can get those out by ACoP.

The JuliaPro registry currently has old versions of the DifferentialEquations packages which introduced some compile time regressions due to interactions with Cassette. These were subsequently made optional, dropping the DelayDiffEq (a dependency of Pumas) compile time from 55 seconds to 0.2 on Win32 platforms. Source: https://github.com/JuliaDiffEq/DelayDiffEq.jl/issues/125 . Additionally, some inference fixes dropped the total test time of DelayDiffEq.jl from 44 minutes on Julia 1.0 to 13 minutes on Julia 1.3, with a good portion of that being compile time (since the delay differential equation solvers are the heaviest on compile time), and we will probably drop the direct dependency on them in the next Pumas release as well (and instead utilize abstract type handling).

There is another compile time issue in Pumas which I think I am tracking down to AdvancedHMC.jl, and we’ll be working that out similarly.

I anticipate finishing the PR by the end of next week: https://github.com/PumasAI/Pumas.jl/pull/607 . Of course, it’ll take longer to release.

Not only are we adding more error messages on the Pumas level, but also looking to change the default error messages of Julia itself to be more friendly to new users. I’ve discussed my thoughts on this in the compiler triage calls, and they were mentioned in Jeff Bezanson’s JuliaCon talk, so it seems they are being taken seriously and will likely get into one of the next Julia versions.

process_nmtran isn’t a function that exists in the released versions. The name change occurred since what we parse to isn’t exactly the NMTRAN format, but we do have things in our parser that allow for NMTRAN files to work with our parser. There are many computational advantages that we get by using our more column-based system.

We used DynamicHMC.jl as the backend for the release version, but this was not producing stable estimates and thus we are changing to AdvancedHMC.jl. Not only is this faster, but it requires less user input and is very robust. We haven’t really detailed the Bayes functionality yet since we knew we were going to swap this out.

I will sharing a whole lot on this quite soon :smile: though it’s blocked right now for various reasons. Expect a very fun blog post.

Thanks for the detailed feedback Vijay and Chris, this is very helpful. Will be following upcoming releases and blog posts closely!

We expect at least one release every month till ACoP. Blog posts will go along with the release. Will keep you posted!