DNA production/consumption in cytoplasm

PkiwiBird commented 2 years ago

Description of the issue:

reaction MAR07160, which represents DNA formation is inside the cytoplasm in the model.
reaction MAR07163, which represents DNA degradation by the activity of DNA exodeoxyribonuclease III is inside the cytoplasm in the model.
reaction MAR08639 transports DNA from the cytoplasm to the nucleus and vice-versa. Note: the three reactions have a low confidence score in vmh (score 1).

Expected feature/value/output:

reaction MAR07160 to happen in the nucleus because it is where DNA replication occurs.
reaction MAR07163 to happen in the nucleus because it is where DNA is, and where the DNA exonucleases correct DNA replication errors.
for the reaction MAR08639 to not exist, as DNA does not cross the nucleus envelope.

What could be added to guarantee flux through the reactions (not be blocked):

nuclear dAMP[n] and dGMP[n] could be added together with corresponding transport reactions from the cytoplasm to the nucleus (both ways) for MAR07163 to not be blocked when moved to the nucleus.

Pros of making these changes:

if MAR07160 was moved to the nucleus (all its metabolites were in the nucleus) then reactions MAR07861 (dATP[c] <=> dATP[n]) and MAR07862 (dGTP[c] <=> dGTP[n]) would not be blocked.
if MAR07163 was moved to the nucleus together with the aforementioned supporting reactions, reaction MAR06613 would not be blocked. Note: there is a bibliographic reference for this reaction in metabolic atlas.
the only reaction that would be blocked when applying the changes would be reaction MAR08639.

Reproducing these results:

Code in python:

import cobra
from cobra import Metabolite, Reaction
from urllib.request import urlretrieve

##### load human1.10 and save:
if not os.path.exists('data/models/human1_10.xml'):
    mdPthSBML='https://raw.githubusercontent.com/SysBioChalmers/Human-GEM/master/model/Human-GEM.xml' # path to generic SBML model
    modelPathSBML, _ = urlretrieve(mdPthSBML)
    m = cobra.io.read_sbml_model(modelPathSBML)
    cobra.io.write_sbml_model(m, 'data/models/human1_10.xml')
else:
    m = cobra.io.read_sbml_model('data/models/human1_10.xml')

##### fva to identify blocked reactions:
blck_rc = cobra.flux_analysis.find_blocked_reactions(m)
len(blck_rc) # 1190
'MAR06613' in blck_rc # True
'MAR07861' in blck_rc # True
'MAR07862' in blck_rc # True

#### move reaction MAR07160 to nucleus:
m.reactions.MAR07160.add_metabolites({
'MAM01642n': -0.3,
'MAM01645n': -0.2,
'MAM01688n': -0.2,
'MAM01753n': -0.3,
'MAM01721n': 1.0,
'MAM02759n': 1.0
})
m.reactions.MAR07160.subtract_metabolites({
'MAM01642c': -0.3,
'MAM01645c': -0.2,
'MAM01688c': -0.2,
'MAM01753c': -0.3,
'MAM01721c': 1.0,
'MAM02759c': 1.0
})

#### move reaction MAR07163 to nucleus:
MAM01639n = Metabolite(
    'MAM01639n',
    formula='C10H12N5O6P',
    name='dAMP',
    compartment='n',
    charge= '-2')
MAM01686n = Metabolite(
    'MAM01686n',
    formula='C10H12N5O7P',
    name='dGMP',
    compartment='n',
    charge= '-2')
m.reactions.MAR07163.subtract_metabolites({
'MAM01721c': -1.0,
'MAM02040c': -1.0,
'MAM01639c': 0.3,
'MAM01644c': 0.2,
'MAM01686c': 0.2,
'MAM01752c': 0.3,
'MAM02039c': 1.0
})
m.reactions.MAR07163.add_metabolites({
'MAM01721n': -1.0,
'MAM02040n': -1.0,
MAM01639n: 0.3,
'MAM01644n': 0.2,
MAM01686n: 0.2,
'MAM01752n': 0.3,
'MAM02039n': 1.0
})
trc_damp_cn = Reaction('trc_damp_cn')
trc_damp_cn.add_metabolites({
    m.metabolites.get_by_id('MAM01639c'):-1.0,
    m.metabolites.get_by_id('MAM01639n'):1.0
})
trc_damp_cn.gene_reaction_rule = ''

trc_dgmp_cn = Reaction('trc_dgmp_cn')
trc_dgmp_cn.add_metabolites({
    m.metabolites.get_by_id('MAM01686c'):-1.0,
    m.metabolites.get_by_id('MAM01686n'):1.0
})
trc_dgmp_cn.gene_reaction_rule = ''

m.add_reactions([trc_damp_cn, trc_dgmp_cn])
m.reactions.trc_damp_cn.bounds = (-1000.0, 1000.0)
m.reactions.trc_dgmp_cn.bounds = (-1000.0, 1000.0)

#### see if reactions are blocked now:
blck_rc_end = cobra.flux_analysis.find_blocked_reactions(m)
len(blck_rc_end) # 1188
'MAR06613' in blck_rc_end # False - not blocked anymore
'MAR07861' in blck_rc_end # False - not blocked anymore
'MAR07862' in blck_rc_end # False - not blocked anymore
'MAR08639' in blck_rc_end # True - DNA transport to nucleus is blocked

'MAR07160' in blck_rc_end # False
'MAR07163' in blck_rc_end # False

I hereby confirm that I have:

[x] Tested my code on my own computer for running the model
[x] Done this analysis in the main branch of the repository
[x] Checked that a similar issue does not exist already

Note: replace [ ] with [X] to check the box. PLEASE DELETE THIS LINE

haowang-bioinfo commented 2 years ago

@PkiwiBird nice job! sounds reasonable to me

JonathanRob commented 2 years ago

@PkiwiBird thanks for your work on this! To move forward for implementation, I would suggest the following:

For newly added reactions, use the MAR##### naming scheme
Update the metabolites.tsv and reactions.tsv annotation files accordingly

Also, although I agree that MAR08639 (DNA transport across the nuclear membrane) should probably be removed, it is associated with many genes, most of which encode the nuclear pore complex. It seems that these genes should instead be associated with a reaction facilitating RNA transport across the nuclear membrane, but that would also involve moving some reactions/metabolites around. So not really necessary to address now, but an idea to keep in mind for future work - maybe I'll add it as a new issue.

mihai-sysbio commented 2 years ago

Thank you @PkiwiBird for the PR - it has now been merged, so I am considering the issue resolved. In case I missed something, please re-open the issue, or create a new one as you see fit.

SysBioChalmers / Human-GEM