Skip to content
GitLab
Commit 8fa83b7f authored by Bentriou Mahmoud's avatar Bentriou Mahmoud
Browse files

Variables, parameters and transitions of models are now represented as 

models. It's a better semantic and improve performance as well as
readability of the code. All the tests passes (except the remark in the
last commit).
parent c3d6a1f7
Hide whitespace changes
Inline Side-by-side
algorithms/_utils_abc.jl
View file @ 8fa83b7f
function _init_abc_automaton!(mat_p_old::Matrix{Float64}, vec_dist::Vector{Float64},
pm::ParametricModel, str_var_aut::String;
pm::ParametricModel, sym_var_aut::VariableAutomaton;
l_obs::Union{Nothing,AbstractVector} = nothing,
func_dist::Union{Nothing,Function} = nothing)
vec_p = zeros(pm.df)
......@@ -9,7 +9,7 @@ function _init_abc_automaton!(mat_p_old::Matrix{Float64}, vec_dist::Vector{Float
mat_p_old[:,i] = vec_p
if l_obs == nothing
S = volatile_simulate(pm, vec_p)
vec_dist[i] = S[str_var_aut]
vec_dist[i] = S[sym_var_aut]
else
σ = simulate(pm, vec_p)
vec_dist[i] = func_dist(σ, l_obs)
......@@ -45,7 +45,7 @@ function _task_worker!(d_mat_p::DArray{Float64,2}, d_vec_dist::DArray{Float64,1}
pm::ParametricModel, epsilon::Float64,
wl_old::Vector{Float64}, mat_p_old::Matrix{Float64},
mat_cov::Matrix{Float64}, tree_mat_p::Union{Nothing,KDTree},
kernel_type::String, str_var_aut::String;
kernel_type::String, sym_var_aut::VariableAutomaton;
l_obs::Union{Nothing,AbstractVector} = nothing,
func_dist::Union{Nothing,Function} = nothing)
mat_p = localpart(d_mat_p)
......@@ -54,7 +54,7 @@ function _task_worker!(d_mat_p::DArray{Float64,2}, d_vec_dist::DArray{Float64,1}
l_nbr_sim = zeros(Int, length(vec_dist))
Threads.@threads for i = eachindex(vec_dist)
_update_param!(mat_p, vec_dist, l_nbr_sim, wl_current, i, pm, epsilon,
wl_old, mat_p_old, mat_cov, tree_mat_p, kernel_type, str_var_aut;
wl_old, mat_p_old, mat_cov, tree_mat_p, kernel_type, sym_var_aut;
l_obs = l_obs, func_dist = func_dist)
end
return sum(l_nbr_sim)
......@@ -86,7 +86,7 @@ function _update_param!(mat_p::Matrix{Float64}, vec_dist::Vector{Float64},
pm::ParametricModel, epsilon::Float64,
wl_old::Vector{Float64}, mat_p_old::Matrix{Float64},
mat_cov::Matrix{Float64}, tree_mat_p::Union{Nothing,KDTree},
kernel_type::String, str_var_aut::String;
kernel_type::String, sym_var_aut::VariableAutomaton;
l_obs::Union{Nothing,AbstractVector} = nothing,
func_dist::Union{Nothing,Function} = nothing)
d_weights = Categorical(wl_old)
......@@ -102,7 +102,7 @@ function _update_param!(mat_p::Matrix{Float64}, vec_dist::Vector{Float64},
end
if l_obs == nothing
S = volatile_simulate(pm, vec_p_prime)
dist_sim = S[str_var_aut]
dist_sim = S[sym_var_aut]
else
σ = simulate(pm, vec_p)
dist_sim = func_dist(σ, l_obs)
......
algorithms/abc_smc.jl
View file @ 8fa83b7f
......@@ -28,7 +28,7 @@ end
function automaton_abc(pm::ParametricModel; nbr_particles::Int = 100, alpha::Float64 = 0.75, kernel_type = "mvnormal",
NT::Float64 = nbr_particles/2, duration_time::Float64 = Inf,
bound_sim::Int = typemax(Int), str_var_aut::String = "d", verbose::Int = 0)
bound_sim::Int = typemax(Int), sym_var_aut::VariableAutomaton = :d, verbose::Int = 0)
@assert typeof(pm.m) <: SynchronizedModel
@assert 0 < nbr_particles
@assert 0.0 < alpha < 1.0
......@@ -41,14 +41,14 @@ function automaton_abc(pm::ParametricModel; nbr_particles::Int = 100, alpha::Flo
write(file_cfg, "kernel type : $(kernel_type) \n")
close(file_cfg)
if nprocs() == 1
return _abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, str_var_aut)
return _abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, sym_var_aut)
end
return _distributed_abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, str_var_aut)
return _distributed_abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, sym_var_aut)
end
"""
`abc_smc(pm::ParametricModel, l_obs, func_dist; nbr_particles, alpha, kernel_type, NT
duration_tiùe, bound_sim, str_var_aut, verbose)`
duration_tiùe, bound_sim, sym_var_aut, verbose)`
Run the ABC-SMC algorithm with the pm parametric model.
`func_dist(l_sim, l_obs)` is the distance function between simulations and observation,
......@@ -66,7 +66,7 @@ If pm is defined on a ContinuousTimeModel, then `T1` should verify `T1 <: Trajec
function abc_smc(pm::ParametricModel, l_obs::AbstractVector, func_dist::Function;
nbr_particles::Int = 100, alpha::Float64 = 0.75, kernel_type = "mvnormal",
NT::Float64 = nbr_particles/2, duration_time::Float64 = Inf,
bound_sim::Int = typemax(Int), str_var_aut::String = "d", verbose::Int = 0)
bound_sim::Int = typemax(Int), sym_var_aut::VariableAutomaton = :d, verbose::Int = 0)
@assert 0 < nbr_particles
@assert 0.0 < alpha < 1.0
@assert kernel_type in ["mvnormal", "knn_mvnormal"]
......@@ -78,9 +78,9 @@ function abc_smc(pm::ParametricModel, l_obs::AbstractVector, func_dist::Function
write(file_cfg, "kernel type : $(kernel_type) \n")
close(file_cfg)
if nprocs() == 1
return _abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, str_var_aut; l_obs = l_obs, func_dist = func_dist)
return _abc_smc(pm, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, sym_var_aut; l_obs = l_obs, func_dist = func_dist)
end
return _distributed_abc_smc(pm, l_obs, dist, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, str_var_aut; l_obs = l_obs, func_dist = func_dist)
return _distributed_abc_smc(pm, l_obs, dist, nbr_particles, alpha, kernel_type, NT, duration_time, bound_sim, dir_res, sym_var_aut; l_obs = l_obs, func_dist = func_dist)
end
......@@ -89,7 +89,7 @@ end
function _abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Float64,
kernel_type::String, NT::Float64, duration_time::Float64,
bound_sim::Int, dir_res::String, str_var_aut::String;
bound_sim::Int, dir_res::String, sym_var_aut::VariableAutomaton;
l_obs::Union{Nothing,AbstractVector} = nothing, func_dist::Union{Nothing,Function} = nothing)
@info "ABC PMC with $(nworkers()) processus and $(Threads.nthreads()) threads"
begin_time = time()
......@@ -102,7 +102,7 @@ function _abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Float64,
vec_dist = zeros(nbr_particles)
wl_old = zeros(nbr_particles)
@info "Step 1 : Init"
_init_abc_automaton!(mat_p_old, vec_dist, pm, str_var_aut; l_obs = l_obs, func_dist = func_dist)
_init_abc_automaton!(mat_p_old, vec_dist, pm, sym_var_aut; l_obs = l_obs, func_dist = func_dist)
prior_pdf!(wl_old, pm, mat_p_old)
normalize!(wl_old, 1)
ess = effective_sample_size(wl_old)
......@@ -141,7 +141,7 @@ function _abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Float64,
end
Threads.@threads for i = eachindex(vec_dist)
_update_param!(mat_p, vec_dist, l_nbr_sim, wl_current, i, pm, epsilon,
wl_old, mat_p_old, mat_cov, tree_mat_p, kernel_type, str_var_aut;
wl_old, mat_p_old, mat_cov, tree_mat_p, kernel_type, sym_var_aut;
l_obs = l_obs, func_dist = func_dist)
end
normalize!(wl_current, 1)
......@@ -182,7 +182,7 @@ end
function _distributed_abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Float64,
kernel_type::String, NT::Float64, duration_time::Float64, bound_sim::Int,
dir_res::String, str_var_aut::String;
dir_res::String, sym_var_aut::VariableAutomaton;
l_obs::Union{Nothing,AbstractVector} = nothing, func_dist::Union{Nothing,Function} = nothing)
@info "Distributed ABC PMC with $(nworkers()) processus and $(Threads.nthreads()) threads"
begin_time = time()
......@@ -195,7 +195,7 @@ function _distributed_abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Fl
vec_dist = zeros(nbr_particles)
wl_old = zeros(nbr_particles)
@info "Step 1 : Init"
_init_abc_automaton!(mat_p_old, vec_dist, pm, str_var_aut; l_obs = l_obs, func_dist = func_dist)
_init_abc_automaton!(mat_p_old, vec_dist, pm, sym_var_aut; l_obs = l_obs, func_dist = func_dist)
prior_pdf!(wl_old, pm, mat_p_old)
normalize!(wl_old, 1)
ess = effective_sample_size(wl_old)
......@@ -240,7 +240,7 @@ function _distributed_abc_smc(pm::ParametricModel, nbr_particles::Int, alpha::Fl
@async l_nbr_sim[t_id_w] =
remotecall_fetch(() -> _task_worker!(d_mat_p, d_vec_dist, d_wl_current,
pm, epsilon, wl_old, mat_p_old, mat_cov, tree_mat_p,
kernel_type, str_var_aut;
kernel_type, sym_var_aut;
l_obs = l_obs, func_dist = func_dist), id_w)
end
wl_current = convert(Array, d_wl_current)
......
automata/automaton_F.jl
View file @ 8fa83b7f
function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, str_obs::String)
@assert str_obs in m.g
function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, sym_obs::VariableModel)
@assert sym_obs in m.g "$(sym_obs) is not observed."
# Locations
locations = [:l0, :l1, :l2, :l3]
......@@ -14,8 +14,8 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
locations_final = [:l2]
#S.n <=> S.values[A.map_var_automaton_idx[:n]]
#P <=> xn[map_var_model_idx[constants[str_O]] with str_O = "P". On stock str_O dans constants
# P = get_value(A, x, str_obs)
#P <=> xn[map_var_model_idx[constants[str_O]] with str_O = :P. On stock str_O dans constants
# P = get_value(A, x, sym_obs)
## Map of automaton variables
map_var_automaton_idx = Dict{VariableAutomaton,Int}(:n => 1, :d => 2, :isabs => 3)
......@@ -38,7 +38,7 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
cc_aut_F_l0l1_1(A::LHA, S::StateLHA) = true
us_aut_F_l0l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1;
S[:n] = get_value(A, x, str_obs);
S[:n] = get_value(A, x, sym_obs);
S[:d] = Inf;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_F_l0l1_1, us_aut_F_l0l1_1!)
......@@ -108,9 +108,9 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
cc_aut_F_l3l1_1(A::LHA, S::StateLHA) = true
us_aut_F_l3l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1;
S[:n] = get_value(A, x, str_obs);
S[:n] = get_value(A, x, sym_obs);
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_F_l3l1_1, us_aut_F_l3l1_1!)
edge1 = Edge([:ALL], cc_aut_F_l3l1_1, us_aut_F_l3l1_1!)
map_edges[:l3][:l1] = [edge1]
......
automata/automaton_G.jl
View file @ 8fa83b7f
function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, str_obs::String)
@assert str_obs in m.g
function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, sym_obs::VariableModel)
@assert sym_obs in m.g
# Locations
locations = [:l0, :l1, :l2, :l3, :l4]
......@@ -38,7 +38,7 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
us_aut_G_l0l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1;
S[:d] = 0;
S[:n] = get_value(A, x, str_obs);
S[:n] = get_value(A, x, sym_obs);
S[:in] = true;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_G_l0l1_1, us_aut_G_l0l1_1!)
......@@ -142,9 +142,9 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
cc_aut_G_l3l1_1(A::LHA, S::StateLHA) = true
us_aut_G_l3l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1;
S[:n] = get_value(A, x, str_obs);
S[:n] = get_value(A, x, sym_obs);
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_G_l3l1_1, us_aut_G_l3l1_1!)
edge1 = Edge([:ALL], cc_aut_G_l3l1_1, us_aut_G_l3l1_1!)
map_edges[:l3][:l1] = [edge1]
......@@ -172,10 +172,10 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
(S.loc = :l1;
S[:d] += S[:tprime] * min(abs(A.constants[:x1] - S[:n]), abs(A.constants[:x2] - S[:n]));
S[:tprime] = 0.0;
S[:n] = get_value(A, x, str_obs);
S[:n] = get_value(A, x, sym_obs);
S[:in] = true;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_G_l4l1_1, us_aut_G_l4l1_1!)
edge1 = Edge([:ALL], cc_aut_G_l4l1_1, us_aut_G_l4l1_1!)
map_edges[:l4][:l1] = [edge1]
......
automata/automaton_G_and_F.jl
View file @ 8fa83b7f
function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, str_obs_G::String,
x3::Float64, x4::Float64, t3::Float64, t4::Float64, str_obs_F::String)
function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, sym_obs_G::VariableModel,
x3::Float64, x4::Float64, t3::Float64, t4::Float64, sym_obs_F::VariableModel)
@assert str_obs_G in m.g
@assert str_obs_F in m.g
@assert sym_obs_G in m.g
@assert sym_obs_F in m.g
# Locations
locations = [:l0G, :l1G, :l2G, :l3G, :l4G,
:l1F, :l2F, :l3F]
......@@ -49,7 +49,7 @@ function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float
us_aut_G_l0Gl1G_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1G;
S[:d] = 0;
S[:n] = get_value(A, x, str_obs_G);
S[:n] = get_value(A, x, sym_obs_G);
S[:in] = true;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_G_l0Gl1G_1, us_aut_G_l0Gl1G_1!)
......@@ -154,9 +154,9 @@ function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float
cc_aut_G_l3Gl1G_1(A::LHA, S::StateLHA) = true
us_aut_G_l3Gl1G_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1G;
S[:n] = get_value(A, x, str_obs_G);
S[:n] = get_value(A, x, sym_obs_G);
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_G_l3Gl1G_1, us_aut_G_l3Gl1G_1!)
edge1 = Edge([:ALL], cc_aut_G_l3Gl1G_1, us_aut_G_l3Gl1G_1!)
map_edges[:l3G][:l1G] = [edge1]
tuple = (:l3G, :l2G)
......@@ -184,10 +184,10 @@ function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float
(S.loc = :l1G;
S[:d] += S[:tprime] * min(abs(A.constants[:x1] - S[:n]), abs(A.constants[:x2] - S[:n]));
S[:tprime] = 0.0;
S[:n] = get_value(A, x, str_obs_G);
S[:n] = get_value(A, x, sym_obs_G);
S[:in] = true;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_G_l4Gl1G_1, us_aut_G_l4Gl1G_1!)
edge1 = Edge([:ALL], cc_aut_G_l4Gl1G_1, us_aut_G_l4Gl1G_1!)
map_edges[:l4G][:l1G] = [edge1]
......@@ -207,7 +207,7 @@ function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float
cc_aut_F_l2Gl1F_1(A::LHA, S::StateLHA) = true
us_aut_F_l2Gl1F_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1F;
S[:n] = get_value(A, x, str_obs_F);
S[:n] = get_value(A, x, sym_obs_F);
S[:dprime] = Inf;
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_F_l2Gl1F_1, us_aut_F_l2Gl1F_1!)
......@@ -280,9 +280,9 @@ function create_automaton_G_and_F(m::ContinuousTimeModel, x1::Float64, x2::Float
cc_aut_F_l3Fl1F_1(A::LHA, S::StateLHA) = true
us_aut_F_l3Fl1F_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = :l1F;
S[:n] = get_value(A, x, str_obs_F);
S[:n] = get_value(A, x, sym_obs_F);
S[:isabs] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_F_l3Fl1F_1, us_aut_F_l3Fl1F_1!)
edge1 = Edge([:ALL], cc_aut_F_l3Fl1F_1, us_aut_F_l3Fl1F_1!)
map_edges[:l3F][:l1F] = [edge1]
......
bench/array_memory_order/access_trajectory.jl
View file @ 8fa83b7f
......@@ -7,7 +7,7 @@ include(get_module_path() * "/core/_tests_simulate.jl")
BenchmarkTools.DEFAULT_PARAMETERS.samples = 20000
if ARGS[1] == "SIR"
bound_time = 200.0
l_var = ["S", "I", "R"]
l_var = [:S, :I, :R]
load_model("_bench_perf_test/SIR_col_buffer")
SIR_col_buffer.time_bound = bound_time
......@@ -19,7 +19,7 @@ if ARGS[1] == "SIR"
model_col_buffer = SIR_col_buffer
model_row_buffer = SIR_row_buffer
elseif ARGS[1] == "ER"
l_var = ["E","S","ES","P"]
l_var = [:E,:S,:ES,:P]
bound_time = 20.0
load_model("_bench_perf_test/ER_col_buffer")
......@@ -42,7 +42,7 @@ function access_trajectory_col(m::Model)
n_sim = 100
for k = 1:n_sim
σ = _simulate_col_buffer(m)
t = _get_values_col(σ, "I")
t = _get_values_col(σ, :I)
res += t[end-1]
end
return res
......@@ -59,7 +59,7 @@ function access_trajectory_row(m::Model)
n_sim = 100
for k = 1:n_sim
σ = _simulate_row_buffer(m)
t = _get_values_row(σ, "I")
t = _get_values_row(σ, :I)
res += t[end-1]
end
return res
......
bench/array_memory_order/multiple_sim.jl
View file @ 8fa83b7f
......@@ -6,7 +6,7 @@ using MarkovProcesses
include(get_module_path() * "/core/_tests_simulate.jl")
if ARGS[1] == "SIR"
l_var = ["S", "I", "R"]
l_var = [:S, :I, :R]
bound_time = 200.0
load_model("_bench_perf_test/SIR_col")
......@@ -23,7 +23,7 @@ if ARGS[1] == "SIR"
model_col_buffer = SIR_col_buffer
model_row_buffer = SIR_row_buffer
elseif ARGS[1] == "ER"
l_var = ["E","S","ES","P"]
l_var = [:E,:S,:ES,:P]
bound_time = 20.0
load_model("_bench_perf_test/ER_col")
......
bench/array_memory_order/read_random_state_trajectory.jl
View file @ 8fa83b7f
......@@ -7,7 +7,7 @@ include(get_module_path() * "/core/_tests_simulate.jl")
BenchmarkTools.DEFAULT_PARAMETERS.samples = 20000
if ARGS[1] == "SIR"
bound_time = 200.0
l_var = ["S", "I", "R"]
l_var = [:S, :I, :R]
load_model("_bench_perf_test/SIR_col_buffer")
SIR_col_buffer.time_bound = bound_time
......@@ -19,7 +19,7 @@ if ARGS[1] == "SIR"
model_col_buffer = SIR_col_buffer
model_row_buffer = SIR_row_buffer
elseif ARGS[1] == "ER"
l_var = ["E","S","ES","P"]
l_var = [:E,:S,:ES,:P]
bound_time = 20.0
load_model("_bench_perf_test/ER_col_buffer")
......
bench/array_memory_order/read_trajectory.jl
View file @ 8fa83b7f
......@@ -7,7 +7,7 @@ include(get_module_path() * "/core/_tests_simulate.jl")
BenchmarkTools.DEFAULT_PARAMETERS.samples = 20000
if ARGS[1] == "SIR"
bound_time = 200.0
l_var = ["S", "I", "R"]
l_var = [:S, :I, :R]
load_model("_bench_perf_test/SIR_col_buffer")
SIR_col_buffer.time_bound = bound_time
......@@ -19,7 +19,7 @@ if ARGS[1] == "SIR"
model_col_buffer = SIR_col_buffer
model_row_buffer = SIR_row_buffer
elseif ARGS[1] == "ER"
l_var = ["E","S","ES","P"]
l_var = [:E,:S,:ES,:P]
bound_time = 20.0
load_model("_bench_perf_test/ER_col_buffer")
......@@ -44,7 +44,7 @@ function read_trajectory_col(m::Model)
n_read = 100000
for k = 1:n_read
for i = 1:n_states
res += _get_value_col(σ, "I", i)
res += _get_value_col(σ, :I, i)
end
end
return res
......@@ -63,7 +63,7 @@ function read_trajectory_row(m::Model)
n_read = 100000
for k = 1:n_read
for i = 1:n_states
res += _get_value_row(σ, "I", i)
res += _get_value_row(σ, :I, i)
end
end
return res
......
bench/array_memory_order/sim.jl
View file @ 8fa83b7f
......@@ -6,7 +6,7 @@ using MarkovProcesses
include(get_module_path() * "/core/_tests_simulate.jl")
if ARGS[1] == "SIR"
l_var = ["S", "I", "R"]
l_var = [:S, :I, :R]
bound_time = 200.0
load_model("_bench_perf_test/SIR_col")
......@@ -23,7 +23,7 @@ if ARGS[1] == "SIR"
model_col_buffer = SIR_col_buffer
model_row_buffer = SIR_row_buffer
elseif ARGS[1] == "ER"
l_var = ["E","S","ES","P"]
l_var = [:E,:S,:ES,:P]
bound_time = 20.0
nbr_sim = 10000
......
bench/pkg/catalyst.jl
View file @ 8fa83b7f
......@@ -5,8 +5,8 @@ using Catalyst
using DiffEqJump
load_model("ER")
set_param!(ER, "k1", 0.2)
set_param!(ER, "k2", 40.0)
set_param!(ER, :k1, 0.2)
set_param!(ER, :k2, 40.0)
ER.buffer_size = 100
ER.estim_min_states = 8000
......
bench/pygmalion/multiple_long_sim_er.jl
View file @ 8fa83b7f
......@@ -29,8 +29,8 @@ println("MarkovProcesses:")
using MarkovProcesses
MarkovProcesses.load_model("ER")
ER.time_bound = 20.0
set_param!(ER, "k1", 0.2)
set_param!(ER, "k2", 40.0)
set_param!(ER, :k1, 0.2)
set_param!(ER, :k2, 40.0)
@timev MarkovProcesses.simulate(ER)
println("Default buffer size=10")
......
bench/pygmalion/read_random_state_trajectory.jl
View file @ 8fa83b7f
......@@ -8,7 +8,7 @@ println("Pygmalion:")
str_m = "enzymatic_reaction"
str_d = "abc_er"
pygmalion.load_model(str_m)
l_var = ["E", "S", "ES", "P"]
l_var = [:E, :S, :ES, :P]
str_oml = "E,S,ES,P,R,time"
ll_om = split(str_oml, ",")
x0 = State(95.0, 95.0, 0.0, 0.0, 0.0, 0.0)
......@@ -18,11 +18,11 @@ tml = 1:400
g_all = create_observation_function([ObserverModel(str_oml, tml)])
so = pygmalion.simulate(f, g_all, x0, u, p_true; on = nothing, full_timeline = true)
function random_trajectory_value_pyg(so::SystemObservation)
n_states = get_number_of_observations(so, "P")
n_states = get_number_of_observations(so, :P)
return to_vec(so, "P", rand(1:n_states))
end
function random_trajectory_state_pyg(so::SystemObservation)
n_states = get_number_of_observations(so, "P")
n_states = get_number_of_observations(so, :P)
return to_vec(so, so.oml, rand(1:n_states))
end
# Bench
......@@ -37,7 +37,7 @@ ER.time_bound = 10.0
σ = MarkovProcesses.simulate(ER)
function random_trajectory_value(σ::AbstractTrajectory)
n_states = length_states(σ)
return σ["P"][rand(1:n_states)]
return σ[:P][rand(1:n_states)]
end
# Bench
@timev random_trajectory_value(σ)
......
bench/pygmalion/read_trajectory_er.jl
View file @ 8fa83b7f
......@@ -26,7 +26,7 @@ function read_trajectory_pyg_v1(so::SystemObservation)
return res
end
function read_trajectory_pyg_v2(so::SystemObservation)
idx_P = om_findfirst("P", so.oml)
idx_P = om_findfirst(:P, so.oml)
n_states = length(so.otll[idx_P])
res = 0.0
for i = 1:n_states
......@@ -51,13 +51,13 @@ function read_trajectory_v1(σ::AbstractTrajectory)
n_states = length_states(σ)
res = 0
for i = 1:n_states
res += σ["P",i]
res += σ[:P,i]
end
return res
end
function read_trajectory_v2(σ::AbstractTrajectory)
n_states = length_states(σ)
vals = σ["P"]
vals = σ[:P]
res = 0
for i = 1:n_states
res += vals[i]
......
core/_tests_simulate.jl
View file @ 8fa83b7f
<
......@@ -3,16 +3,16 @@ struct OldTrajectory <: AbstractTrajectory
m::ContinuousTimeModel
values::Matrix{Int}
times::Vector{Float64}
transitions::Vector{Union{String,Nothing}}
transitions::Vector{Union{Symbol,Nothing}}
end
# File for benchmarking simulation and memory access of the package.
# Trajectories
_get_values_col(σ::OldTrajectory, var::String) =
_get_values_col(σ::OldTrajectory, var::Symbol) =
@view σ.values[(σ.m)._map_obs_var_idx[var],:]
_get_values_row(σ::OldTrajectory, var::String) =
_get_values_row(σ::OldTrajectory, var::Symbol) =
@view σ.values[:,(σ.m)._map_obs_var_idx[var]]
_get_state_col(σ::OldTrajectory, idx::Int) =
......@@ -20,9 +20,9 @@ _get_state_col(σ::OldTrajectory, idx::Int) =
_get_state_row(σ::OldTrajectory, idx::Int) =
@view σ.values[idx,:]
_get_value_col(σ::OldTrajectory, var::String, idx::Int) =
_get_value_col(σ::OldTrajectory, var::Symbol, idx::Int) =
σ.values[(σ.m)._map_obs_var_idx[var],idx]
_get_value_row(σ::OldTrajectory, var::String, idx::Int) =
_get_value_row(σ::OldTrajectory, var::Symbol, idx::Int) =
σ.values[idx,(σ.m)._map_obs_var_idx[var]]
# Model
......@@ -31,7 +31,7 @@ function _simulate_col(m::ContinuousTimeModel)
# trajectory fields
full_values = zeros(m.dim_state, 0)
times = zeros(0)
transitions = Vector{Union{String,Nothing}}(undef,0)
transitions = Vector{Union{Symbol,Nothing}}(undef,0)
# values at time n
n = 0
xn = @view m.x0[:]
......@@ -65,7 +65,7 @@ function _simulate_row(m::ContinuousTimeModel)
# trajectory fields
full_values = zeros(m.dim_state, 0)
times = zeros(0)
transitions = Vector{Union{String,Nothing}}(undef,0)
transitions = Vector{Union{Symbol,Nothing}}(undef,0)
# values at time n
n = 0
xn = @view m.x0[:]
......@@ -100,7 +100,7 @@ function _simulate_col_buffer(m::ContinuousTimeModel; buffer_size::Int = 5)
# trajectory fields
full_values = zeros(m.dim_state, 0)
times = zeros(0)
transitions = Vector{Union{String,Nothing}}(undef,0)
transitions = Vector{Union{Symbol,Nothing}}(undef,0)
# values at time n
n = 0
xn = @view m.x0[:]
......@@ -108,7 +108,7 @@ function _simulate_col_buffer(m::ContinuousTimeModel; buffer_size::Int = 5)
# at time n+1
mat_x = zeros(Int, m.dim_state, buffer_size)
l_t = zeros(Float64, buffer_size)
l_tr = Vector{Union{String,Nothing}}(undef, buffer_size)
l_tr = Vector{Union{Symbol,Nothing}}(undef, buffer_size)
isabsorbing = m.isabsorbing(m.p,xn)::Bool
while !isabsorbing && (tn <= m.time_bound)
i = 0
......@@ -140,7 +140,7 @@ function _simulate_row_buffer(m::ContinuousTimeModel; buffer_size::Int = 5)
# trajectory fields
full_values = zeros(0, m.dim_state)
times = zeros(0)
transitions = Vector{Union{String,Nothing}}(undef,0)
transitions = Vector{Union{Symbol,Nothing}}(undef,0)
# values at time n
n = 0
xn = @view m.x0[:]
......