Merge pull request #76 from LCSB-BioCore/sew-small-fixes

Small cosmetic fixes for Memote.jl

Project.toml

@@ -1,7 +1,7 @@
 name = "FBCModelTests"
 uuid = "1060cd45-1c33-4a14-be81-8fa38fdd82bf"
 authors = ["The developers of FBCModelTests.jl"]
-version = "0.2.5"
+version = "0.3.0"
 
 [deps]
 COBREXA = "babc4406-5200-4a30-9033-bf5ae714c842"

src/memote/checks/Annotation.jl

@@ -19,31 +19,31 @@ is_annotated_gene(model, gid) = !isempty(gene_annotations(model, gid))
 """
 $(TYPEDSIGNATURES)
 
-Helper function to find all the annotations that are missing for a component in a model.
+Helper function to find all the annotations for a component in a model.
 """
-function _find_unannotated_components(
+function _find_annotated_components(
     model::MetabolicModel,
     id_accessor::String,
     annotation_accessor,
     annotation_kws,
 )
-    missing_annos = String[]
+    found_annos = String[]
     for anno_kw in annotation_kws
         annos = parse_annotations(annotation_accessor(model, id_accessor))
-        !haskey(annos, anno_kw) && push!(missing_annos, anno_kw)
+        haskey(annos, anno_kw) && push!(found_annos, anno_kw)
     end
-    return missing_annos
+    return found_annos
 end
 
 """
 $(TYPEDSIGNATURES)
 
 Checks if the databases listed in `config.annotation.gene_annotation_keywords`
 are present in the gene annotations. Returns a vector of annotation keywords
-that were not found.
+that were found.
 """
-findall_unannotated_gene_databases(model, gid::String; config = Config.memote_config) =
-    _find_unannotated_components(
+findall_annotated_gene_databases(model, gid::String; config = Config.memote_config) =
+    _find_annotated_components(
         model,
         gid,
         gene_annotations,
@@ -55,28 +55,25 @@ $(TYPEDSIGNATURES)
 
 Checks if the databases listed in `config.annotation.metabolite_annotation_keywords`
 are present in the metabolite annotations. Returns a vector of annotation keywords
-that were not found.
+that were found.
 """
-findall_unannotated_metabolite_databases(
-    model,
-    mid::String;
-    config = Config.memote_config,
-) = _find_unannotated_components(
-    model,
-    mid,
-    metabolite_annotations,
-    config.annotation.metabolite_annotation_keywords,
-)
+findall_annotated_metabolite_databases(model, mid::String; config = Config.memote_config) =
+    _find_annotated_components(
+        model,
+        mid,
+        metabolite_annotations,
+        config.annotation.metabolite_annotation_keywords,
+    )
 
 """
 $(TYPEDSIGNATURES)
 
 Checks if the databases listed in `config.annotation.reaction_annotation_keywords`
 are present in the reaction annotations. Returns a vector of annotation keywords
-that were not found.
+that were found.
 """
-findall_unannotated_reaction_databases(model, rid::String; config = Config.memote_config) =
-    _find_unannotated_components(
+findall_annotated_reaction_databases(model, rid::String; config = Config.memote_config) =
+    _find_annotated_components(
         model,
         rid,
         reaction_annotations,
@@ -86,23 +83,24 @@ findall_unannotated_reaction_databases(model, rid::String; config = Config.memot
 """
 $(TYPEDSIGNATURES)
 
-Helper function to find all the annotations for a component that do not conform in the model.
+Helper function to find all the annotations for a component that conform to
+recognized patterns in the model.
 """
-function _find_nonconformal_annotations(
+function _find_conformal_annotations(
     model::MetabolicModel,
     id_accessor::String,
     annotation_accessor,
     annotation_regex,
 )
-    no_annos = String[]
+    conform_annos = String[]
     for (anno_kw, anno_regex) in annotation_regex
         annos = parse_annotations(annotation_accessor(model, id_accessor))
         if haskey(annos, anno_kw)
-            in.(nothing, Ref(match.(anno_regex, annos[anno_kw]))) &&
-                push!(no_annos, anno_kw)
+            in.(nothing, Ref(match.(anno_regex, annos[anno_kw]))) ||
+                push!(conform_annos, anno_kw)
         end
     end
-    return no_annos
+    return conform_annos
 end
 
 """
@@ -111,10 +109,10 @@ $(TYPEDSIGNATURES)
 Check if the gene annotation entry conforms to commonly recognized formats of
 annotation database using regex patterns. Uses the database formats listed in
 `config.annotation.gene_annotation_regexes` to test the conformity. Returns a
-string vector of database ids do not conform.
+string vector of database ids that do conform.
 """
-findall_nonconformal_gene_annotations(model, gid::String; config = Config.memote_config) =
-    _find_nonconformal_annotations(
+findall_conformal_gene_annotations(model, gid::String; config = Config.memote_config) =
+    _find_conformal_annotations(
         model,
         gid,
         gene_annotations,
@@ -127,13 +125,13 @@ $(TYPEDSIGNATURES)
 Check if the metabolite annotation entry conforms to commonly recognized formats
 of annotation database using regex patterns. Uses the database formats listed in
 `config.annotation.metabolite_annotation_regexes` to test the conformity.
-Returns a string vector of database ids do not conform.
+Returns a string vector of database ids that do conform.
 """
-findall_nonconformal_metabolite_annotations(
+findall_conformal_metabolite_annotations(
     model,
     mid::String;
     config = Config.memote_config,
-) = _find_nonconformal_annotations(
+) = _find_conformal_annotations(
     model,
     mid,
     metabolite_annotations,
@@ -146,17 +144,14 @@ $(TYPEDSIGNATURES)
 Check if the reaction annotation entry conforms to commonly recognized formats
 of annotation database using regex patterns. Uses the database formats listed in
 `config.annotation.reaction_annotation_regexes` to test the conformity. Returns
-a string vector of database ids do not conform.
+a string vector of database ids that do conform.
 """
-findall_nonconformal_reaction_annotations(
-    model,
-    rid::String;
-    config = Config.memote_config,
-) = _find_nonconformal_annotations(
-    model,
-    rid,
-    reaction_annotations,
-    config.annotation.reaction_annotation_regexes,
-)
+findall_conformal_reaction_annotations(model, rid::String; config = Config.memote_config) =
+    _find_conformal_annotations(
+        model,
+        rid,
+        reaction_annotations,
+        config.annotation.reaction_annotation_regexes,
+    )
 
 end # module

src/memote/checks/Metabolites.jl

@@ -18,33 +18,24 @@ $(TYPEDSIGNATURES)
 Test if metabolites `m1` and `m2` are different by comparing their
 `config.metabolite.test_annotation` field in the annotations of each
 metabolite. Note, if no annotations are present for one or both of the
-metabolites, then return `true`.
+metabolites, then return `false`.
 """
-function metabolites_are_duplicated(
-    model::MetabolicModel,
-    m1,
-    m2;
-    config = Config.memote_config,
-)
-    k1s = get(
-        parse_annotations(metabolite_annotations(model, m1)),
-        config.metabolite.test_annotation,
-        nothing,
-    )
-    isnothing(k1s) && return true
-    k2s = get(
-        parse_annotations(metabolite_annotations(model, m2)),
-        config.metabolite.test_annotation,
-        nothing,
-    )
-    isnothing(k2s) && return true
+function metabolites_are_duplicated(model::MetabolicModel, m1, m2, test_annotation)
+    k1s =
+        get(parse_annotations(metabolite_annotations(model, m1)), test_annotation, nothing)
+    isnothing(k1s) && return false
+    k2s =
+        get(parse_annotations(metabolite_annotations(model, m2)), test_annotation, nothing)
+    isnothing(k2s) && return false
     any(in.(k1s, Ref(k2s)))
 end
 
 """
 $(TYPEDSIGNATURES)
 
-Return a dictionary of metabolites that are duplicated in their compartment.
+Return a dictionary of metabolites that are duplicated in their compartment. If
+any of the test annotations, stored in `config.metabolite.test_annotations`, are
+repeated, then the metabolite is counted as duplicated. Missing annotations are ignored.
 """
 function metabolites_duplicated_in_compartment(
     model::MetabolicModel;
@@ -55,7 +46,12 @@ function metabolites_duplicated_in_compartment(
         c1 = metabolite_compartment(model, m1)
         for m2 in metabolites(model)
             c2 = metabolite_compartment(model, m2)
-            if c1 == c2 && m1 != m2 && metabolites_are_duplicated(model, m1, m2; config)
+            if c1 == c2 &&
+               m1 != m2 &&
+               any(
+                   metabolites_are_duplicated(model, m1, m2, test_annotation) for
+                   test_annotation in config.metabolite.test_annotations
+               )
                 push!(unique_metabolites, m1)
             end
         end
@@ -66,28 +62,34 @@ end
 """
 $(TYPEDSIGNATURES)
 
-Helper function to find orphan or deadend metabolites. Specify `consumed=true`
-to consider orphan metabolites or `false` to consider deadend metabolites. Set
-`complete_medium=true` to open all boundary reactions to simulate a complete
-medium.
+Helper function to find orphan or deadend metabolites. Specify
+`only_consumed=true` to consider orphan metabolites or `false` to consider
+deadend metabolites.
 """
-function _find_orphan_or_deadend_metabolites(model::MetabolicModel; consumed = true)
+function _find_orphan_or_deadend_metabolites(model::MetabolicModel; only_consumed = true)
     mids = metabolites(model)
+    rids = reactions(model)
     mets = String[]
     S = stoichiometry(model)
     lbs, ubs = bounds(model)
     for idx in axes(S, 1)
-        rids, vals = findnz(S[idx, :])
-        if length(vals) == 1 && (consumed ? first(vals) < 0 : first(vals) > 0)
-            ridx = first(rids)
-            rid = reactions(model)[ridx]
-            met = mids[idx]
-            v = reaction_stoichiometry(model, rid)[met]
-            # check if reaction can actually make or consume the metabolite
-            lbs[ridx] < 0 < ubs[ridx] && continue # ignore reversible reactions
-            consumed && lbs[ridx] * v <= 0 && push!(mets, met)
-            !consumed && ubs[ridx] * v >= 0 && push!(mets, met)
+        mid = mids[idx]
+        one_sided_metabolite = true
+        for (ridx, _) in zip(findnz(S[idx, :])...)
+            lbs[ridx] < 0 < ubs[ridx] && (one_sided_metabolite = false; break) # can be produced or consumed
+            rid = rids[ridx]
+            rs = reaction_stoichiometry(model, rid)[mid]
+            if only_consumed # consumed only
+                (rs < 0 && lbs[ridx] >= 0) ||
+                    (rs > 0 && ubs[ridx] <= 0) ||
+                    (one_sided_metabolite = false; break)
+            else # produced only
+                (rs > 0 && lbs[ridx] >= 0) ||
+                    (rs < 0 && ubs[ridx] <= 0) ||
+                    (one_sided_metabolite = false; break)
+            end
         end
+        one_sided_metabolite && push!(mets, mid)
     end
     return mets
 end
@@ -96,20 +98,19 @@ end
 $(TYPEDSIGNATURES)
 
 Find all metabolites that can only (excludes reversible reactions) be consumed
-in the `model` by inspecting the stoichiometric matrix.
+in the `model` by inspecting the stoichiometric matrix and reaction bounds.
 """
 find_orphan_metabolites(model::MetabolicModel) =
-    _find_orphan_or_deadend_metabolites(model, consumed = true)
+    _find_orphan_or_deadend_metabolites(model, only_consumed = true)
 
 """
 $(TYPEDSIGNATURES)
 
 Find all metabolites that can only (excludes reversible reactions) be produced
-in the `model` by inspecting the stoichiometric matrix.
+in the `model` by inspecting the stoichiometric matrix and reaction bounds.
 """
 find_deadend_metabolites(model::MetabolicModel) =
-    _find_orphan_or_deadend_metabolites(model, consumed = false)
-
+    _find_orphan_or_deadend_metabolites(model, only_consumed = false)
 
 """
 $(TYPEDSIGNATURES)
@@ -119,5 +120,4 @@ Returns a list of all metabolites that aren't part of any reactions.
 find_disconnected_metabolites(model::MetabolicModel) =
     metabolites(model)[all(stoichiometry(model) .== 0, dims = 2)[:]]
 
-
 end # module

src/memote/checks/Network.jl

@@ -108,7 +108,7 @@ function find_all_universally_blocked_reactions(
 
     fvas = flux_variability_analysis(
         stdmodel,
-        collect(1:n_reactions(model)),
+        collect(1:n_reactions(stdmodel)),
         optimizer;
         modifications = config.network.optimizer_modifications,
         workers,
@@ -121,7 +121,7 @@ function find_all_universally_blocked_reactions(
         isnothing(vs[2]) && continue
         abs(vs[1]) <= config.network.blocked_tol &&
             abs(vs[2]) <= config.network.blocked_tol &&
-            (push!(blocked_rxns, rid); break)
+            push!(blocked_rxns, rid)
     end
 
     blocked_rxns

src/memote/config.jl

@@ -108,8 +108,8 @@ Base.@kwdef mutable struct AnnotationConfig
     metabolite_annotation_regexes::Dict{String,Regex}
     reaction_annotation_keywords::Vector{String}
     reaction_annotation_regexes::Dict{String,Regex}
-    maximum_nonconformal_references::Int64
-    maximum_missing_databases::Int64
+    minimum_conformal_crossreferences::Int64
+    minimum_crossreferences::Int64
 end
 
 annotation_config = AnnotationConfig(
@@ -190,8 +190,8 @@ annotation_config = AnnotationConfig(
             r"^\d+\.-\.-\.-|\d+\.\d+\.-\.-|\d+\.\d+\.\d+\.-|\d+\.\d+\.\d+\.(n)?\d+$",
         "biocyc" => r"^[A-Z-0-9]+(?<!CHEBI)(\:)?[A-Za-z0-9+_.%-]+$",
     ),
-    maximum_nonconformal_references = 5,
-    maximum_missing_databases = 6,
+    minimum_conformal_crossreferences = 3,
+    minimum_crossreferences = 4,
 )
 
 """
@@ -319,10 +319,11 @@ $(TYPEDEF)
 $(TYPEDFIELDS)
 """
 Base.@kwdef mutable struct MetaboliteConfig
-    test_annotation::String
+    test_annotations::Vector{String}
 end
 
-metabolite_config = MetaboliteConfig(test_annotation = "inchi_key")
+metabolite_config =
+    MetaboliteConfig(test_annotations = ["inchi_key", "chebi", "metanetx.chemical"])
 
 """
 $(TYPEDEF)
@@ -340,7 +341,7 @@ end
 network_config = NetworkConfig(
     condition_number = 1e9,
     cycle_tol = 1e-3,
-    blocked_tol = 1e-3,
+    blocked_tol = 1e-8,
     optimizer_modifications = [silence],
 )