Feature request: Or is it already there?

Basic Usage
1

Hi folks,

I was thinking about this implementation to sample $\epsilon$ with multivariate distribution. Similar to the way we handle variables in random block. Here is an example: (See comments below the code):

const N_THETA = 3 
const N_OMEGA = 3 
const N_SIGMA = 2

model = @model begin
    @param begin
        θ ∈ VectorDomain(N_THETA)
        Ω ∈ PDiagDomain(N_OMEGA)
        Σ ∈ PDiagDomain(N_SIGMA)
    end
    @random begin
        η ~ MvNormal(Ω)
        ϵ ~ MvNormal(Σ)
    end
    @pre begin
        CLᵢ = exp(log(θ[1]) + η[1])
        Vᵢ = exp(log(θ[2]) + η[2])  
        kaᵢ = exp(log(θ[3]) + η[3])
        keᵢ = CLᵢ/Vᵢ 
    end
    @dynamics begin
        y1' = - kaᵢ * y1
        y2' = kaᵢ* y1 - keᵢ * y2 
    end
    @derived begin
        cp = @. y2/Vᵢ 
        dv = @. cp + cp * ϵ[1] + ϵ[2]
    end
end

COMMENT: Currently, epsilon will have one value per subject and will not change which is expected from @random block. Is there a way that one can sample within the same subject?

Here is a dummy setup:

event = DosageRegimen(100)
sub = Subject(;id=1, events=event)
par_values = (θ=[8, 25, 1.5], Ω = Diagonal([0.09, 0.09, 0.09]), Σ = Diagonal([0.015, 0.0625]), )

sim = simobs(model, sub, par_values) |> DataFrame

You can see the values of epsilons from dataframe.

2

Before trying to come up with a technical solution, I would be curious to know more about what exactly you’re trying to do here. Can you explain the motivation and what kind of properties you’re looking for? This way of writing it would mean that dv would not carry a distribution and thus it’s excluded as a possible dependent variable in estimation for example.

3

Hi Patrick -

The ultimate goal is to be able to implement Berkson term e.g., add variability term to independent variables such as time, dose, and others that change within the subject. Using this implementation can could give me more flexibility to work mathematically with error terms. I also would like to look at these sampled epsilons in the dataframe. Does any of that make sense?

Best,
Mutaz

4

I think I understand where you’re trying to go, but I also think it would be nice to maybe do a brainstorming/write-up to get a better feel of what is actually needed and if it will work as expected/correctly if done as shown at the top. Variability to time (I suppose you mean that it’s not clear when the observation as actually made?) could be hard/impossible to do with the current interface for example. I could show you how to add time-varying per subject epsilons, but it would not work as you may expect in estimation for example.