<commit_before># Copyright 2015 Novo Nordisk Foundation Center for Biosustainability, DTU.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""Methods for stoichiometric analaysis"""

from copy import copy

import sympy
from cameo.core import SolverBasedModel
from cameo import Reaction
from cobra import Metabolite
import optlang
import six
import logging

logger = logging.getLogger(__name__)
logger.setLevel(logging.DEBUG)

from cameo.exceptions import SolveError, Infeasible
from cameo.util import TimeMachine


def find_dead_end_reactions(model):
    """
    Identify reactions that are structurally prevented from carrying flux (dead ends).
    """
    stoichiometries = {}
    for reaction in model.reactions:
        for met, coef in reaction.metabolites.items():
            stoichiometries.setdefault(met.id, {})[reaction.id] = coef

    blocked_reactions = set()
    while True:
        new_blocked = set()
        for met_id, stoichiometry in stoichiometries.items():
            if len(stoichiometry) == 1:  # Metabolite is only associated with 1 reaction, which can thus not be active
                new_blocked.add(list(stoichiometry)[0])
        if len(new_blocked) == 0:
            break  # No more blocked reactions

        # Remove blocked reactions from stoichiometries
        stoichiometries = {
            met_id: {reac_id: coef for reac_id, coef in stoichiometry.items() if reac_id not in new_blocked}
            for met_id, stoichiometry in stoichiometries.items()
        }
        blocked_reactions.update(new_blocked)

    return blocked_reactions


class ShortestElementaryFluxModes(six.Iterator):
    def __init__(self, model, reactions=None, c=1e-5, copy=True, change_bounds=True):
        self._indicator_variables = None
        if copy:
            self._model = model.copy()
        else:
            self._model = model
        if reactions is None:
            self._reactions = self.model.reactions
        else:
            self._reactions = []
            for reaction in reactions:
                if isinstance(reaction, six.string_types):
                    self._reactions.append(self.model.reactions.get_by_id(reaction))
                else:
                    self._reactions.append(reaction)
        if change_bounds:
            self.__set_exchange_bounds()<commit_msg>chore: add functions for finding dead ends and coupled reactions<commit_after># Copyright 2015 Novo Nordisk Foundation Center for Biosustainability, DTU.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""Methods for stoichiometric analaysis"""

from copy import copy

import sympy
from cameo.core import SolverBasedModel
from cameo import Reaction
from cobra import Metabolite
import optlang
import six
import logging

logger = logging.getLogger(__name__)
logger.setLevel(logging.DEBUG)

from cameo.exceptions import SolveError, Infeasible
from cameo.util import TimeMachine


def find_dead_end_reactions(model):
    """
    Identify reactions that are structurally prevented from carrying flux (dead ends).
    """
    stoichiometries = {}
    for reaction in model.reactions:
        for met, coef in reaction.metabolites.items():
            stoichiometries.setdefault(met.id, {})[reaction.id] = coef

    blocked_reactions = set()
    while True:
        new_blocked = set()
        for met_id, stoichiometry in stoichiometries.items():
            if len(stoichiometry) == 1:  # Metabolite is only associated with 1 reaction, which can thus not be active
                new_blocked.add(list(stoichiometry)[0])
        if len(new_blocked) == 0:
            break  # No more blocked reactions

        # Remove blocked reactions from stoichiometries
        stoichiometries = {
            met_id: {reac_id: coef for reac_id, coef in stoichiometry.items() if reac_id not in new_blocked}
            for met_id, stoichiometry in stoichiometries.items()
        }
        blocked_reactions.update(new_blocked)

    return blocked_reactions


def find_coupled_reactions(model, return_dead_ends=False):
    """Find reaction sets that are structurally forced to carry equal flux"""
    blocked = find_dead_end_reactions(model)
    stoichiometries = {}
    for reaction in model.reactions:
        if reaction.id in blocked:
            continue
        for met, coef in reaction.metabolites.items():
            stoichiometries.setdefault(met.id, {})[reaction.id] = coef

    # Find reaction pairs that are constrained to carry equal flux
    couples = []
    for met_id, stoichiometry in stoichiometries.items():
        if len(stoichiometry) == 2 and set(stoichiometry.values()) == {1, -1}:
            couples.append(set(stoichiometry.keys()))

    # Aggregate the pairs into groups
    coupled_groups = []
    for couple in couples:
        for group in coupled_groups:
            if len(couple & group) != 0:
                group.update(couple)
                break
        else:
            coupled_groups.append(couple)

    if return_dead_ends:
        return coupled_groups, blocked
    else:
        return coupled_groups


class ShortestElementaryFluxModes(six.Iterator):
    def __init__(self, model, reactions=None, c=1e-5, copy=True, change_bounds=True):
        self._indicator_variables = None
        if copy:
            self._model = model.copy()
        else:
            self._model = model
        if reactions is None:
            self._reactions = self.model.reactions
        else:
            self._reactions = []
            for reaction in reactions:
                if isinstance(reaction, six.string_types):
                    self._reactions.append(self.model.reactions.get_by_id(reaction))
                else:
                    self._reactions.append(reaction)
        if change_bounds:
            self.__set_exchange_bounds()
{'fraction_non_alphanumeric': 0.07482692590955958, 'fraction_numerical': 0.004713507143909265, 'mean_word_length': 2.8210467079347215, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '501534', 'n_tokens_mistral': 1889, 'n_tokens_neox': 1827, 'n_words': 608}