Fix parameters during full Bayesian estimate

1

Dear Pumas team,

I had a question related to fix some parameters in full Bayesian approach for popPK params estimate.
I tried for using: constantcoef = (ωKa = 0.2, ωCL = 0.2, ωVc = 0.2) in theopmodel_bayes_v2 in this tutorial but it had an error:

ERROR: MethodError: no method matching fit(::PumasModel{ParamSet{NamedTuple{(:θ, :ωKa, :ωCL, :ωVc, :σ), Tuple{Constrained{DiagNormal, VectorDomain{Vector{Float64}, Vector{Float64}, Vector{Float64}}}, Gamma{Float64}, Gamma{Float64}, Gamma{Float64}, Gamma{Float64}}}}, var"#5#19", var"#6#20", var"#8#22", var"#10#24", Depots1Central1, var"#11#25", var"#15#29", Nothing}, ::Vector{Subject{NamedTuple{(:DV,), Tuple{Vector{Union{Missing, Float64}}}}, Pumas.ConstantCovar{NamedTuple{(:WT,), Tuple{Float64}}}, Vector{Pumas.Event{Float64, Float64, Float64, Float64, Float64, Float64, Int64}}, Vector{Float64}}}, ::NamedTuple{(:θ, :ωKa, :ωCL, :ωVc, :σ), Tuple{Vector{Float64}, Float64, Float64, Float64, Float64}}, ::Pumas.BayesMCMC; constantcoef=(ωKa = 0.2, ωCL = 0.2, ωVc = 0.2), nsamples=2000, nadapts=1000)
Closest candidates are:
  fit(::PumasModel, ::AbstractVector{T} where T<:Subject, ::NamedTuple, ::Pumas.BayesMCMC; target_accept, nadapts, nsamples, progress, nchains, rng, diffeq_options) at /builds/PumasAI/PumasSystemImages-jl/.julia/packages/Pumas/XHa5R/src/estimation/bayes.jl:109 got unsupported keyword argument "constantcoef"
  fit(::PumasModel, ::AbstractVector{T} where T<:Subject, ::NamedTuple, ::Union{Pumas.LikelihoodApproximation, Pumas.MAP}; optimize_fn, constantcoef, omegas, ensemblealg, checkidentification, diffeq_options) at /builds/PumasAI/PumasSystemImages-jl/.julia/packages/Pumas/XHa5R/src/estimation/likelihoods.jl:2372 got unsupported keyword arguments "nsamples", "nadapts"
  fit(::Pumas.PumasEMModel, ::AbstractVector{var"#s35"} where var"#s35"<:Pumas.AbstractSubject, ::NamedTuple, ::Union{Pumas.LikelihoodApproximation, Pumas.MAP}; optimize_fn, constantcoef, ensemblealg, checkidentification, diffeq_options) at /builds/PumasAI/PumasSystemImages-jl/.julia/packages/Pumas/XHa5R/src/estimation/likelihoods.jl:2625 got unsupported keyword arguments "nsamples", "nadapts"
  ...
Stacktrace:
 [1] kwerr(::NamedTuple{(:constantcoef, :nsamples, :nadapts), Tuple{NamedTuple{(:ωKa, :ωCL, :ωVc), Tuple{Float64, Float64, Float64}}, Int64, Int64}}, ::Function, ::PumasModel{ParamSet{NamedTuple{(:θ, :ωKa, :ωCL, :ωVc, :σ), Tuple{Constrained{DiagNormal, VectorDomain{Vector{Float64}, Vector{Float64}, Vector{Float64}}}, Gamma{Float64}, Gamma{Float64}, Gamma{Float64}, Gamma{Float64}}}}, var"#5#19", var"#6#20", var"#8#22", var"#10#24", Depots1Central1, var"#11#25", var"#15#29", Nothing}, ::Vector{Subject{NamedTuple{(:DV,), Tuple{Vector{Union{Missing, Float64}}}}, Pumas.ConstantCovar{NamedTuple{(:WT,), Tuple{Float64}}}, Vector{Pumas.Event{Float64, Float64, Float64, Float64, Float64, Float64, Int64}}, Vector{Float64}}}, ::NamedTuple{(:θ, :ωKa, :ωCL, :ωVc, :σ), Tuple{Vector{Float64}, Float64, Float64, Float64, Float64}}, ::Pumas.BayesMCMC)
   @ Base ./error.jl:157
 [2] top-level scope
   @ REPL[46]:1

Thank you very much for your help!
Tien Nguyen.

2

This is currently not supported but it will be in the next release. For now, please just hard code the values in the model and redefine the model.

3

Thank you for answering, but unfortunately, I don’t understand what exactly you mean:

For now, please just hard code the values in the model and redefine the model.

I tried this code, but it did not work because I didn’t define ω in @param block. I don’t know how to fix it in @model code.
Thank you for your help!

theopmodel_bayes = @model begin
    @param begin
      # Mode at [2.0, 0.2, 0.8]
      θ ~ Constrained(
        MvNormal(
          [2.0, 0.2, 0.8],
          Diagonal(ones(3))
        ),
        lower = zeros(3),
        upper = fill(10.0, 3),
        init  = [2.0, 0.2, 0.8])
    end
  
    @random begin
      ηKa ~ Normal(0.0, ωKa)
      ηCL ~ Normal(0.0, ωCL)
      ηVc ~ Normal(0.0, ωVc)
    end
  
    @pre begin
      ωKa = 0.2 
      ωCL = 0.2
      ωVc = 0.2 
      σ = 0.5 

      Ka = θ[1]        *exp(ηKa)
      CL = θ[2]*(WT/70)*exp(ηCL)
      Vc = θ[3]        *exp(ηVc)
    end
  
    @covariates WT
  
    @dynamics Depots1Central1
  
    @derived begin
      μ := @. Central / Vc
      DV ~ @. Normal(μ, σ)    
    end
  end
  
  
# Read data 
theo_sd = CSV.read("theo_sd.csv", DataFrame)
data = read_pumas(
  theo_sd;
  id = :ID, time = :TIME, amt = :AMT, observations = [:DV], evid = :EVID, cmt = :CMT,
  covariates = [:WT]
)

param_v2 = init_params(theopmodel_bayes)

result_v2 = fit(theopmodel_bayes, data, param_v2, Pumas.BayesMCMC(), 
  #constantcoef = (ωKa = 0.2, ωCL = 0.2, ωVc = 0.2),
  nsamples = 2000, nadapts = 1000)

ERROR: UndefVarError: ωKa not defined

And this code also didn’t work!

theopmodel_bayes = @model begin
    @param begin
      # Mode at [2.0, 0.2, 0.8]
      θ ~ Constrained(
        MvNormal(
          [2.0, 0.2, 0.8],
          Diagonal(ones(3))
        ),
        lower = zeros(3),
        upper = fill(10.0, 3),
        init  = [2.0, 0.2, 0.8])
  
      # Mean at 0.2 and positive density at 0.0
      ωKa ~ Gamma(1.0, 0.2)
      ωCL ~ Gamma(1.0, 0.2)
      ωVc ~ Gamma(1.0, 0.2)
  
      # Mean at 0.5 and positive density at 0.0
      σ ~ Gamma(1.0, 0.5)
    end
  
    @random begin
      ηKa ~ Normal(0.0, ωKa)
      ηCL ~ Normal(0.0, ωCL)
      ηVc ~ Normal(0.0, ωVc)
    end
  
    @pre begin
      ωKa = 0.2 #~ Gamma(1.0, 0.2)
      ωCL = 0.2 #~ Gamma(1.0, 0.2)
      ωVc = 0.2 #~ Gamma(1.0, 0.2)
      σ = 0.5 #~ Gamma(1.0, 0.5)

      Ka = θ[1]        *exp(ηKa)
      CL = θ[2]*(WT/70)*exp(ηCL)
      Vc = θ[3]        *exp(ηVc)
    end
  
    @covariates WT
  
    @dynamics Depots1Central1
  
    @derived begin
      # The conditional mean
      μ := @. Central / Vc
      # Additive error model
      DV ~ @. Normal(μ, σ)    
    end
  end

ERROR: Variable ωKa already defined
4

Replace the above in the first attempt by:

@random begin
   ηKa ~ Normal(0.0, 0.2)
   ηCL ~ Normal(0.0, 0.2)
   ηVc ~ Normal(0.0, 0.2)
end

and remove the ω definitions in the @pre block. Then do the same thing for σ. This is often called “hard coding” the values in the model.

Additionally, since Pumas gives you access to a whole programming language (Julia), one can go a step further and define a function that outputs models with specific values fixed. For example, I can use the @model macro inside a function as follows:

function define_model(ωKa, ωCL, ωVc, σ)
  return @model begin
      @param begin
        # Mode at [2.0, 0.2, 0.8]
        θ ~ Constrained(
          MvNormal(
            [2.0, 0.2, 0.8],
            Diagonal(ones(3))
          ),
          lower = zeros(3),
          upper = fill(10.0, 3),
          init  = [2.0, 0.2, 0.8])
      end
  
      @random begin
        ηKa ~ Normal(0.0, ωKa)
        ηCL ~ Normal(0.0, ωCL)
        ηVc ~ Normal(0.0, ωVc)
      end
  
      @pre begin
        Ka = θ[1]        *exp(ηKa)
        CL = θ[2]*(WT/70)*exp(ηCL)
        Vc = θ[3]        *exp(ηVc)
      end
  
      @covariates WT
  
      @dynamics Depots1Central1
  
      @derived begin
        # The conditional mean
        μ := @. Central / Vc
        # Additive error model
        DV ~ @. Normal(μ, σ)    
      end
  end
end

Now I can use this function to define a model as follows:

theopmodel_bayes = define_model(0.2, 0.2, 0.2, 0.5)

and the values of the ωs and σ will come from the inputs to the function.

5

Oh, thank you very much for kindly support.

Tien Nguyen.

1 Like