# MonteCarloMeasurements & MTK, again
**URL:** https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205
**Category:** Modelling & Simulations
**Created:** [November 14, 2023, 12:48pm UTC](https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205 "2023-11-14T12:48:25Z")
**Posts on this page:** 4
**Page:** 1
### Author:  [@BLI](https://discourse.julialang.org/u/BLI)
#### Post date: [November 14, 2023, 12:48pm UTC](https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205/1 "2023-11-14T12:48:25Z")
I’m trying to use MonteCarloMeasurements with MTK. I have some problems wrt. (i) the number of particles, and (ii) interpretation/presentation of the results.
Here is a simple RC circuit example:
```julia
using Plots, LaTeXStrings
using ModelingToolkit, DifferentialEquations
using MonteCarloMeasurements
# 1. Independent variables
@variables t
# 2. Differential operator
Dt = Differential(t)
# 3. Parameters
pars_RC = @parameters begin
C=8, [description="Capacitance, C"]
R=0.1, [description="Resistance, Ω"]
end
# 4. Dependent variables
vars_RC = @variables begin
i(t), [description="Total current, A"]
i_c(t), [description="Capacitor current, A"]
i_r(t), [description="Resistor current, A"]
(u(t)=0), [description="Voltage, V"]
end
# RC circuit
eqs_a = [Dt(u) ~ i_c/C,
i_r ~ u/R,
i ~ i_r + i_c]
@named mod_RC_a = ODESystem(eqs_a,t,vars_RC,pars_RC)
# Specifying inputs
@register_symbolic i_u(t)
eqs_i = [i ~ i_u(t)]
@named i_input = ODESystem(eqs_i,t)
# Causal system
mod_RC = extend(i_input, mod_RC_a) |> structural_simplify
i_u(t) = 10
tspan=(0,5)
prob = ODEProblem(mod_RC, [],tspan)
sol = solve(prob)
plot(sol, lw=2.5)
```
This produces the following plot:
Next, using `MonteCarloMeasurements`:
```julia
prob = remake(prob; u0=[u => -0.2..0.2], p=[C=> 8±1])
sol = solve(prob)
plot(sol)
```
produces:
[`Warning: Using arrays or dicts to store parameters of different types can hurt performance. Consider using tuples instead.` – I neglected this advice…]
**Question 1** : In order for this to work, it is required that _at least_ one _parameter_ **and** _at least_ one _state/initial value_ simultaneously are specified with uncertainty, right? […to make sure that the quantities are properly promoted to the new type…]
According to the documentations, the above syntax for specifying uncertainty produces 2000 _particles_.
The following _alternative_ description also works:
```julia
Np = 10
prob = remake(prob; u0=[u => 0.2*Particles(Np)], p=[C=> 8 + Particles(Np)])
sol = solve(prob)
plot(sol)
```
with result:
**Question 2:** : This only works if I specify the same number of particles in every uncertainty, right? [`Np` above is used in both cases of calling `Particles`…]
**Question 3:** Clearly, the blue region in the two plots above are not the “shadows” of plotting multiple parameters – in the last plot, one may discern a “small” number of particle trajectories, probably 10.
- Does the blue region represent a ribbon plot of the computed mean +/- the computed standard deviation of the particles?
- Does the thicker blue line in the center of the ribbon plot represent the mean?
I observe that I can change the color of the filled area (“ribbon”?) with setting fill color (`fc=...`), the line color by specifying `lc`, and the line thickness by specifying `lw`.
**Question 4:** Regarding solutions…
- Is it possible to specify the line width, line color, line alpha value, and line style individually for the mean and the particles?
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [November 14, 2023, 1:40pm UTC](https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205/2 "2023-11-14T13:40:54Z")
> [@BLI](#):
>
> **Question 1** : In order for this to work, it is required that _at least_ one _parameter_ **and** _at least_ one _state/initial value_ simultaneously are specified with uncertainty, right?
In DiffEq it’s enough that the initial condition has the correct type, I’m not sure if MTK requires something different here.
> [@BLI](#):
>
> **Question 2:** : This only works if I specify the same number of particles in every uncertainty, right? [`Np` above is used in both cases of calling `Particles`…]
Correct
> [@BLI](#):
>
> Does the blue region represent a ribbon plot of the computed mean +/- the computed standard deviation of the particles?
almost, they are quantiles, see [plotting](https://baggepinnen.github.io/MonteCarloMeasurements.jl/latest/#Plotting-1) for more info.
> [@BLI](#):
>
> Does the thicker blue line in the center of the ribbon plot represent the mean?
yes
> [@BLI](#):
>
> Is it possible to specify the line width, line color, line alpha value, and line style individually for the mean and the particles?
No. There are also two plot recipes that are interacting here, the recipe for solutions and the one for particles. I know that this interaction is not always working perfectly. If you want full control, you can issue two plot commands, one without sample trajectories and one with.
* * *
If you want to use a custom number of particles by default, you can define your own operator, e.g.,
```julia
julia> import MonteCarloMeasurements # Don't use using to avoid conflicting import of ..
julia> (..)(a, b) = MonteCarloMeasurements.Particles(10, Uniform(a,b))
.. (generic function with 1 method)
julia> 1 .. 2
1.5 ± 0.303 Particles{Float64, 10}
```
---
### Author:  [@BLI](https://discourse.julialang.org/u/BLI)
#### Post date: [November 14, 2023, 2:20pm UTC](https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205/3 "2023-11-14T14:20:53Z")
> [@baggepinnen](#):
>
> lmost, they are quantiles, see [plotting](https://baggepinnen.github.io/MonteCarloMeasurements.jl/latest/#Plotting-1) for more info
Ah. OK – so presumably defaults to 0.025/0.975 quantiles for ribbon?
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [November 14, 2023, 2:45pm UTC](https://discourse.julialang.org/t/montecarlomeasurements-mtk-again/106205/4 "2023-11-14T14:45:42Z")
yes