# Defining different error models for data from different sources
**URL:** https://discourse.pumas.ai/t/defining-different-error-models-for-data-from-different-sources/838
**Category:** Estimation
**Created:** [February 24, 2023, 4:27am UTC](https://discourse.pumas.ai/t/defining-different-error-models-for-data-from-different-sources/838 "2023-02-24T04:27:54Z")
**Posts on this page:** 3
**Page:** 1
### Author:  [@ahmed.salem](https://discourse.pumas.ai/u/ahmed.salem)
#### Post date: [February 24, 2023, 4:27am UTC](https://discourse.pumas.ai/t/defining-different-error-models-for-data-from-different-sources/838/1 "2023-02-24T04:27:54Z")
Hello,
I have datasets collected from multiple clinical sites and I would like to define different error model based on the clinical site. I wrote the code as following
```auto
mod = @model begin
@param begin
TVCL ∈ RealDomain(lower= 0.0001, upper = 6)
tvKTV ∈ RealDomain(lower = 0.001, upper = 2)
TVVc ∈ RealDomain(lower = 0.001, upper = 20)
ecmo ∈ RealDomain(lower = -3,upper = 3)
cl ∈ RealDomain(lower = -3,upper = 3)
v ∈ RealDomain(lower = -3,upper = 3)
circuit ∈ RealDomain(lower = 0.001, upper = 5)
Ω ∈ PDiagDomain(3)
σ_prop_1 ∈ RealDomain(lower = 0.00001, upper = 1)
σ_prop_2 ∈ RealDomain(lower = 0.00001, upper = 1)
σ_add ∈ RealDomain(lower = 1)
end
@random begin
η ~ MvNormal(Ω)
end
@covariates WT ECMO_DAYS IS_CIRCUIT_CHANGE site
@pre begin
_site = site
KTV = tvKTV*exp(η[1])*(ECMO_DAYS/6)^ecmo
CL = TVCL * exp(η[2])*(1-exp(-KTV*t))*(WT/10)^cl*circuit^IS_CIRCUIT_CHANGE
Vc = TVVc * exp(η[3])*(WT/10)^v
end
@dynamics begin
Central' = -(CL/Vc)*Central
end
@derived begin
CP = @. Central / Vc
CONC = if _site ==2
@. Normal(CP, sqrt(CP^2*σ_prop_1 + σ_add))
else
@. Normal(CP, sqrt(CP^2*σ_prop_2))
end
end
end
```
When I run the `fit` function I get the following error
```auto
julia> modfit = fit(mod, estimation, param, Pumas.FOCEI())
[ Info: Checking the initial parameter values.
[ Info: The initial negative log likelihood and its gradient are finite. Check passed.
Iter Function value Gradient norm
0 1.897660e+04 1.198511e+05
* time: 3.0040740966796875e-5
ERROR: gradient is exactly zero in σ_prop_1. This indicates that σ_prop_1 isn't identified.
Stacktrace:
[1] _unidentified_exception(#unused#::Float64, i::Int64, j::Int64, d::Int64, k::Symbol)
@ Pumas /build/_work/PumasSystemImages/PumasSystemImages/julia_depot/packages/Pumas/Td3Jp/src/estimation/likelihoods.jl:3484
[2] _check_zero_gradient(vparam::Vector{Float64}, g::Vector{Float64}, fixedtrf::TransformVariables.TransformTuple{NamedTuple{(:TVCL, :tvKTV, :TVVc, :ecmo, :cl, :v, :circuit, :Ω, :σ_prop_1, :σ_prop_2, :σ_add), Tuple{TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Float64}, Pumas.PDiagTransform, TransformVariables.ShiftedExp{true, Float64}, TransformVariables.Identity, TransformVariables.ShiftedExp{true, Float64}}}})
@ Pumas /build/_work/PumasSystemImages/PumasSystemImages/julia_depot/packages/Pumas/Td3Jp/src/estimation/likelihoods.jl:3475
[3] (::Pumas.var"#614#616"{FOCE, Bool, NLSolversBase.OnceDifferentiable{Float64, Vector{Float64}, Vector{Float64}}, TransformVariables.TransformTuple{NamedTuple{(:TVCL, :tvKTV, :TVVc, :ecmo, :cl, :v, :circuit, :Ω, :σ_prop_1, :σ_prop_2, :σ_add), Tuple{TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Float64}, Pumas.PDiagTransform, TransformVariables.ShiftedExp{true, Float64}, TransformVariables.Identity, TransformVariables.ShiftedExp{true, Float64}}}}, Vector{Vector{Float64}}, Vector{Vector{Float64}}})(state::Vector{Optim.OptimizationState{Float64, Optim.BFGS{LineSearches.InitialStatic{Float64}, LineSearches.BackTracking{Float64, Int64}, Nothing, Float64, Optim.Flat}}})
@ Pumas /build/_work/PumasSystemImages/PumasSystemImages/julia_depot/packages/Pumas/Td3Jp/src/estimation/likelihoods.jl:4047
[4] #585
@ /build/_work/PumasSystemImages/PumasSystemImages/julia_depot/packages/Pumas/Td3Jp/src/estimation/likelihoods.jl:3000 [inlined]
[5] update!(tr::Vector{Optim.OptimizationState{Float64, Optim.BFGS{LineSearches.InitialStatic{Float64}, LineSearches.BackTracking{Float64, Int64}, Nothing, Float64, Optim.Flat}}}, iteration::Int64, f_x::Float64, grnorm::Float64, dt::Dict{Any, Any}, store_trace::Bool, show_trace::Bool, show_every::Int64, callback::Pumas.var"#585#587"{Pumas.DefaultOptimizeFN{Nothing, NamedTuple{(:show_trace, :store_trace, :extended_trace, :g_tol, :allow_f_increases), Tuple{Bool, Bool, Bool, Float64, Bool}}}, Pumas.var"#614#616"{FOCE, Bool, NLSolversBase.OnceDifferentiable{Float64, Vector{Float64}, Vector{Float64}}, TransformVariables.TransformTuple{NamedTuple{(:TVCL, :tvKTV, :TVVc, :ecmo, :cl, :v, :circuit, :Ω, :σ_prop_1, :σ_prop_2, :σ_add), Tuple{TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Float64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Int64}, TransformVariables.ScaledShiftedLogistic{Float64}, Pumas.PDiagTransform, TransformVariables.ShiftedExp{true, Float64}, TransformVariables.Identity, TransformVariables.ShiftedExp{true, Float64}}}}, Vector{Vector{Float64}}, Vector{Vector{Float64}}}}, trace_simplex::Bool)
```
When I fit each dataset separately, the `fit` function works fine. I am not sure if the syntax in the `@derived` block should be modified. Any help would be appreciated.
Thanks
---
### Author:  [@vijay](https://discourse.pumas.ai/u/vijay)
#### Post date: [February 24, 2023, 5:30am UTC](https://discourse.pumas.ai/t/defining-different-error-models-for-data-from-different-sources/838/2 "2023-02-24T05:30:10Z")
1. What’s your initial estimate for sigma\_prop\_1?
2. You have an upper bound of 1 on that parameter, is that intentional?
3. Could you rewrite the model by constructing the covariate effect in the pre block and then passing it into the derived block? e.g. like this
`siteeffect = site==1 ? sigma_prop_1 : sigma_prop_2`
Then pass that `siteeffect` into the derived block
Best,
Vijay
---
### Author:  [@ahmed.salem](https://discourse.pumas.ai/u/ahmed.salem)
#### Post date: [February 24, 2023, 4:32pm UTC](https://discourse.pumas.ai/t/defining-different-error-models-for-data-from-different-sources/838/3 "2023-02-24T16:32:04Z")
> [@vijay](#):
>
> siteeffect = site==1 ? sigma\_prop\_1 : sigma\_prop\_2
Hi @vijay,
Thanks so much for your input, it was very helpful. Here is the modified code if a future user wants to implement it
```auto
mod = @model begin
@param begin
TVCL ∈ RealDomain(lower= 0.0001, upper = 6)
tvKTV ∈ RealDomain(lower = 0.001, upper = 2)
TVVc ∈ RealDomain(lower = 0.001, upper = 20)
ecmo ∈ RealDomain(lower = -3,upper = 3)
cl ∈ RealDomain(lower = -3,upper = 3)
v ∈ RealDomain(lower = -3,upper = 3)
circuit ∈ RealDomain(lower = 0.001, upper = 5)
Ω ∈ PDiagDomain(3)
σ_prop_1 ∈ RealDomain(lower = 0.00001)
σ_prop_2 ∈ RealDomain(lower = 0.00001)
σ_add ∈ RealDomain(lower = 1)
end
@random begin
η ~ MvNormal(Ω)
end
@covariates WT ECMO_DAYS IS_CIRCUIT_CHANGE site
@pre begin
siteeffect_1 = if site==2
σ_prop_1
else
σ_prop_2
end
siteeffect_2 = site==2 ? σ_add : 0.0
KTV = tvKTV*exp(η[1])*(ECMO_DAYS/6)^ecmo
CL = TVCL * exp(η[2])*(1-exp(-KTV*t))*(WT/10)^cl*circuit^IS_CIRCUIT_CHANGE
Vc = TVVc * exp(η[3])*(WT/10)^v
end
@dynamics begin
Central' = -(CL/Vc)*Central
end
@derived begin
CP = @. Central / Vc
CONC = @. Normal(CP, sqrt(CP^2*siteeffect_1 + siteeffect_2))
end
end
```