fba.jl 3.45 KB
Newer Older
1
"""
2
    flux_balance_analysis(model::M, optimizer) where {M<:MetabolicModel}
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
3
4

Flux balance analysis solves the following problem for the input `model`:
Sylvain Arreckx's avatar
Sylvain Arreckx committed
5
6
7
8
9
```
max cᵀx
s.t. S x = b
     xₗ ≤ x ≤ xᵤ
```
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
10

11
Returns a solved model from [`optimize_model`](@ref).
Sylvain Arreckx's avatar
Sylvain Arreckx committed
12
"""
13
14
flux_balance_analysis(model::M, optimizer) where {M<:MetabolicModel} =
    optimize_model(model, optimizer; sense = MOI.MAX_SENSE)
15

16
"""
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
17
    flux_balance_analysis_vec(args...)::Union{Vector{Float64},Nothing}
18
19
20
21
22

A variant of FBA that returns a vector of fluxes in the same order as reactions
of the model, if the solution is found.
Arguments are passed to [`flux_balance_analysis`](@ref).
"""
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
23
function flux_balance_analysis_vec(args...)::Union{Vector{Float64},Nothing}
St. Elmo's avatar
St. Elmo committed
24
    optmodel = flux_balance_analysis(args...)
St. Elmo's avatar
St. Elmo committed
25
    
26
    termination_status(optmodel) in [MOI.OPTIMAL, MOI.LOCALLY_SOLVED] || return nothing
St. Elmo's avatar
St. Elmo committed
27
    value.(optmodel[:x])
28
29
30
end

"""
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
31
    flux_balance_analysis_dict(model::M, args...)::Union{Dict{String, Float64},Nothing} where {M <: MetabolicModel}
32
33
34
35
36
37
38

A variant of FBA that returns a dictionary assigning fluxes to reactions, if
the solution is found. Arguments are passed to [`flux_balance_analysis`](@ref).
"""
function flux_balance_analysis_dict(
    model::M,
    args...,
Miroslav Kratochvil's avatar
Miroslav Kratochvil committed
39
)::Union{Dict{String,Float64},Nothing} where {M<:MetabolicModel}
40
41
42
43
44
    v = flux_balance_analysis_vec(model, args...)
    isnothing(v) && return nothing
    Dict(zip(reactions(model), v))
end

45
"""
St. Elmo's avatar
St. Elmo committed
46
    flux_balance_analysis(model::CobraModel, optimizer; objective_func::Union{Reaction, Array{Reaction, 1}}=Reaction[], weights=Float64[], solver_attributes=Dict{Any, Any}(), constraints=Dict{String, Tuple{Float64,Float64}}())
47

St. Elmo's avatar
St. Elmo committed
48
Run flux balance analysis (FBA) on the `model` optionally specifying `objective_rxn(s)` and their `weights` (empty `weights` mean equal weighting per reaction).
49
Optionally also specify any additional flux constraints with `constraints`, a dictionary mapping reaction `id`s to tuples of (lb, ub) flux constraints.
50
Note, the `optimizer` must be set to perform the analysis, any JuMP solver will work.
51
52
53
The `solver_attributes` can also be specified in the form of a dictionary where each (key, value) pair will be passed to `set_optimizer_attribute(cbmodel, key, value)`.
This function builds the optimization problem from the model, and hence uses the constraints implied by the model object.
Returns a dictionary of reaction `id`s mapped to fluxes if solved successfully, otherwise an empty dictionary.
54

55
56
57
58
59
60
61
62
# Example
```
optimizer = Gurobi.Optimizer
atts = Dict("OutputFlag" => 0)
model = CobraTools.read_model("iJO1366.json")
biomass = findfirst(model.reactions, "BIOMASS_Ec_iJO1366_WT_53p95M")
sol = fba(model, biomass, optimizer; solver_attributes=atts)
```
63
"""
St. Elmo's avatar
St. Elmo committed
64
function flux_balance_analysis(
65
66
    model::CobraModel,
    optimizer;
St. Elmo's avatar
St. Elmo committed
67
    modifications::Union{Vector, Function} = [(model, opt_model)->nothing]
68
)
St. Elmo's avatar
St. Elmo committed
69
    # get core optimization problem
St. Elmo's avatar
St. Elmo committed
70
    cbm = make_optimization_model(model, optimizer)
St. Elmo's avatar
St. Elmo committed
71
    
St. Elmo's avatar
St. Elmo committed
72
73
74
75
76
77
78
    # apply callbacks - user can also just put in a function
    if typeof(modifications) <: Vector
        for mod in modifications
            mod(model, cbm)
        end    
    else
        modifications(model, cbm)
St. Elmo's avatar
St. Elmo committed
79
    end
St. Elmo's avatar
St. Elmo committed
80

81
    optimize!(cbm)
St. Elmo's avatar
St. Elmo committed
82
83
84
    
    return cbm
end
85

St. Elmo's avatar
St. Elmo committed
86
function flux_balance_analysis_vec()
87

St. Elmo's avatar
St. Elmo committed
88
89
90
91
92
93
94
95
96
97
98
99
end

status = (
    termination_status(cbm) == MOI.OPTIMAL ||
    termination_status(cbm) == MOI.LOCALLY_SOLVED
)

if status
    return map_fluxes(cbm[:v], model)
else
    @warn "Optimization issues occurred."
    return Dict{String,Float64}()
100
end