Unverified Commit 4fae8d8d authored by Miroslav Kratochvil's avatar Miroslav Kratochvil 🚴
Browse files

The Big Rename of LinearModel

parent 38ae8761
......@@ -128,7 +128,7 @@ Note, this function only runs serially.
Consider using a different model type for parallel implementations.
Returns two dictionaries (`fva_max` and `fva_min`) that map each reaction `id` to dictionaries of the resultant flux distributions when that `id` is optimized.
See also: [`LinearModel`](@ref)
See also: [`CoreModel`](@ref)
# Example
```
......
"""
Convert LinearModel to the JuMP model, place objectives and the equality
Convert CoreModel to the JuMP model, place objectives and the equality
constraint.
"""
function make_optimization_model(
......@@ -21,7 +21,7 @@ function make_optimization_model(
end
"""
Use JuMP to solve an instance of LinearModel
Use JuMP to solve an instance of CoreModel
"""
function optimize_model(
model::LM,
......
"""
Convert a dictionary read from a MAT file to LinearModel
Convert a dictionary read from a MAT file to CoreModel
"""
function _convert_m_dict_to_linear_model(model::Dict)
model_keys = ["S", "b", "c", "ub", "lb"]
......@@ -18,14 +18,14 @@ function _convert_m_dict_to_linear_model(model::Dict)
rxns = vec(string.(model["rxns"]))
mets = vec(string.(model["mets"]))
return LinearModel(S, b, c, lb, ub, rxns, mets)
return CoreModel(S, b, c, lb, ub, rxns, mets)
end
"""
Convert a LinearModel to exportable format
Convert a CoreModel to exportable format
SparseVectors are not written and read properly, SparseMatrix is okay
"""
function _convert_to_m_exportable_dict(model::LinearModel)
function _convert_to_m_exportable_dict(model::CoreModel)
xl, xu = bounds(model)
return Dict(
"S" => stoichiometry(model),
......
......@@ -140,7 +140,7 @@ function _read_model(file_location::String, ::Type{MFile}, ::Type{StandardModel}
return StandardModel(model_id, rxns, mets, genes)
end
function _read_model(file_location::String, ::Type{MFile}, ::Type{LinearModel})
function _read_model(file_location::String, ::Type{MFile}, ::Type{CoreModel})
matfile = matread(file_location)
model_name = collect(keys(matfile))[1] # assume only one model per m-file
# of not then need to make this more resilient, maybe keyword args.. ?
......
......@@ -56,7 +56,7 @@ NB: Does NOT export general inequality constraints (eg coupling)
See also: `MAT.jl`
"""
function _write_model(model::LinearModel, ::Type{MFile}, file_path::String)
function _write_model(model::CoreModel, ::Type{MFile}, file_path::String)
var_name = "model" # maybe make a field for this in the model?
matwrite(file_path, Dict(var_name => _convert_to_m_exportable_dict(model)))
end
......@@ -2,9 +2,9 @@ function _read_model(file_location::String, ::Type{SBMLFile}, ::Type{SBMLModel})
return SBMLModel(SBML.readSBML(file_location)) # direct import
end
function _read_model(file_location::String, ::Type{SBMLFile}, ::Type{LinearModel})
function _read_model(file_location::String, ::Type{SBMLFile}, ::Type{CoreModel})
sbmlmodel = _read_model(file_location, SBMLFile, SBMLModel) # use other import
return convert(LinearModel, sbmlmodel)
return convert(CoreModel, sbmlmodel)
end
function _read_model(file_location::String, ::Type{SBMLFile}, ::Type{StandardModel})
......
"""
Add constraints of the following form to a CoupledLinearModel and return a modified one.
Add constraints of the following form to a CoreCoupledModel and return a modified one.
The arguments are same as for in-place `add_coupling_constraints!`.
"""
function add_coupling_constraints(m::CoupledLinearModel, args...)
function add_coupling_constraints(m::CoreCoupledModel, args...)
new_lp = deepcopy(m)
add_coupling_constraints!(new_lp, args...)
return new_lp
end
"""
Add constraints to a plain `LinearModel` (converts it to the coupled type)
Add constraints to a plain `CoreModel` (converts it to the coupled type)
"""
add_coupling_constraints(m::LinearModel, args...) =
add_coupling_constraints(convert(CoupledLinearModel, m), args...)
add_coupling_constraints(m::CoreModel, args...) =
add_coupling_constraints(convert(CoreCoupledModel, m), args...)
"""
In-place add coupling constraints in form
......@@ -22,7 +22,7 @@ In-place add coupling constraints in form
```
"""
function add_coupling_constraints!(
m::CoupledLinearModel,
m::CoreCoupledModel,
c::V,
cl::AbstractFloat,
cu::AbstractFloat,
......@@ -37,7 +37,7 @@ end
function add_coupling_constraints!(
m::CoupledLinearModel,
m::CoreCoupledModel,
C::M,
cl::V,
cu::V,
......@@ -59,7 +59,7 @@ Remove coupling constraints from the linear model and return the modified model.
Arguments are the same as for in-place version `remove_coupling_constraints!`.
"""
function remove_coupling_constraints(m::CoupledLinearModel, args...)
function remove_coupling_constraints(m::CoreCoupledModel, args...)
new_model = deepcopy(m)
remove_coupling_constraints!(new_model, args...)
return new_model
......@@ -67,14 +67,14 @@ end
"""
Removes a set of coupling constraints from a CoupledLinearModel in-place.
Removes a set of coupling constraints from a CoreCoupledModel in-place.
"""
function remove_coupling_constraints!(m::CoupledLinearModel, constraint::Int)
function remove_coupling_constraints!(m::CoreCoupledModel, constraint::Int)
remove_coupling_constraints!(m, [constraint])
end
function remove_coupling_constraints!(m::CoupledLinearModel, constraints::Vector{Int})
function remove_coupling_constraints!(m::CoreCoupledModel, constraints::Vector{Int})
to_be_kept = filter(e -> e constraints, 1:n_coupling_constraints(m))
m.C = m.C[to_be_kept, :]
m.cl = m.cl[to_be_kept]
......@@ -86,7 +86,7 @@ end
Change the lower and/or upper bounds ('cl' and 'cu') for given coupling constraints
"""
function change_coupling_bounds!(
model::CoupledLinearModel,
model::CoreCoupledModel,
constraints::Vector{Int};
cl::V = Float64[],
cu::V = Float64[],
......
......@@ -2,7 +2,7 @@
Adds reactions to the model `m`
"""
function add_reactions(
m::LinearModel,
m::CoreModel,
s::V1,
b::V2,
c::AbstractFloat,
......@@ -23,7 +23,7 @@ end
function add_reactions(
m::LinearModel,
m::CoreModel,
s::V1,
b::V2,
c::AbstractFloat,
......@@ -47,7 +47,7 @@ function add_reactions(
end
function add_reactions(
m::LinearModel,
m::CoreModel,
Sp::M,
b::V,
c::V,
......@@ -71,7 +71,7 @@ function add_reactions(
end
function add_reactions(m1::LinearModel, m2::LinearModel; check_consistency = false)
function add_reactions(m1::CoreModel, m2::CoreModel; check_consistency = false)
return add_reactions(
m1,
m2.S,
......@@ -87,7 +87,7 @@ end
function add_reactions(
m::LinearModel,
m::CoreModel,
Sp::M,
b::V,
c::V,
......@@ -147,8 +147,8 @@ function add_reactions(
newxl = vcat(m.xl, xl[new_reactions])
newxu = vcat(m.xu, xu[new_reactions])
new_rxns = vcat(m.rxns, rxns[new_reactions])
#new_lp = LinearModel(new_s, newb, new_c, m.cl, m.cu, newc, newxl, newxu, new_rxns, new_mets)
new_lp = LinearModel(new_s, newb, newc, newxl, newxu, new_rxns, new_mets)
#new_lp = CoreModel(new_s, newb, new_c, m.cl, m.cu, newc, newxl, newxu, new_rxns, new_mets)
new_lp = CoreModel(new_s, newb, newc, newxl, newxu, new_rxns, new_mets)
if check_consistency
return (new_lp, new_reactions, new_metabolites)
......@@ -161,7 +161,7 @@ end
Verifies the consistency of a given model
"""
function verify_consistency(
m::LinearModel,
m::CoreModel,
Sp::M,
b::V,
c::V,
......@@ -201,10 +201,10 @@ function verify_consistency(
end
"""
Removes a set of reactions from a LinearModel.
Removes a set of reactions from a CoreModel.
Also removes the metabolites not involved in any reaction.
"""
function remove_reactions(m::LinearModel, rxns::Vector{Int})
function remove_reactions(m::CoreModel, rxns::Vector{Int})
rxns_to_keep = filter(e -> e rxns, 1:n_reactions(m))
temp_s = m.S[:, rxns_to_keep]
......@@ -220,22 +220,22 @@ function remove_reactions(m::LinearModel, rxns::Vector{Int})
new_rxns = m.rxns[rxns_to_keep]
new_mets = m.mets[mets_to_keep]
new_model =
#LinearModel(new_s, newb, new_c, m.cl, m.cu, newc, newxl, newxu, new_rxns, new_mets)
LinearModel(new_s, newb, newc, newxl, newxu, new_rxns, new_mets)
#CoreModel(new_s, newb, new_c, m.cl, m.cu, newc, newxl, newxu, new_rxns, new_mets)
CoreModel(new_s, newb, newc, newxl, newxu, new_rxns, new_mets)
return new_model
end
function remove_reactions(m::LinearModel, rxn::Integer)
function remove_reactions(m::CoreModel, rxn::Integer)
return remove_reactions(m, [rxn])
end
function remove_reactions(m::LinearModel, rxn::String)
function remove_reactions(m::CoreModel, rxn::String)
return remove_reactions(m, [rxn])
end
function remove_reactions(m::LinearModel, rxns::Vector{String})
function remove_reactions(m::CoreModel, rxns::Vector{String})
rxn_indices = [findfirst(isequal(name), m.rxns) for name in intersect(rxns, m.rxns)]
if isempty(rxn_indices)
return m
......@@ -250,7 +250,7 @@ Returns indices of exchange reactions.
Exchange reactions are identified based on most commonly used prefixes.
"""
function find_exchange_reactions(
model::LinearModel;
model::CoreModel;
exclude_biomass = false,
biomass_str::String = "biomass",
exc_prefs = ["EX_"; "Exch_"; "Ex_"],
......@@ -273,7 +273,7 @@ In practice returns the metabolites consumed by the reactions given by `find_exc
and if called with the same arguments, the two outputs correspond.
"""
function find_exchange_metabolites(
model::LinearModel;
model::CoreModel;
exclude_biomass = false,
biomass_str::String = "biomass",
exc_prefs = ["EX_"; "Exch_"; "Ex_"],
......@@ -293,7 +293,7 @@ end
Change the lower and/or upper bounds ('xl' and 'xu') for given reactions
"""
function change_bounds!(
model::LinearModel,
model::CoreModel,
rxns::Vector{Int};
xl::V = Float64[],
xu::V = Float64[],
......@@ -316,7 +316,7 @@ end
function change_bounds!(
model::LinearModel,
model::CoreModel,
rxns::Vector{String};
xl::V = Float64[],
xu::V = Float64[],
......
"""
struct CoupledLinearModel <: MetabolicModel
struct CoreCoupledModel <: MetabolicModel
The linear model with additional coupling constraints in the form
```
cₗ ≤ C x ≤ cᵤ
```
"""
mutable struct CoupledLinearModel <: MetabolicModel
lm::LinearModel
mutable struct CoreCoupledModel <: MetabolicModel
lm::CoreModel
C::SparseMat
cl::SparseVec
cu::SparseVec
function CoupledLinearModel(
function CoreCoupledModel(
lm::MetabolicModel,
C::M,
cl::V,
......@@ -25,50 +25,50 @@ mutable struct CoupledLinearModel <: MetabolicModel
size(C) == (length(cu), n_reactions(lm)) ||
throw(DimensionMismatch("wrong dimensions of `C`"))
new(convert(LinearModel, lm), sparse(C), sparse(cl), sparse(cu))
new(convert(CoreModel, lm), sparse(C), sparse(cl), sparse(cu))
end
end
reactions(a::CoupledLinearModel) = reactions(a.lm)
metabolites(a::CoupledLinearModel) = metabolites(a.lm)
stoichiometry(a::CoupledLinearModel) = stoichiometry(a.lm)
bounds(a::CoupledLinearModel) = bounds(a.lm)
balance(a::CoupledLinearModel) = balance(a.lm)
objective(a::CoupledLinearModel) = objective(a.lm)
reactions(a::CoreCoupledModel) = reactions(a.lm)
metabolites(a::CoreCoupledModel) = metabolites(a.lm)
stoichiometry(a::CoreCoupledModel) = stoichiometry(a.lm)
bounds(a::CoreCoupledModel) = bounds(a.lm)
balance(a::CoreCoupledModel) = balance(a.lm)
objective(a::CoreCoupledModel) = objective(a.lm)
"""
coupling(a::CoupledLinearModel)::SparseMat
coupling(a::CoreCoupledModel)::SparseMat
Coupling constraint matrix for a `CoupledLinearModel`.
Coupling constraint matrix for a `CoreCoupledModel`.
"""
function coupling(a::CoupledLinearModel)::SparseMat
function coupling(a::CoreCoupledModel)::SparseMat
a.C
end
"""
n_coupling_constraints(a::CoupledLinearModel)::Int
n_coupling_constraints(a::CoreCoupledModel)::Int
The number of coupling constraints in a `CoupledLinearModel`.
The number of coupling constraints in a `CoreCoupledModel`.
"""
function n_coupling_constraints(a::CoupledLinearModel)::Int
function n_coupling_constraints(a::CoreCoupledModel)::Int
return size(a.C, 1)
end
"""
coupling_bounds(a::CoupledLinearModel)::Tuple{SparseVec,SparseVec}
coupling_bounds(a::CoreCoupledModel)::Tuple{SparseVec,SparseVec}
Coupling bounds for a `CoupledLinearModel`.
Coupling bounds for a `CoreCoupledModel`.
"""
function coupling_bounds(a::CoupledLinearModel)::Tuple{SparseVec,SparseVec}
function coupling_bounds(a::CoreCoupledModel)::Tuple{SparseVec,SparseVec}
(a.cl, a.cu)
end
"""
Base.convert(::Type{CoupledLinearModel}, m::M) where {M <: MetabolicModel}
Base.convert(::Type{CoreCoupledModel}, m::M) where {M <: MetabolicModel}
Make a `CoupledLinearModel` out of any compatible model type.
Make a `CoreCoupledModel` out of any compatible model type.
"""
function Base.convert(::Type{CoupledLinearModel}, m::M) where {M<:MetabolicModel}
function Base.convert(::Type{CoreCoupledModel}, m::M) where {M<:MetabolicModel}
(cl, cu) = coupling_bounds(m)
CoupledLinearModel(convert(LinearModel, m), coupling(m), cl, cu)
CoreCoupledModel(convert(CoreModel, m), coupling(m), cl, cu)
end
"""
struct LinearModel <: MetabolicModel
struct CoreModel <: MetabolicModel
A concrete linear optimization problem of the form:
A "bare bones" core linear optimization problem of the form, with reaction and
metabolite names.
```
min c^T x
s.t. S x = b
xₗ ≤ x ≤ xᵤ
```
"""
mutable struct LinearModel <: MetabolicModel
mutable struct CoreModel <: MetabolicModel
S::SparseMat
b::SparseVec
c::SparseVec
......@@ -18,7 +19,7 @@ mutable struct LinearModel <: MetabolicModel
rxns::Vector{String}
mets::Vector{String}
function LinearModel(
function CoreModel(
S::M,
b::V,
c::V,
......@@ -39,86 +40,86 @@ mutable struct LinearModel <: MetabolicModel
end
"""
reactions(a::LinearModel)::Vector{String}
reactions(a::CoreModel)::Vector{String}
Get the reactions in a `LinearModel`.
Get the reactions in a `CoreModel`.
"""
function reactions(a::LinearModel)::Vector{String}
function reactions(a::CoreModel)::Vector{String}
a.rxns
end
"""
metabolites(a::LinearModel)::Vector{String}
metabolites(a::CoreModel)::Vector{String}
Metabolites in a `LinearModel`.
Metabolites in a `CoreModel`.
"""
function metabolites(a::LinearModel)::Vector{String}
function metabolites(a::CoreModel)::Vector{String}
a.mets
end
"""
stoichiometry(a::LinearModel)::SparseMat
stoichiometry(a::CoreModel)::SparseMat
`LinearModel` stoichiometry matrix.
`CoreModel` stoichiometry matrix.
"""
function stoichiometry(a::LinearModel)::SparseMat
function stoichiometry(a::CoreModel)::SparseMat
a.S
end
"""
bounds(a::LinearModel)::Tuple{SparseVec,SparseVec}
bounds(a::CoreModel)::Tuple{SparseVec,SparseVec}
`LinearModel` flux bounds.
`CoreModel` flux bounds.
"""
function bounds(a::LinearModel)::Tuple{SparseVec,SparseVec}
function bounds(a::CoreModel)::Tuple{SparseVec,SparseVec}
(a.xl, a.xu)
end
"""
balance(a::LinearModel)::SparseVec
balance(a::CoreModel)::SparseVec
`LinearModel` target flux balance.
`CoreModel` target flux balance.
"""
function balance(a::LinearModel)::SparseVec
function balance(a::CoreModel)::SparseVec
a.b
end
"""
objective(a::LinearModel)::SparseVec
objective(a::CoreModel)::SparseVec
`LinearModel` objective vector.
`CoreModel` objective vector.
"""
function objective(a::LinearModel)::SparseVec
function objective(a::CoreModel)::SparseVec
a.c
end
"""
coupling(a::LinearModel)::SparseMat
coupling(a::CoreModel)::SparseMat
Coupling constraint matrix for a `LinearModel`, actually empty.
Coupling constraint matrix for a `CoreModel`, actually empty.
"""
function coupling(a::LinearModel)::SparseMat
function coupling(a::CoreModel)::SparseMat
spzeros(0, n_reactions(a))
end
"""
coupling_bounds(a::LinearModel)::Tuple{SparseVec,SparseVec}
coupling_bounds(a::CoreModel)::Tuple{SparseVec,SparseVec}
Coupling bounds for a `LinearModel`, in fact always empty.
Coupling bounds for a `CoreModel`, in fact always empty.
"""
function coupling_bounds(a::LinearModel)::Tuple{SparseVec,SparseVec}
function coupling_bounds(a::CoreModel)::Tuple{SparseVec,SparseVec}
(spzeros(0), spzeros(0))
end
"""
Base.convert(::Type{LinearModel}, m::M) where {M <: MetabolicModel}
Base.convert(::Type{CoreModel}, m::M) where {M <: MetabolicModel}
Make a `LinearModel` out of any compatible model type.
Make a `CoreModel` out of any compatible model type.
"""
function Base.convert(::Type{LinearModel}, m::M) where {M<:MetabolicModel}
function Base.convert(::Type{CoreModel}, m::M) where {M<:MetabolicModel}
(xl, xu) = bounds(m)
LinearModel(
CoreModel(
stoichiometry(m),
balance(m),
objective(m),
......
Base.isequal(model1::LinearModel, model2::LinearModel) =
Base.isequal(model1::CoreModel, model2::CoreModel) =
isequal(model1.S, model2.S) &&
isequal(model1.b, model2.b) &&
isequal(model1.c, model2.c) &&
......@@ -7,14 +7,13 @@ Base.isequal(model1::LinearModel, model2::LinearModel) =
isequal(model1.rxns, model2.rxns) &&
isequal(model1.mets, model2.mets)
Base.copy(model::LinearModel) =
LinearModel(model.S, model.b, model.c, model.xl, model.xu, model.rxns, model.mets)
Base.copy(model::CoreModel) =
CoreModel(model.S, model.b, model.c, model.xl, model.xu, model.rxns, model.mets)
Base.isequal(model1::CoupledLinearModel, model2::CoupledLinearModel) =
Base.isequal(model1::CoreCoupledModel, model2::CoreCoupledModel) =
isequal(model1.lm, model2.lm) &&
isequal(model1.C, model2.C) &&
isequal(model1.cl, model2.cl) &&
isequal(model1.cu, model2.cu)
Base.copy(model::CoupledLinearModel) =
CoupledLinearModel(model.lm, model.C, model.cl, model.cu)
Base.copy(model::CoreCoupledModel) = CoreCoupledModel(model.lm, model.C, model.cl, model.cu)
@testset "Flux balance analysis with LinearModel" begin
@testset "Flux balance analysis with CoreModel" begin
cp = test_simpleLP()
lp = flux_balance_analysis(cp, GLPK.Optimizer)
@test termination_status(lp) === MOI.OPTIMAL
......@@ -24,7 +24,7 @@
model_path,
"d17be86293d4caafc32b829da4e2d0d76eb45e1bb837e0138327043a83e20c6e",
)
cp = read_model(model_path, LinearModel)
cp = read_model(model_path, CoreModel)
expected_optimum = 0.9219480950504393
lp = flux_balance_analysis(cp, GLPK.Optimizer)
......
......@@ -15,7 +15,7 @@
@test fluxes == Matrix{Float64}([2 2])
# a special testcase for slightly sub-optimal FVA (gamma<1)
cp = LinearModel(
cp = CoreModel(
[-1.0 -1.0 -1.0],
[0.0],
[1.0, 0.0, 0.0],
......
@testset "LinearModel type" begin
@testset "CoreModel type" begin
cp = test_LP()
@test cp isa LinearModel
@test cp isa CoreModel
end
@testset "LinearModel simple functions" begin
@testset "CoreModel simple functions" begin
cp = test_LP()
@test n_reactions(cp) == 3
@test n_metabolites(cp) == 4
......
test_LP() = LinearModel(
test_LP() = CoreModel(
zeros(4, 3),
zeros(4),
ones(3),
......@@ -8,7 +8,7 @@ test_LP() = LinearModel(
["m$x" for x = 1:4],
)
test_simpleLP() = LinearModel(
test_simpleLP() = CoreModel(
[
1.0 1.0
-1.0 1.0
......@@ -21,7 +21,7 @@ test_simpleLP() = LinearModel(
["m$x" for x = 1:2],
)