Unverified Commit 334d9140 authored by Miroslav Kratochvil's avatar Miroslav Kratochvil Committed by GitHub
Browse files

Merge pull request #328 from LCSB-BioCore/mo-summarize-results

馃毟馃拕 Added summarize function
parents 97b7f639 4718fbb5
Pipeline #43080 passed with stages
in 30 minutes and 29 seconds
"""
FluxSummary
A struct used to store summary information about the solution
of a constraint based analysis result.
"""
struct FluxSummary
biomass_fluxes::OrderedDict{String,Float64}
import_fluxes::OrderedDict{String,Float64}
export_fluxes::OrderedDict{String,Float64}
unbounded_fluxes::OrderedDict{String,Float64}
end
"""
flux_summary(flux_result::Dict{String, Float64};
exclude_exchanges = false,
exchange_prefixes = _constants.exchange_prefixes,
biomass_strings = _constants.biomass_strings,
exclude_biomass = false,
small_flux_bound = 1.0/_constants.default_reaction_bound^2,
large_flux_bound = _constants.default_reaction_bound,
keep_unbounded = false,
)::FluxSummary
Summarize a dictionary of fluxes into small, useful representation of the most
important information contained. Useful for pretty-printing and quickly
exploring the results. Internally this function uses
[`looks_like_biomass_reaction`](@ref) and
[`looks_like_exchange_reaction`](@ref). The corresponding keyword arguments
passed to these functions. Use this if your model has non-standard ids for
reactions. Fluxes smaller than `small_flux_bound` are not stored, while fluxes
larger than `large_flux_bound` are only stored if `keep_unbounded` is `true`.
# Example
```
julia> sol = flux_balance_analysis_dict(model, Tulip.Optimizer)
julia> fr = flux_summary(sol)
Biomass:
BIOMASS_Ecoli_core_w_GAM: 0.8739
Import:
EX_o2_e: -21.7995
EX_glc__D_e: -10.0
EX_nh4_e: -4.7653
EX_pi_e: -3.2149
Export:
EX_h_e: 17.5309
EX_co2_e: 22.8098
EX_h2o_e: 29.1758
```
"""
function flux_summary(
flux_result::Dict{String,Float64};
exclude_exchanges = false,
exchange_prefixes = _constants.exchange_prefixes,
biomass_strings = _constants.biomass_strings,
exclude_biomass = false,
small_flux_bound = 1.0 / _constants.default_reaction_bound^2,
large_flux_bound = _constants.default_reaction_bound,
keep_unbounded = false,
)
rxn_ids = collect(keys(flux_result))
ex_rxns = filter(
x -> looks_like_exchange_reaction(
x,
exclude_biomass = exclude_biomass,
biomass_strings = biomass_strings,
exchange_prefixes = exchange_prefixes,
),
rxn_ids,
)
bmasses = filter(
x -> looks_like_biomass_reaction(
x;
exclude_exchanges = exclude_exchanges,
exchange_prefixes = exchange_prefixes,
biomass_strings = biomass_strings,
),
rxn_ids,
)
ex_fluxes = [flux_result[k] for k in ex_rxns]
bmass_fluxes = [flux_result[k] for k in bmasses]
idx_srt_fluxes = sortperm(ex_fluxes)
import_fluxes = [
idx for
idx in idx_srt_fluxes if -large_flux_bound < ex_fluxes[idx] <= -small_flux_bound
]
export_fluxes = [
idx for
idx in idx_srt_fluxes if small_flux_bound < ex_fluxes[idx] <= large_flux_bound
]
if keep_unbounded
lower_unbounded =
[idx for idx in idx_srt_fluxes if ex_fluxes[idx] <= -large_flux_bound]
upper_unbounded =
[idx for idx in idx_srt_fluxes if ex_fluxes[idx] >= large_flux_bound]
return FluxSummary(
OrderedDict(k => v for (k, v) in zip(bmasses, bmass_fluxes)),
OrderedDict(
k => v for (k, v) in zip(ex_rxns[import_fluxes], ex_fluxes[import_fluxes])
),
OrderedDict(
k => v for (k, v) in zip(ex_rxns[export_fluxes], ex_fluxes[export_fluxes])
),
OrderedDict(
k => v for (k, v) in zip(
[ex_rxns[lower_unbounded]; ex_rxns[upper_unbounded]],
[ex_fluxes[lower_unbounded]; ex_fluxes[upper_unbounded]],
)
),
)
else
return FluxSummary(
OrderedDict(k => v for (k, v) in zip(bmasses, bmass_fluxes)),
OrderedDict(
k => v for (k, v) in zip(ex_rxns[import_fluxes], ex_fluxes[import_fluxes])
),
OrderedDict(
k => v for (k, v) in zip(ex_rxns[export_fluxes], ex_fluxes[export_fluxes])
),
OrderedDict{String,Float64}(),
)
end
end
"""
FluxVariabilitySummary
A struct used to store summary information about the
solution of a flux variability analysis result.
"""
struct FluxVariabilitySummary
biomass_fluxes::Dict{String,Vector{Union{Float64,Nothing}}}
exchange_fluxes::Dict{String,Vector{Union{Float64,Nothing}}}
end
"""
flux_variability_summary(flux_result::Tuple{Dict{String, Dict{String, Float64}}, Dict{String, Dict{String, Float64}}};
exclude_exchanges = false,
exchange_prefixes = _constants.exchange_prefixes,
biomass_strings = _constants.biomass_strings,
exclude_biomass = false,
)::FluxVariabilitySummary
Summarize a dictionary of flux dictionaries obtained eg. from
flux_variability_analysis_dict. The simplified summary representation is useful
for pretty-printing and easily showing the most important results. Internally
this function uses [`looks_like_biomass_reaction`](@ref) and
[`looks_like_exchange_reaction`](@ref). The corresponding keyword arguments
passed to these functions. Use this if your model has non-standard ids for
reactions.
# Example
```
julia> sol = flux_variability_analysis_dict(model, Gurobi.Optimizer; bounds = objective_bounds(0.99))
julia> flux_res = flux_variability_summary(sol)
Biomass Lower bound Upper bound
BIOMASS_Ecoli_core_w_GAM: 0.8652 0.8652
Exchange
EX_h2o_e: 28.34 28.34
EX_co2_e: 22.0377 22.0377
EX_o2_e: -22.1815 -22.1815
EX_h_e: 17.3556 17.3556
EX_glc__D_e: -10.0 -10.0
EX_nh4_e: -4.8448 -4.8448
EX_pi_e: -3.2149 -3.2149
EX_for_e: 0.0 0.0
... ... ...
```
"""
function flux_variability_summary(
flux_result::Tuple{Dict{String,Dict{String,Float64}},Dict{String,Dict{String,Float64}}};
exclude_exchanges = false,
exchange_prefixes = _constants.exchange_prefixes,
biomass_strings = _constants.biomass_strings,
exclude_biomass = false,
)
rxn_ids = keys(flux_result[1])
ex_rxns = filter(
x -> looks_like_exchange_reaction(
x,
exclude_biomass = exclude_biomass,
biomass_strings = biomass_strings,
exchange_prefixes = exchange_prefixes,
),
rxn_ids,
)
bmasses = filter(
x -> looks_like_biomass_reaction(
x;
exclude_exchanges = exclude_exchanges,
exchange_prefixes = exchange_prefixes,
biomass_strings = biomass_strings,
),
rxn_ids,
)
biomass_fluxes = Dict{String,Vector{Union{Float64,Nothing}}}()
for rxn_id in bmasses
lb = isnothing(flux_result[1][rxn_id]) ? nothing : flux_result[1][rxn_id][rxn_id]
ub = isnothing(flux_result[2][rxn_id]) ? nothing : flux_result[2][rxn_id][rxn_id]
biomass_fluxes[rxn_id] = [lb, ub]
end
ex_rxn_fluxes = Dict{String,Vector{Union{Float64,Nothing}}}()
for rxn_id in ex_rxns
lb = isnothing(flux_result[1][rxn_id]) ? nothing : flux_result[1][rxn_id][rxn_id]
ub = isnothing(flux_result[2][rxn_id]) ? nothing : flux_result[2][rxn_id][rxn_id]
ex_rxn_fluxes[rxn_id] = [lb, ub]
end
return FluxVariabilitySummary(biomass_fluxes, ex_rxn_fluxes)
end
function _pad_spaces(slen::Int, maxlen::Int)
" "^max(0, maxlen - slen + 1)
end
_pad_spaces(str::String, maxlen::Int) = _pad_spaces(length(str), maxlen)
function Base.show(io::IO, ::MIME"text/plain", flux_res::FluxSummary)
longest_biomass_len = maximum([length(k) for k in keys(flux_res.biomass_fluxes)])
longest_import_len = maximum([length(k) for k in keys(flux_res.import_fluxes)])
longest_export_len = maximum([length(k) for k in keys(flux_res.export_fluxes)])
if !isempty(flux_res.unbounded_fluxes)
longest_unbounded_len =
maximum([length(k) for k in keys(flux_res.unbounded_fluxes)])
word_pad = max(
longest_biomass_len,
longest_export_len,
longest_import_len,
longest_unbounded_len,
)
else
word_pad = max(longest_biomass_len, longest_export_len, longest_import_len)
end
println(io, "Biomass")
for (k, v) in flux_res.biomass_fluxes
println(io, " ", k, ":", _pad_spaces(k, word_pad), round(v, digits = 4))
end
println(io, "Import")
for (k, v) in flux_res.import_fluxes
println(io, " ", k, ":", _pad_spaces(k, word_pad), round(v, digits = 4))
end
println(io, "Export")
for (k, v) in flux_res.export_fluxes
println(io, " ", k, ":", _pad_spaces(k, word_pad), round(v, digits = 4))
end
if !isempty(flux_res.unbounded_fluxes)
println(io, "Unbounded")
for (k, v) in flux_res.unbounded_fluxes
println(io, " ", k, ":", _pad_spaces(k, word_pad), round(v, digits = 4))
end
end
return nothing
end
function Base.show(io::IO, ::MIME"text/plain", flux_res::FluxVariabilitySummary)
longest_biomass_len = maximum([length(k) for k in keys(flux_res.biomass_fluxes)])
longest_exchange_len = maximum([length(k) for k in keys(flux_res.exchange_fluxes)])
word_pad_len = max(longest_biomass_len, longest_exchange_len)
longest_biomass_flux_len = maximum([
length(string(round(v[1], digits = 4))) for v in values(flux_res.biomass_fluxes)
])
longest_exchange_flux_len = max(
maximum([
length(string(round(v[1], digits = 4))) for
v in values(flux_res.exchange_fluxes)
]),
length("Lower bound "),
)
number_pad_len = max(longest_biomass_flux_len, longest_exchange_flux_len)
println(
io,
"Biomass",
_pad_spaces(length("Biomass"), word_pad_len + 3),
"Lower bound ",
_pad_spaces("Lower bound ", number_pad_len),
"Upper bound",
)
for (k, v) in flux_res.biomass_fluxes
lb = isnothing(v[1]) ? " " : string(round(v[1], digits = 4))
ub = isnothing(v[2]) ? " " : string(round(v[1], digits = 4))
println(
io,
" ",
k,
":",
_pad_spaces(k, word_pad_len),
lb,
_pad_spaces(lb, number_pad_len),
ub,
)
end
println(io, "Exchange")
ex_ids = collect(keys(flux_res.exchange_fluxes))
vs = [
abs(flux_res.exchange_fluxes[k][1]) + abs(flux_res.exchange_fluxes[k][2]) for
k in ex_ids
]
idxs = sortperm(vs, rev = true)
for k in ex_ids[idxs]
v = flux_res.exchange_fluxes[k]
lb = isnothing(v[1]) ? " " : string(round(v[1], digits = 4))
ub = isnothing(v[2]) ? " " : string(round(v[1], digits = 4))
println(
io,
" ",
k,
":",
_pad_spaces(k, word_pad_len),
lb,
_pad_spaces(lb, number_pad_len),
ub,
)
end
end
@testset "Flux summary" begin
model = load_model(model_paths["e_coli_core.json"])
sol = flux_balance_analysis_dict(
model,
Tulip.Optimizer;
modifications = [change_optimizer_attribute("IPM_IterationsLimit", 200)],
)
fr = flux_summary(sol; keep_unbounded = true, large_flux_bound = 25)
@test isapprox(
fr.biomass_fluxes["BIOMASS_Ecoli_core_w_GAM"],
0.8739215022690006;
atol = TEST_TOLERANCE,
)
@test isapprox(fr.export_fluxes["EX_co2_e"], 22.80983339307183; atol = TEST_TOLERANCE)
@test isapprox(fr.import_fluxes["EX_o2_e"], -21.799492758430517; atol = TEST_TOLERANCE)
@test isapprox(
fr.unbounded_fluxes["EX_h2o_e"],
29.175827202663395;
atol = TEST_TOLERANCE,
)
end
@testset "Flux variability summary" begin
model = load_model(model_paths["e_coli_core.json"])
sol = flux_variability_analysis_dict(
model,
Tulip.Optimizer;
bounds = objective_bounds(0.90),
modifications = [change_optimizer_attribute("IPM_IterationsLimit", 2000)],
)
fr = flux_variability_summary(sol)
@test isapprox(
fr.biomass_fluxes["BIOMASS_Ecoli_core_w_GAM"][1],
0.7865293520891825;
atol = TEST_TOLERANCE,
)
@test isapprox(
fr.exchange_fluxes["EX_for_e"][2],
11.322324494491848;
atol = TEST_TOLERANCE,
)
end
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment