Merge pull request #433 from LCSB-BioCore/mk-remove-stuff

start the grand unification of remove_whatever

docs/src/notebooks/7_community_model.jl

@@ -26,15 +26,15 @@ using Tulip
 
 base_model = load_model(CoreModel, "e_coli_core.json") # base from from which the knockouts will be constructed
 
-cytbd_knockout_model = remove_reactions(base_model, "CYTBD") # knockout the CYTBD (cytochrome oxidase) reaction
+cytbd_knockout_model = remove_reaction(base_model, "CYTBD") # knockout the CYTBD (cytochrome oxidase) reaction
 sol = flux_balance_analysis_dict(cytbd_knockout_model, Tulip.Optimizer)
 sol["BIOMASS_Ecoli_core_w_GAM"] # Cytochrome oxidase knockout μ (growth rate)
 #
-atps4r_knockout_model = remove_reactions(base_model, "ATPS4r") # knockout the ATP synthase reaction
+atps4r_knockout_model = remove_reaction(base_model, "ATPS4r") # knockout the ATP synthase reaction
 sol = flux_balance_analysis_dict(atps4r_knockout_model, Tulip.Optimizer)
 sol["BIOMASS_Ecoli_core_w_GAM"] # ATP synthase knockout μ
 #
-eno_knockout_model = remove_reactions(base_model, "ENO") # knockout the enolase reaction
+eno_knockout_model = remove_reaction(base_model, "ENO") # knockout the enolase reaction
 sol = flux_balance_analysis_dict(eno_knockout_model, Tulip.Optimizer)
 sol["BIOMASS_Ecoli_core_w_GAM"] # Enolase knockout μ, cannot grow by itself
 

src/base/macros/remove_item.jl

+macro _remove_fn(objname, model_type, idx_type, args...)
+    body = last(args)
+    typeof(body) == Expr || throw(DomainError(body, "missing function body"))
+    plural = :plural in args
+    plural_s = plural ? "s" : ""
+    inplace = :inplace in args
+    fname = Symbol(:remove_, objname, plural_s, inplace ? "!" : "")
+    idx_var = Symbol(
+        objname,
+        idx_type == :Int ? "_idx" :
+        idx_type == :String ? "_id" :
+        throw(DomainError(idx_type, "unsupported index type for _remove_fn macro")),
+        plural_s,
+    )
+
+    if plural
+        idx_type = AbstractVector{eval(idx_type)}
+    end
+
+    docstring = """
+        $fname(model::$model_type, $idx_var::$idx_type)
+
+    Remove $objname$plural_s from the model of type `$model_type`
+    $(inplace ? "in-place" : "and return the modified model").
+    """
+
+    return :(@doc $docstring $fname(model::$model_type, $idx_var::$idx_type) = $body)
+end

src/reconstruction/CoreModel.jl

@@ -232,12 +232,10 @@ function add_reactions(
 
     new_s = hcat(ext_s, ext_sp)
     newb = vcat(m.b, b[new_metabolites])
-    #new_c = hcat(m.C, spzeros(size(m.C, 1), length(new_reactions)))
     newc = vcat(m.c, c[new_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 = 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
@@ -305,110 +303,6 @@ function verify_consistency(
     return (new_reactions, new_metabolites)
 end
 
-"""
-    remove_metabolites(model::CoreModel, mets)
-
-Removes a set of `metabolites` from the `model` of type `CoreModel` and returns
-a new `CoreModel` without those metabolites. Here, `metabolites` can be either a
-string, a vector of strings, an index or a vector of indices. Also removes any
-reactions that have no associated metabolites after the metabolites have been
-removed.
-
-# Example
-```
-model = load_model(CoreModel, "e_coli_core.json")
-
-m1 = remove_metabolites(model, ["glc__D_e", "for_c"])
-m2 = remove_metabolites(model, "glc__D_e")
-m3 = remove_metabolites(model, indexin(["glc__D_e", "for_c"], metabolites(model)))
-m4 = remove_metabolites(model, first(indexin(["glc__D_e"], metabolites(model))))
-```
-"""
-function remove_metabolites(model::CoreModel, mets)
-    mets_to_keep = filter(!in(mets), 1:n_metabolites(model))
-    temp_S = model.S[mets_to_keep, :]
-
-    (I, rxns_to_keep, val) = findnz(temp_S)
-    sort!(rxns_to_keep)
-    unique!(rxns_to_keep)
-    new_S = model.S[mets_to_keep, rxns_to_keep]
-    new_b = model.b[mets_to_keep]
-    new_c = model.c[rxns_to_keep]
-    new_lbs = model.xl[rxns_to_keep]
-    new_ubs = model.xu[rxns_to_keep]
-    new_rxns = model.rxns[rxns_to_keep]
-    new_mets = model.mets[mets_to_keep]
-
-    return CoreModel(new_S, new_b, new_c, new_lbs, new_ubs, new_rxns, new_mets)
-end
-
-function remove_metabolites(model::CoreModel, met::Int)
-    return remove_metabolites(model, [met])
-end
-
-function remove_metabolites(model::CoreModel, met::String)
-    return remove_metabolites(model, [met])
-end
-
-function remove_metabolites(model::CoreModel, mets::Vector{String})
-    met_indices = filter(!isnothing, indexin(mets, metabolites(model)))
-    if isempty(met_indices)
-        return model
-    else
-        return remove_metabolites(model, met_indices)
-    end
-end
-
-"""
-    remove_reactions(m::CoreModel, rxns::Vector{Int})
-
-Removes a set of reactions from a CoreModel.
-Also removes the metabolites not involved in any reaction.
-"""
-function remove_reactions(m::CoreModel, rxns::Vector{Int})
-    rxns_to_keep = filter(!in(rxns), 1:n_reactions(m))
-    temp_s = m.S[:, rxns_to_keep]
-
-    (mets_to_keep, J, val) = findnz(temp_s)
-    sort!(mets_to_keep)
-    unique!(mets_to_keep)
-    new_s = m.S[mets_to_keep, rxns_to_keep]
-    newb = m.b[mets_to_keep]
-    newc = m.c[rxns_to_keep]
-    newxl = m.xl[rxns_to_keep]
-    newxu = m.xu[rxns_to_keep]
-    new_rxns = m.rxns[rxns_to_keep]
-    new_mets = m.mets[mets_to_keep]
-    new_model = CoreModel(new_s, newb, newc, newxl, newxu, new_rxns, new_mets)
-    return new_model
-end
-
-"""
-    remove_reactions(m::CoreModel, rxn::Int)
-"""
-function remove_reactions(m::CoreModel, rxn::Int)
-    return remove_reactions(m, [rxn])
-end
-
-"""
-    remove_reactions(m::CoreModel, rxn::String)
-"""
-function remove_reactions(m::CoreModel, rxn::String)
-    return remove_reactions(m, [rxn])
-end
-
-"""
-    remove_reactions(m::CoreModel, 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
-    else
-        return remove_reactions(m, rxn_indices)
-    end
-end
-
 @_change_bounds_fn CoreModel Int inplace begin
     isnothing(lower) || (model.xl[rxn_idx] = lower)
     isnothing(upper) || (model.xu[rxn_idx] = upper)
@@ -453,8 +347,93 @@ end
 @_change_bounds_fn CoreModel String plural begin
     change_bounds(
         model,
-        Vector{Int}(indexin(rxn_ids, reactions(model))),
+        Int.(indexin(rxn_ids, reactions(model))),
         lower = lower,
         upper = upper,
     )
 end
+
+@_remove_fn reaction CoreModel Int inplace begin
+    remove_reactions!(model, [reaction_idx])
+end
+
+@_remove_fn reaction CoreModel Int inplace plural begin
+    mask = .!in.(1:n_reactions(model), Ref(reaction_idxs))
+    model.S = model.S[:, mask]
+    model.c = model.c[mask]
+    model.xl = model.xl[mask]
+    model.xu = model.xu[mask]
+    model.rxns = model.rxns[mask]
+    return nothing
+end
+
+@_remove_fn reaction CoreModel Int begin
+    remove_reactions(model, [reaction_idx])
+end
+
+@_remove_fn reaction CoreModel Int plural begin
+    n = copy(model)
+    remove_reactions!(n, reaction_idxs)
+    return n
+end
+
+@_remove_fn reaction CoreModel String inplace begin
+    remove_reactions!(model, [reaction_id])
+end
+
+@_remove_fn reaction CoreModel String inplace plural begin
+    remove_reactions!(model, Int.(indexin(reaction_ids, reactions(model))))
+end
+
+@_remove_fn reaction CoreModel String begin
+    remove_reactions(model, [reaction_id])
+end
+
+@_remove_fn reaction CoreModel String plural begin
+    remove_reactions(model, Int.(indexin(reaction_ids, reactions(model))))
+end
+
+@_remove_fn metabolite CoreModel Int inplace begin
+    remove_metabolites!(model, [metabolite_idx])
+end
+
+@_remove_fn metabolite CoreModel Int plural inplace begin
+    remove_reactions!(
+        model,
+        [
+            ridx for ridx in 1:n_reactions(model) if
+            any(in.(findnz(model.S[:, ridx])[2], Ref(metabolite_idxs)))
+        ],
+    )
+    mask = .!in.(1:n_metabolites(model), Ref(metabolite_idxs))
+    model.S = model.S[mask, :]
+    model.b = model.b[mask]
+    model.mets = model.mets[mask]
+    return nothing
+end
+
+@_remove_fn metabolite CoreModel Int begin
+    remove_metabolites(model, [metabolite_idx])
+end
+
+@_remove_fn metabolite CoreModel Int plural begin
+    n = deepcopy(model) #everything gets changed anyway
+    remove_metabolites!(n, metabolite_idxs)
+    return n
+end
+
+@_remove_fn metabolite CoreModel String inplace begin
+    remove_metabolites!(model, [metabolite_id])
+end
+
+@_remove_fn metabolite CoreModel String inplace plural begin
+    remove_metabolites!(model, Int.(indexin(metabolite_ids, metabolites(model))))
+end
+
+@_remove_fn metabolite CoreModel String begin
+    remove_metabolites(model, [metabolite_id])
+end
+
+@_remove_fn metabolite CoreModel String plural begin
+    remove_metabolites(model, Int.(indexin(metabolite_ids, metabolites(model))))
+end

src/reconstruction/CoreModelCoupled.jl

@@ -164,52 +164,6 @@ function add_reactions(
     )
 end
 
-"""
-    remove_reactions(m::CoreModelCoupled, rxns::Vector{Int})
-
-Remove reaction(s) from a `CoreModelCoupled`.
-
-Also removes any metabolites not involved in any reaction after the deletion.
-
-"""
-function remove_reactions(m::CoreModelCoupled, rxns::Vector{Int})
-    return CoreModelCoupled(
-        remove_reactions(m.lm, rxns),
-        m.C[:, filter(!in(rxns), 1:n_reactions(m))],
-        m.cl,
-        m.cu,
-    )
-end
-
-"""
-    remove_reactions(m::CoreModelCoupled, rxn::Int)
-
-"""
-function remove_reactions(m::CoreModelCoupled, rxn::Int)
-    return remove_reactions(m, [rxn])
-end
-
-"""
-    remove_reactions(m::CoreModelCoupled, rxn::String)
-
-"""
-function remove_reactions(m::CoreModelCoupled, rxn::String)
-    return remove_reactions(m, [rxn])
-end
-
-"""
-    remove_reactions(m::CoreModelCoupled, rxns::Vector{String})
-
-"""
-function remove_reactions(m::CoreModelCoupled, rxns::Vector{String})
-    rxn_indices = [findfirst(isequal(rid), m.lm.rxns) for rid in intersect(rxns, m.lm.rxns)]
-    if isempty(rxn_indices)
-        return m
-    else
-        return remove_reactions(m, rxn_indices)
-    end
-end
-
 """
 Add constraints of the following form to a CoreModelCoupled and return a modified one.
 
@@ -400,3 +354,78 @@ end
     n.lm = change_bounds(model.lm, rxn_ids, lower = lower, upper = upper)
     n
 end
+
+@_remove_fn reaction CoreModelCoupled Int inplace begin
+    remove_reactions!(model, [reaction_idx])
+end
+
+@_remove_fn reaction CoreModelCoupled Int inplace plural begin
+    orig_rxns = reactions(model.lm)
+    remove_reactions!(model.lm, reaction_idxs)
+    model.C = model.C[:, in.(orig_rxns, Ref(Set(reactions(model.lm))))]
+    return nothing
+end
+
+@_remove_fn reaction CoreModelCoupled Int begin
+    remove_reactions(model, [reaction_idx])
+end
+
+@_remove_fn reaction CoreModelCoupled Int plural begin
+    n = copy(model)
+    n.lm = remove_reactions(n.lm, reaction_idxs)
+    n.C = n.C[:, in.(reactions(model.lm), Ref(Set(reactions(n.lm))))]
+    return n
+end
+
+@_remove_fn reaction CoreModelCoupled String inplace begin
+    remove_reactions!(model, [reaction_id])
+end
+
+@_remove_fn reaction CoreModelCoupled String inplace plural begin
+    remove_reactions!(model, Int.(indexin(reaction_ids, reactions(model))))
+end
+
+@_remove_fn reaction CoreModelCoupled String begin
+    remove_reactions(model, [reaction_id])
+end
+
+@_remove_fn reaction CoreModelCoupled String plural begin
+    remove_reactions(model, Int.(indexin(reaction_ids, reactions(model))))
+end
+
+@_remove_fn metabolite CoreModelCoupled Int inplace begin
+    remove_metabolites!(model, [metabolite_idx])
+end
+
+@_remove_fn metabolite CoreModelCoupled Int plural inplace begin
+    orig_rxns = reactions(model.lm)
+    model.lm = remove_metabolites(model.lm, metabolite_idxs)
+    model.C = model.C[:, in.(orig_rxns, Ref(Set(reactions(model.lm))))]
+    return nothing
+end
+
+@_remove_fn metabolite CoreModelCoupled Int begin
+    remove_metabolites(model, [metabolite_idx])
+end
+
+@_remove_fn metabolite CoreModelCoupled Int plural begin
+    n = deepcopy(model) #almost everything gets changed anyway
+    remove_metabolites!(n, metabolite_idxs)
+    return n
+end
+
+@_remove_fn metabolite CoreModelCoupled String inplace begin
+    remove_metabolites!(model, [metabolite_id])
+end
+
+@_remove_fn metabolite CoreModelCoupled String inplace plural begin
+    remove_metabolites!(model, Int.(indexin(metabolite_ids, metabolites(model))))
+end
+
+@_remove_fn metabolite CoreModelCoupled String begin
+    remove_metabolites(model, [metabolite_id])
+end
+
+@_remove_fn metabolite CoreModelCoupled String plural begin
+    remove_metabolites(model, Int.(indexin(metabolite_ids, metabolites(model))))
+end

src/reconstruction/StandardModel.jl

@@ -106,62 +106,6 @@ macro add_reactions!(model::Symbol, ex::Expr)
     return all_reactions
 end
 
-"""
-    remove_reactions!(model::StandardModel, ids::Vector{String})
-
-Remove all reactions with `ids` from `model`. Note, may result in orphan metabolites.
-
-# Example
-```
-remove_reactions!(model, ["EX_glc__D_e", "fba"])
-```
-"""
-function remove_reactions!(model::StandardModel, ids::Vector{String})
-    pop!.(Ref(model.reactions), ids)
-end
-
-"""
-    remove_reaction!(model::StandardModel, id::String)
-
-Remove reaction with `id` from `model`. Note, may result in orphan metabolites.
-
-# Example
-```
-remove_reaction!(model, "EX_glc__D_e")
-```
-"""
-remove_reaction!(model::StandardModel, id::String) = remove_reactions!(model, [id])
-
-"""
-    remove_metabolites!(model::StandardModel, ids::Vector{String})
-
-Remove all metabolites with `ids` from `model`.
-Warning, this could leave the model inconsistent, e.g. a reaction might
-require the deleted metabolite, in which case analysis functions will error.
-
-# Example
-```
-remove_metabolites!(model, ["atp_c", "adp_c"])
-```
-"""
-function remove_metabolites!(model::StandardModel, ids::Vector{String})
-    pop!.(Ref(model.metabolites), ids)
-end
-
-"""
-    remove_metabolite!(model::StandardModel, id::String)
-
-Remove metabolite with `id` from `model`.
-Warning, this could leave the model inconsistent, e.g. a reaction might
-require the deleted metabolite, in which case analysis functions will error.
-
-# Example
-```
-remove_metabolite!(model, "atp_c")
-```
-"""
-remove_metabolite!(model::StandardModel, id::String) = remove_metabolites!(model, [id])
-
 """
     remove_genes!(
         model::StandardModel,
@@ -217,7 +161,7 @@ remove_gene!(model::StandardModel, gid::String; knockout_reactions::Bool = false
 @_change_bounds_fn StandardModel String inplace begin
     isnothing(lower) || (model.reactions[rxn_id].lb = lower)
     isnothing(upper) || (model.reactions[rxn_id].ub = upper)
-    nothing
+    return nothing
 end
 
 @_change_bounds_fn StandardModel String inplace plural begin
@@ -237,5 +181,52 @@ end
         n.reactions[i] = copy(n.reactions[i])
     end
     change_bounds!(n, rxn_ids, lower = lower, upper = upper)
-    n
+    return n
+end
+
+@_remove_fn reaction StandardModel String inplace begin
+    delete!(model.reactions, reaction_id)
+end
+
+@_remove_fn reaction StandardModel String inplace plural begin
+    remove_reaction!.(Ref(model), reaction_ids)
+end
+
+@_remove_fn reaction StandardModel String begin
+    remove_reactions(model, [reaction_id])
+end
+
+@_remove_fn reaction StandardModel String plural begin
+    n = copy(model)
+    n.reactions = copy(model.reactions)
+    remove_reactions!(n, reaction_ids)
+    return n
+end
+
+@_remove_fn metabolite StandardModel String inplace begin
+    remove_metabolites!(model, [metabolite_id])
+end
+
+@_remove_fn metabolite StandardModel String inplace plural begin
+    remove_reactions!(
+        model,
+        [
+            rid for (rid, rn) in model.reactions if
+            any(haskey.(Ref(rn.metabolites), metabolite_ids))
+        ],
+    )
+    delete!.(Ref(model.metabolites), metabolite_ids)
+    return nothing
+end
+
+@_remove_fn metabolite StandardModel String begin
+    remove_metabolites(model, [metabolite_id])
+end
+
+@_remove_fn metabolite StandardModel String plural begin
+    n = copy(model)
+    n.reactions = copy(model.reactions)
+    n.metabolites = copy(model.metabolites)
+    remove_metabolites!(n, metabolite_ids)
+    return n
 end

src/reconstruction/modifications/generic.jl

@@ -14,16 +14,23 @@ Specifies a model variant that has new bounds set. Forwards arguments to
 """
 with_changed_bounds(args...; kwargs...) = m -> change_bounds(m, args...; kwargs...)
 
+"""
+    with_removed_metabolite(args...; kwargs...)
+
+Specifies a model variant without a certain metabolite. Forwards arguments to
+[`remove_metabolite`](@ref). Intended to be used with [`screen`](@ref).
+"""
+with_removed_metabolite(args...; kwargs...) = m -> remove_metabolite(m, args...; kwargs...)
+
 """
     with_removed_metabolites(args...; kwargs...)
 
-Specifies a model variant without specified metabolites. Forwards arguments to
-[`remove_metabolites`](@ref). Intended to be used with [`screen`](@ref).
+Plural version of [`with_removed_metabolite`](@ref), calls
+[`remove_metabolites`](@ref) internally.
 """
 with_removed_metabolites(args...; kwargs...) =
     m -> remove_metabolites(m, args...; kwargs...)
 
-
 """
     with_added_reactions(args...; kwargs...)
 
@@ -32,10 +39,18 @@ Specifies a model variant with reactions added. Forwards the arguments to
 """
 with_added_reactions(args...; kwargs...) = m -> add_reactions(m, args...; kwargs...)
 
+"""
+    with_removed_reaction(args...; kwargs...)
+
+Specifies a model variant without a certain reaction. Forwards arguments to
+[`remove_reaction`](@ref). Intended to be used with [`screen`](@ref).
+"""
+with_removed_reaction(args...; kwargs...) = m -> remove_reaction(m, args...; kwargs...)
+
 """
     with_removed_reactions(args...; kwargs...)
 
-Specifies a model variant with specified reactions removed. Forwards arguments
-to [`remove_reactions`](@ref). Intended to be used with [`screen`](@ref).
+Plural version of [`with_removed_reaction`](@ref), calls
+[`remove_reactions`](@ref) internally.
 """
 with_removed_reactions(args...; kwargs...) = m -> remove_reactions(m, args...; kwargs...)

test/reconstruction/CoreModel.jl

@@ -43,7 +43,6 @@ end
 
 @testset "Verify consistency" begin
     cp = test_LP()
-    @test size(cp.S) == (4, 3)
     (new_reactions, new_mets) = verify_consistency(
         cp,
         reshape(cp.S[:, end], :, 1),
@@ -77,7 +76,6 @@ end
 
 @testset "Add reactions (checking existence and consistency)" begin
     cp = test_LP()
-    @test size(cp.S) == (4, 3)
     (new_cp, new_reactions, new_mets) = add_reactions(
         cp,
         cp.S[:, end],
@@ -106,7 +104,6 @@ end
 
 @testset "Add reactions" begin
     cp = test_LP()
-    @test size(cp.S) == (4, 3)
     cp = add_reactions(cp, 2.0 * ones(4), 3 .* ones(4), 2.0, -1.0, 1.0)
     @test size(cp.S) == (8, 4)
     cp = add_reactions(cp, 2.0 * ones(4, 1), 3 .* ones(4), 2 .* ones(1), -ones(1), ones(1))
@@ -184,20 +181,16 @@ end
 
 @testset "Remove reactions" begin
     cp = test_LP()
-    @test size(cp.S) == (4, 3)
-    cp = remove_reactions(cp, 2)
-    @test size(cp.S) == (0, 2)
+    cp = remove_reaction(cp, 2)
+    @test size(cp.S) == (4, 2