# Unable to specify a correct prior distribution

**URL:** <https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371>\
**Category:** bayesian\
**Created:** [July 14, 2020, 8:08pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371 "2020-07-14T20:08:41Z")\
**Posts on this page:** 16\
**Page:** 1

<div class="post-metadata">

**Author:** ![anqipan](https://avatars.discourse-cdn.com/v4/letter/a/ecae2f/32.png) [@anqipan](https://discourse.pumas.ai/u/anqipan)\
**Post date:** [July 14, 2020, 8:08pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/1 "2020-07-14T20:08:41Z")

</div>

Hi,  
I am having difficulty specifying a prior distribution in Pumas with lower bound and initials, my codes is:

```auto
    @param begin

      θ ~ MvNormal([16.6,13.1,4.05,4.89,0.623,0,0,0],
      Matrix{Float64}([0.0895	0.0167	0.00318	0.0033	0.000907	0	0	0;
                       0.0167	0.0655	-0.00765	0.00333	-0.000917	0	0	0;
                       0.00318	-0.00765	0.0238	0.0147	-0.000225	0	0	0;
                       0.0033	0.00333	0.0147	0.0161	-0.000243	0	0	0;
                       0.000907	-0.000917	-0.000225	-0.000243	0.0016	0	0	0;
                       0	0	0	0	0	1000000	0	0;
                       0	0	0	0	0	0	1000000	0;
                       0	0	0	0	0	0	0	1000000]),lower=[0,0,0,0,-Inf,-Inf,-Inf,-Inf],init=[18,14.2,2.13,4.29,1,1,1,1])

      Ω ~ Wishart(3, Matrix{Float64}([0.282	0	0; 0	0.149	0; 0	0	0.0377]))

      σ_prop ∈ RealDomain(lower = 0)

    end

```

I mainly had 2 error message:  
1.

```auto
LoadError: MethodError: no method matching MvNormal(::Array{Float64,1}, ::Array{Float64,2}; lower=[0.0, 0.0, 0.0, 0.0, -Inf, -Inf, -Inf, -Inf])

```

1. 

```auto
LoadError: ArgumentError: no Domain(::Wishart{Float64,PDMats.PDMat{Float64,Array{Float64,2}}}) constructor defined

```

I have google both of them and didn’t find anything matched…

Basically I am trying to specify a multivariate normal distribution for the Thetas with lower bound and initials, a Wishart distribution for the Omegas with initials.

May I have any recommendations on which part I missed in my codes?

Thanks!  
Anqi

---

<div class="post-metadata">

**Author:** ![anqipan](https://avatars.discourse-cdn.com/v4/letter/a/ecae2f/32.png) [@anqipan](https://discourse.pumas.ai/u/anqipan)\
**Post date:** [July 17, 2020, 4:51pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/2 "2020-07-17T16:51:47Z")

</div>

Hello,  
Just follow up with the complete coding for this topic:

```auto
@time using Distributions, Pumas, CSV, TableView, StatsPlots, DataFrames, LinearAlgebra

ped_mod = @model begin

    @param begin

      θ ~ Constrained(MvNormal([16.6,13.1,4.05,4.89,0.623,0,0,0],
      Matrix{Float64}([0.0895	0.0167	0.00318	0.0033	0.000907	0	0	0;
                       0.0167	0.0655	-0.00765	0.00333	-0.000917	0	0	0;
                       0.00318	-0.00765	0.0238	0.0147	-0.000225	0	0	0;
                       0.0033	0.00333	0.0147	0.0161	-0.000243	0	0	0;
                       0.000907	-0.000917	-0.000225	-0.000243	0.0016	0	0	0;
                       0	0	0	0	0	1000000	0	0;
                       0	0	0	0	0	0	1000000	0;
                       0	0	0	0	0	0	0	1000000])),lower=[0,0,0,0,-Inf,-Inf,-Inf,-Inf])

      Ω ~ InverseWishart(3, diagm([0.282, 0.149, 0.0377]))

      σ_prop ∈ RealDomain(lower = 0)

    end

    @random begin
      η ~ MvNormal(Ω)
    end

    @covariates WTKG CRCL

    @pre begin
    #--INSERT COVARIATE EFFECTS
    COV1=(CRCL/100)^θ[5]
    COV2=(WTKG/70)^θ[6]
    COV3=(WTKG/70)^θ[7]
    COV4=(WTKG/70)^θ[8]

    TVCLI = θ[1]*COV1*COV4
    TVCL = TVCLI

    TVV1I = θ[2]*COV2
    TVV1 = TVV1I

    TVQI = θ[3]
    TVQ = TVQI

    TVV2I = θ[4]*COV3
    TVV2 = TVV2I

    CL = TVCL*exp(η[1])
    V1 = TVV1*exp(η[2])
    Q = TVQ
    V2 = TVV2*exp(η[3])

    Vc = V1
    Vp = V2

    S1 = V1

    #CALCULATION OF SECONDARY PARAMETERS
    KE = CL/V1
    K12 = Q/V1
    K21 = Q/V2
    AA = KE+K12+K21
    ALPH = (AA+sqrt(AA*AA-4*KE*K21))/2
    BETA = (AA-sqrt(AA*AA-4*KE*K21))/2

    end

    @dynamics Central1Periph1 #a two compartment model

    @derived begin
        cp = @. 1000*(Central/Vc)
        DV ~ @. Normal(cp, sqrt((cp^2*σ_prop)))
      end
  end

inputDataset = CSV.read("dat.csv")
df = read_pumas(inputDataset, id = :ID, dvs =[:DV], cvs=[:WTKG, :CRCL], evid=:EVID, amt=:AMT,cmt=:CMT, rate=:RATE, time=:TIME)

param = (
    θ=([18,14.2,2.13,4.29,1,1,1,1]),
    Ω = Diagonal([0.25,0.25,0.25]),
    σ_prop = 0.04)

 pkres = fit(ped_mod, df, param, Pumas.FOCEI(),
              optimize_fn=Pumas.DefaultOptimizeFN(show_trace=true, extended_trace=false))

```

The error is saying no method matching inverse:

```auto
LoadError: MethodError: no method matching inverse!(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Pumas.PSDTransform, ::Diagonal{Float64,Array{Float64,1}})

Closest candidates are:

inverse!(::AbstractArray{T,1} where T<:Real, !Matched::TransformVariables.ArrayTransform, ::AbstractArray) at C:\Users\.juliapro\JuliaPro_v1.4.2-2\packages\TransformVariables\a4AMY\src\aggregation.jl:79

inverse!(::AbstractArray, !Matched::Pumas.ElementArrayTransform, ::AbstractArray) at C:\Users\.juliapro\JuliaPro_v1.4.2-2\packages\Pumas\NCmSe\src\estimation\transforms.jl:25

inverse!(::AbstractArray, !Matched::Pumas.ConstantTransform, ::Any) at C:\Users\juliapro\JuliaPro_v1.4.2-2\packages\Pumas\NCmSe\src\estimation\transforms.jl:54

...

in expression starting at C:\Users\118

_inverse!_tuple(::Array{Float64,1}, ::Tuple{Pumas.ElementArrayTransform{TransformVariables.ShiftedExp{true,Float64},1},Pumas.PSDTransform,TransformVariables.ShiftedExp{true,Int64}}, ::Tuple{Array{Float64,1},Diagonal{Float64,Array{Float64,1}},Float64}) at aggregation.jl:196

inverse!(::Array{Float64,1}, ::TransformVariables.TransformTuple{NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Pumas.ElementArrayTransform{TransformVariables.ShiftedExp{true,Float64},1},Pumas.PSDTransform,TransformVariables.ShiftedExp{true,Int64}}}}, ::NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Array{Float64,1},Diagonal{Float64,Array{Float64,1}},Float64}}) at aggregation.jl:237

inverse(::TransformVariables.TransformTuple{NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Pumas.ElementArrayTransform{TransformVariables.ShiftedExp{true,Float64},1},Pumas.PSDTransform,TransformVariables.ShiftedExp{true,Int64}}}}, ::NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Array{Float64,1},Diagonal{Float64,Array{Float64,1}},Float64}}) at generic.jl:206

fit(::PumasModel{ParamSet{NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Constrained{MvNormal{Float64,PDMats.PDMat{Float64,Array{Float64,2}},Array{Float64,1}},VectorDomain{Array{Float64,1},Array{TransformVariables.Infinity{true},1},Array{Float64,1}}},InverseWishart{Float64,PDMats.PDMat{Float64,Array{Float64,2}}},RealDomain{Int64,TransformVariables.Infinity{true},Float64}}}},var"#548#603",var"#549#604",var"#550#606",Central1Periph1,var"#551#607",var"#577#633"}, ::Array{Subject{NamedTuple{(:DV,),Tuple{Array{Union{Missing, Float64},1}}},Pumas.ConstantCovar{NamedTuple{(:WTKG, :CRCL),Tuple{Float64,Float64}}},Array{Pumas.Event{Float64,Float64,Float64,Float64,Float64,Float64,Int64},1},Array{Float64,1},Float64},1}, ::NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Array{Float64,1},Diagonal{Float64,Array{Float64,1}},Float64}}, ::Pumas.FOCEI; optimize_fn::Pumas.DefaultOptimizeFN{Nothing,NamedTuple{(:show_trace, :store_trace, :extended_trace, :g_tol, :allow_f_increases),Tuple{Bool,Bool,Bool,Float64,Bool}}}, constantcoef::NamedTuple{,Tuple{}}, omegas::Tuple{}, ensemblealg::EnsembleSerial, kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{,Tuple{}}}) at likelihoods.jl:1330

(::StatsBase.var"#fit##kw")(::NamedTuple{(:optimize_fn,),Tuple{Pumas.DefaultOptimizeFN{Nothing,NamedTuple{(:show_trace, :store_trace, :extended_trace, :g_tol, :allow_f_increases),Tuple{Bool,Bool,Bool,Float64,Bool}}}}}, ::typeof(fit), ::PumasModel{ParamSet{NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Constrained{MvNormal{Float64,PDMats.PDMat{Float64,Array{Float64,2}},Array{Float64,1}},VectorDomain{Array{Float64,1},Array{TransformVariables.Infinity{true},1},Array{Float64,1}}},InverseWishart{Float64,PDMats.PDMat{Float64,Array{Float64,2}}},RealDomain{Int64,TransformVariables.Infinity{true},Float64}}}},var"#548#603",var"#549#604",var"#550#606",Central1Periph1,var"#551#607",var"#577#633"}, ::Array{Subject{NamedTuple{(:DV,),Tuple{Array{Union{Missing, Float64},1}}},Pumas.ConstantCovar{NamedTuple{(:WTKG, :CRCL),Tuple{Float64,Float64}}},Array{Pumas.Event{Float64,Float64,Float64,Float64,Float64,Float64,Int64},1},Array{Float64,1},Float64},1}, ::NamedTuple{(:θ, :Ω, :σ_prop),Tuple{Array{Float64,1},Diagonal{Float64,Array{Float64,1}},Float64}}, ::Pumas.FOCEI) at likelihoods.jl:1328

top-level scope at tp4_5.jl:118

```

Thanks!  
Anqi

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/andreasnoack/32/22_2.png) [@andreasnoack](https://discourse.pumas.ai/u/andreasnoack)\
**Post date:** [July 20, 2020, 7:09pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/3 "2020-07-20T19:09:09Z")

</div>

Error 2. in your first post is a missing feature in Pumas but it’s tiny so we might be able to consider it a bugfix and include it Pumas 1.0. Right now, we only support the `InverseWishart` but extending the support to cover `Wishart` should be very easy.

The problem in your followup post is `Ω = Diagonal([0.25,0.25,0.25])` because it restricts `Ω` do a diagonal matrix but the `InverseWishart` isn’t restricted to diagonal matrices so you should be able to make this work by setting `Ω = diagm([0.25,0.25,0.25])` instead. We might be able to make your version to work as well. I’ve filed an issue.

---

<div class="post-metadata">

**Author:** ![anqipan](https://avatars.discourse-cdn.com/v4/letter/a/ecae2f/32.png) [@anqipan](https://discourse.pumas.ai/u/anqipan)\
**Post date:** [July 20, 2020, 8:00pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/4 "2020-07-20T20:00:17Z")

</div>

Thank you! The error is fixed after I changed the initial to a diagm.

And may I just confirm that, earlier this morning you mentioned the prior is just a domain(at least for the current version), which sounds it won’t enter the posterior distribution like a real prior.  
Is this the case for both FOCE and HMC-NUTS?  
Or the distributional domain(prior) does enter the posterior as a real prior in HMC-NUTS?

Thanks!  
Anqi

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/andreasnoack/32/22_2.png) [@andreasnoack](https://discourse.pumas.ai/u/andreasnoack)\
**Post date:** [July 21, 2020, 6:55am UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/5 "2020-07-21T06:55:59Z")

</div>

> And may I just confirm that, earlier this morning you mentioned the prior is just a domain(at least for the current version), which sounds it won’t enter the posterior distribution like a real prior.  
> Is this the case for both FOCE and HMC-NUTS?

Yes. That is indeed the case for FOCE(I) (and `LaplaceI`) but for the NUTS sampler aka `BayesMCMC` the likelihood contribution of the prior is applied as you would expect.

---

<div class="post-metadata">

**Author:** ![anqipan](https://avatars.discourse-cdn.com/v4/letter/a/ecae2f/32.png) [@anqipan](https://discourse.pumas.ai/u/anqipan)\
**Post date:** [July 21, 2020, 4:20pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/6 "2020-07-21T16:20:27Z")

</div>

I see, thank you for the confirmation!

---

<div class="post-metadata">

**Author:** ![anqipan](https://avatars.discourse-cdn.com/v4/letter/a/ecae2f/32.png) [@anqipan](https://discourse.pumas.ai/u/anqipan)\
**Post date:** [July 22, 2020, 6:19pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/7 "2020-07-22T18:19:48Z")

</div>

I do have a further question regard the same model. Is there anyway that I could fix the off-diagonal elements of OMEGA matrix to 0?  
I tried to use Diagonal() to force it as diagonal matrix both in @param and the initials, but it doesn’t run through.

Thank you!  
Anqi

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/andreasnoack/32/22_2.png) [@andreasnoack](https://discourse.pumas.ai/u/andreasnoack)\
**Post date:** [July 22, 2020, 8:49pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/8 "2020-07-22T20:49:57Z")

</div>

I think the only way to do that would be to specify each of the diagonal elements separately, i.e. something like

```auto
@random begin
  ηCL ~ InverseGamma(...)
  ηV1 ~ InverseGamma(...)
  ηV2 ~ InverseGamma(...)
end

```

I chose `InverseGamma` here since it’s equivalent to a scalar `InverseWishart` but it might be better to use the `Gamma` since it has support at `zero`.

---

<div class="post-metadata">

**Author:** ![Vaibhavdixit02](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/vaibhavdixit02/32/86_2.png) [@Vaibhavdixit02](https://discourse.pumas.ai/u/Vaibhavdixit02)\
**Post date:** [August 14, 2020, 1:36pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/9 "2020-08-14T13:36:08Z")

</div>

This should be fixed on the new release @anqipan

---

<div class="post-metadata">

**Author:** ![NamTienNguyen](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/namtiennguyen/32/143_2.png) [@NamTienNguyen](https://discourse.pumas.ai/u/NamTienNguyen)\
**Post date:** [November 17, 2021, 12:03pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/10 "2021-11-17T12:03:22Z")

</div>

@andreasnoack By the way for this issue, I had a confusion in [Bayesian tutorial](https://tutorials.pumas.ai/html/pkpd/bayesianestimation.html) that:

> Since the mode of `InverseWishart` is `Ψ/(ν + p + 1)` , the mode of our prior is `diagm(fill(0.2, 3))` which corresponds to an inter-subject variability peaking at 20%.

In this, the IIV is 20% or `sqrt(0.2)*100` % because the OMEGA is in the matrix.

Thank for your help.

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/andreasnoack/32/22_2.png) [@andreasnoack](https://discourse.pumas.ai/u/andreasnoack)\
**Post date:** [November 17, 2021, 12:27pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/11 "2021-11-17T12:27:09Z")

</div>

Thanks. That is indeed correct. I should probably have been `fill(0.2^2, 3)`. We’ll have it fixed.

---

<div class="post-metadata">

**Author:** ![NamTienNguyen](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/namtiennguyen/32/143_2.png) [@NamTienNguyen](https://discourse.pumas.ai/u/NamTienNguyen)\
**Post date:** [November 17, 2021, 5:14pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/12 "2021-11-17T17:14:42Z")

</div>

@andreasnoack In this case, could the code write that, the InverseGamma() distribution is for Omega or for the ETA?

And because I didn’t set the matrix for omega, then the IIV = 20℅ was represented as in SD = 0.2 rather than 0.2^2 in code block?

Thank you again for help!

```auto
param begin
    # Mode at [2.0, 0.2, 0.8, 2.0]
    θ ~ Constrained(
      MvNormal(
        [2.0, 0.2, 0.8, 2.0],
        Diagonal(ones(4))
      ),
      lower = zeros(4),
      upper = fill(10.0, 4),
      init = [2.0, 0.2, 0.8, 2.0])

    # Mean at 0.2 and positive density at 0.0
    ωKa ~ InverseGamma(1.0, 0.2)
    ωCL ~ InverseGamma(1.0, 0.2)
    ωVc ~ InverseGamma(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

```

---

<div class="post-metadata">

**Author:** ![NamTienNguyen](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/namtiennguyen/32/143_2.png) [@NamTienNguyen](https://discourse.pumas.ai/u/NamTienNguyen)\
**Post date:** [November 18, 2021, 2:47pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/13 "2021-11-18T14:47:13Z")

</div>

Dear `Pumas` team. if I want to impose a diagonal structure on `Ω` then should I use `InverseGama` rather than `Gamma` for each element of `Ω` because it’s equivalent to `InverseWishart` as mentioned above by @andreasnoack?

And here, should I use 0.2 (IIV = 20% for instance) rather than 0.2^2 in the `InverseGamma` because it’s not in the matrix structure?

Thank you very much for helping!

```auto
   ωKa ~ InverseGamma(1.0, 0.2)
    ωCL ~ InverseGamma(1.0, 0.2)
    ωVc ~ InverseGamma(1.0, 0.2)

Then: 
@random begin
    ηKa ~ Normal(0.0, ωKa)
    ηCL ~ Normal(0.0, ωCL)
    ηVc ~ Normal(0.0, ωVc)

```

---

<div class="post-metadata">

**Author:** ![andreasnoack](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/andreasnoack/32/22_2.png) [@andreasnoack](https://discourse.pumas.ai/u/andreasnoack)\
**Post date:** [November 19, 2021, 8:07am UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/14 "2021-11-19T08:07:31Z")

</div>

If you want to impose a diagonal structure on `Ω` then, indeed, you should define each scalar component the way you do it here.

Also, you are right that the `InverseGamma` is the scalar version of the `InverseWishart`.

However, that doesn’t mean that you should use that distribution. Actually, there are [good reasons to avoid it](http://www.stat.columbia.edu/~gelman/research/published/taumain.pdf) so using `Gamma` with shape parameter equal to one might be a better choice since it has positive density at zero. In contrast, the `InverseGamma` forces the `ω` away from zero which might not always be reasonable.

Finally, if you decide to use the `InverseGamma` anyway then you’d have to choose the parameters differently. You can always check the mode and mean easily in Pumas (or Julia) with

```auto
julia> mean(InverseGamma(1.0, 0.2))
Inf

julia> mode(InverseGamma(1.0, 0.2))
0.1

```

I usually, look up the relationship between the parameters and the mean/variance on [Wikipedia](https://en.wikipedia.org/wiki/Inverse-gamma_distribution). Then you pick a mean value and a shape parameter value and then solve for the scale.

```auto
julia> mean(InverseGamma(2.0, 0.2))
0.2

julia> mode(InverseGamma(2.0, 0.2*3))
0.20000000000000004

```

---

<div class="post-metadata">

**Author:** ![storopoli](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/storopoli/32/190_2.png) [@storopoli](https://discourse.pumas.ai/u/storopoli)\
**Post date:** [November 19, 2021, 1:20pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/15 "2021-11-19T13:20:33Z")

</div>

I am new to Pumas but I would like to point out that while using NUTS in a MCMC sampler a much faster and efficient prior would be a LKJ (2009) prior:

This is just a benchmark for random sampling using `Distributions.jl`:

```auto
using Distributions
using LinearAlgebra

function rand_invwishart()
           d = InverseWishart(6, diagm(fill(0.2, 3)))
           return rand(d)
end

function rand_lkj()
           d = LKJ(3, 1)
           return rand(d) * diagm(fill(0.2, 3))
end

function rand_lkj_corr()
           d = LKJCholesky(3, 1)
           A = rand(d)
           return A.L * A.U * diagm(fill(0.2, 3))
end

julia> @btime rand_invwishart()
@btime 1.718 μs (16 allocations: 1.38 KiB)
3×3 Matrix{Float64}:
 0.11615 0.0392431 0.0272846
 0.0392431 0.0425069 0.00645005
 0.0272846 0.00645005 0.0218819

julia> @btime rand_lkj()
  797.143 ns (15 allocations: 1.50 KiB)
3×3 Matrix{Float64}:
  0.2 0.046439 -0.150757
  0.046439 0.2 0.0195246
 -0.150757 0.0195246 0.2

julia> @btime rand_lkj_corr()
  529.849 ns (8 allocations: 944 bytes)
3×3 Matrix{Float64}:
  0.2 -0.0523386 0.136233
 -0.0523386 0.2 0.0646473
  0.136233 0.0646473 0.2

```

This is also a great paper that deals with the several options to specifying priors for covariance/correlation matrices in a Bayesian setting: [New Prairie Press - Conference on Applied Statistics in Agriculture: BAYESIAN INFERENCE FOR A COVARIANCE MATRIX](https://newprairiepress.org/agstatconference/2014/proceedings/8/)

Lewandowski, Daniel, Dorota Kurowicka, and Harry Joe. “Generating random correlation matrices based on vines and extended onion method.” _Journal of multivariate analysis_ 100.9 (2009): 1989-2001.

> One small caveat is that LKJCorr inside a Turing models will have to await for this issue to solve: [Feature request: Add LKJCholesky · Issue #1629 · TuringLang/Turing.jl · GitHub](https://github.com/TuringLang/Turing.jl/issues/1629#issuecomment-962920241)

---

<div class="post-metadata">

**Author:** ![mohamed82008](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/mohamed82008/32/127_2.png) [@mohamed82008](https://discourse.pumas.ai/u/mohamed82008)\
**Post date:** [November 19, 2021, 3:12pm UTC](https://discourse.pumas.ai/t/unable-to-specify-a-correct-prior-distribution/371/16 "2021-11-19T15:12:52Z")

</div>

@storopoli the LKJ prior is currently supported in the released Pumas. Here is an example of its usage:

```auto
theopmodel_bayes = @model begin
  @param begin
    θ ∈ VectorDomain(5,
      lower=[0.1,0.008,0.0004,0.1,0.0001],
      upper=[5,0.5,0.09,5,1.5],
      init=[1.9,0.0781,0.0463,1.5,0.4]
    )
    Ω ~ LKJ(1, 2.0)
    σ² ~ Uniform(0.0, 2*0.388)
  end

  @random begin
    η ~ MvNormal(Ω)
  end

  @covariates SEX WT

  @pre begin
    Ka = (SEX == 1 ? θ[1] : θ[4]) + η[1]
    K = θ[2]
    CL = θ[3]*(WT/70)^θ[5]
    Vc = CL/K
    SC = Vc/(WT/70)
  end

  @vars begin
    conc = Central / SC
  end

  @dynamics Depots1Central1

  @derived begin
    dv ~ @. Normal(conc, sqrt(σ²))
  end
end

```

For those new to Pumas/Julia, to know what each parameter passed to `LKJ` above means, write `?LKJ` in the Pumas/Julia REPL then press enter.
