# How to handle several random effects with pairwise correlation between some

**URL:** <https://discourse.pumas.ai/t/how-to-handle-several-random-effects-with-pairwise-correlation-between-some/933>\
**Category:** How-to\
**Created:** [September 17, 2023, 7:54pm UTC](https://discourse.pumas.ai/t/how-to-handle-several-random-effects-with-pairwise-correlation-between-some/933 "2023-09-17T19:54:47Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![patrick](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.pumas.ai/patrick/32/206_2.png) [@patrick](https://discourse.pumas.ai/u/patrick)\
**Post date:** [September 17, 2023, 7:54pm UTC](https://discourse.pumas.ai/t/how-to-handle-several-random-effects-with-pairwise-correlation-between-some/933/1 "2023-09-17T19:54:47Z")

</div>

Sometimes it is useful to have some random effects correlated and others uncorrelated with each other. In Pumas, this is possible by specifying separate random effects for those that are correlated and those that are not correlated. In the below example, we have that the clearance and volume of distribution parameters are correlated through the first and second component of the `η_vc_cl` random effect, while the other random effect on the rate of absorption `η_ka` is uncorrelated with the other two.

```auto
model = @model begin
    @param begin
        tvvc ∈ RealDomain(lower = 0)
        tvcl ∈ RealDomain(lower = 0)
        tvka ∈ RealDomain(lower = 0)
        Ωvc_cl ∈ PSDDomain(2)
        ωka ∈ RealDomain(lower = 1e-5)
        σ_pk ∈ RealDomain(lower = 1e-5)
    end
    @random begin
        η_vc_cl ~ MvNormal(Ωvc_cl)
        η_ka ~ Normal(0, ωka)
    end
    @pre begin
        Ka = tvka*exp(η_ka)
        CL = tvcl*exp(η_vc_cl[1])
        Vc = tvvc*exp(η_vc_cl[2])
    end
    @dynamics Depots1Central1
    @derived begin
        conc_pk ~ @. LogNormal(log(Central/Vc), σ_pk)
    end
end

doses = DosageRegimen(10, addl = 3, ii = 12)

raw_subject = Subject(id="1", events = doses)

param = (tvvc = 2.1, tvcl = 2.1, tvka = 0.9, Ωvc_cl = [0.4 0.2; 0.2 0.4], ωka = 0.3, σ_pk = 0.2)

subject = Subject(simobs(model, raw_subject, param; obstimes = 0:0.5:72))

using PumasUtilities
sim_plot(subject)

```

The simulation gives a trajectory as follows

 ![Skærmbillede 2023-09-17 kl. 21.30.02](https://canada1.discourse-cdn.com/flex030/uploads/pumas/original/1X/d13e7765d1c585b0fece9ded7f32d66e276911eb.png)

To get a confirmation that the two components are indeed correlated, we can simulate many subjects as follows

```auto
raw_population = [Subject(id=string(id), events = doses) for id = 1:1000]
sim_population = simobs(model, raw_population, param; obstimes = 0:0.5:72)
population = Subject.(sim_population)
sim_plot(population)

```

and this gives trajectories as follows

 ![Skærmbillede 2023-09-17 kl. 21.32.42](https://canada1.discourse-cdn.com/flex030/uploads/pumas/original/1X/a6011585146d951d39a37c15797b1843811bfd8c.png)

Then, we can look at their tabular representation using the `DataFrame` constructor

```auto
sim_population_df = DataFrame(sim_population)
# First the first row of each subject
sim_popuplation_df_first = unique(sim_population_df, [:id])

```

and from this we can look at correlations between the different random effects using an AlgebraOfGraphics plot

```auto
draw(data(sim_population_df) * (visual(AlgebraOfGraphics.Scatter) + AlgebraOfGraphics.linear()*visual(color=:red) * mapping([:η_vc_cl_1, :η_vc_cl_2,:η_ka_1], [:η_vc_cl_1 :η_vc_cl_2 :η_ka_1], col = dims(1), row = dims(2)))

```

 ![Skærmbillede 2023-09-17 kl. 21.48.40](https://canada1.discourse-cdn.com/flex030/uploads/pumas/original/1X/e64ed148d381a871e9b347c8967cd52ae830dd05.png)  
and we see that the correlations are as expected.

This example had a pair with correlations and a third random effect that was uncorrelated with the other two, but the scalar random effect for rate of absorption that had a `Normal` distribution could just as well have been multivariate, or we could have had ever more random effects each with their own correlation structure between the elements such as the following pseudo-example

```auto
@model begin
    @param begin
        [...]
         Ω1 ∈ PSDDomain(2)
         Ω2 ∈ PSDDomain(4)
         Ω3 ∈ PSDDomain(3)
        [...]
    end
    @random begin
        η1 ~ MvNormal(Ω1)
        η2 ~ MvNormal(Ω2)
        η3 ~ MvNormal(Ω3)
    end
    [...]
end

```

Hope this helps!