Commit da8d0dc1 authored by eladnoor's avatar eladnoor

fixing command line scripts to new syntax

parent 55f7d771
This diff is collapsed.
......@@ -25,11 +25,11 @@
# THE SOFTWARE.
from equilibrator_api import ComponentContribution, Reaction, ReactionMatcher
from equilibrator_api import ComponentContribution, Reaction, ureg, kJ_per_mol
import logging
import argparse
import csv
from numpy import sqrt, nan
import pandas as pd
import numpy as np
import sys
if __name__ == '__main__':
......@@ -53,37 +53,44 @@ if __name__ == '__main__':
args = parser.parse_args()
sys.stderr.write('pH = %.1f\n' % args.ph)
sys.stderr.write('I = %.1f M\n' % args.i)
p_h = args.ph
ionic_strength = args.i * ureg.molar
sys.stderr.write('pH = %.1f\n' % p_h)
sys.stderr.write('I = %s\n' % ionic_strength)
infile_lines = list(filter(None, map(str.strip, args.infile.readlines())))
reaction_matcher = ReactionMatcher()
if args.plaintext:
reactions = list(map(reaction_matcher.match, infile_lines))
raise Exception("ReactionMatcher needs to be reimplemented")
#reaction_matcher = ReactionMatcher()
#reactions = list(map(reaction_matcher.match, infile_lines))
else:
reactions = list(map(Reaction.parse_formula, infile_lines))
equilibrator = ComponentContribution(pH=args.ph, ionic_strength=args.i)
dG0_prime, U = equilibrator.dG0_prime_multi(reactions)
writer = csv.writer(args.outfile)
header = ['reaction (text)', 'reaction (KEGG)', 'pH', 'ionic strength [M]', 'dG\'0 [kJ/mol]',
'uncertainty [kJ/mol]', 'ln(Reversibility Index)', 'comment']
writer.writerow(header)
for s, r, dg0, u in zip(infile_lines, reactions,
dG0_prime.flat, U.diagonal().flat):
if args.plaintext:
row = [s, r.write_formula()]
else:
row = [reaction_matcher.write_text_formula(r), s]
row += [args.ph, args.i]
if r.check_full_reaction_balancing():
ln_RI = r.calculate_reversibility_index_from_dG0_prime(dg0)
row += ['%.2f' % dg0, '%.2f' % sqrt(u), '%.2f' % ln_RI, '']
else:
row += [nan, nan, nan, 'reaction is not chemically balanced']
writer.writerow(row)
eq_api = ComponentContribution(p_h=p_h, ionic_strength=ionic_strength)
result_df = pd.DataFrame(data=list(zip(infile_lines, reactions)),
columns=["reaction (original)",
"reaction (parsed)"])
result_df['pH'] = p_h
result_df['I [M]'] = float(ionic_strength/ureg.molar)
dG0_prime, U = eq_api.standard_dg_prime_multi(reactions)
result_df["ΔG'° [kJ/mol]"] = [float(dg0/kJ_per_mol) for dg0 in dG0_prime]
result_df["σ(ΔG'°) [kJ/mol]"] = [
float(np.sqrt(u)/kJ_per_mol) for u in U.diagonal().flat]
result_df["ln(Reversibility index)"] = [
float(eq_api.ln_reversibility_index(r)) for r in reactions]
result_df["comment"] = ""
non_balanced = result_df["reaction (parsed)"].apply(
lambda r: r.is_balanced(raise_exception=False) == False
)
result_df.loc[non_balanced, "ln(Reversibility index)"] = np.nan
result_df.loc[non_balanced, "comment"] = "reaction is not chemically balanced"
result_df.to_csv(args.outfile)
args.outfile.flush()
This diff is collapsed.
......@@ -32,15 +32,15 @@
import argparse
import logging
import sys
from equilibrator_api import Reaction, ComponentContribution
from equilibrator_api import Reaction, ComponentContribution, ureg
def MakeParser():
parser = argparse.ArgumentParser(
description=('Estimate the Gibbs energy of a reaction. For example,'
'the following calculates dGr0 for ATP hydrolysis '
'at pH 6: calc_dGr0.py --ph 6 "C00002 + C00001 = '
'C00008 + C00009"'))
'at pH 6: equlibrator_cmd.py --ph 6 "KEGG:C00002 + '
'KEGG:C00001 = KEGG:C00008 + KEGG:C00009"'))
parser.add_argument('--ph', type=float, help='pH level', default=7.0)
parser.add_argument('--i', type=float,
help='ionic strength in M',
......@@ -55,8 +55,11 @@ args = parser.parse_args()
logging.getLogger().setLevel(logging.WARNING)
sys.stderr.write('pH = %.1f\n' % args.ph)
sys.stderr.write('I = %.1f M\n' % args.i)
p_h = args.ph
ionic_strength = args.i * ureg.molar
sys.stderr.write('pH = %.1f\n' % p_h)
sys.stderr.write('I = %s\n' % ionic_strength)
sys.stderr.write('Reaction: %s\n' % args.reaction)
sys.stderr.flush()
......@@ -67,26 +70,26 @@ except ValueError as e:
logging.error(str(e))
sys.exit(-1)
equilibrator = ComponentContribution(pH=args.ph, ionic_strength=args.i)
equilibrator = ComponentContribution(p_h=p_h, ionic_strength=ionic_strength)
n_e = reaction.check_half_reaction_balancing()
if n_e is None:
logging.error('reaction is not chemically balanced')
sys.exit(-1)
elif n_e == 0:
dG0_prime, dG0_uncertainty = equilibrator.dG0_prime(reaction)
sys.stdout.write(u'\u0394G\'\u00B0 = %.1f \u00B1 %.1f kJ/mol\n' %
dG0_prime, dG0_uncertainty = equilibrator.standard_dg_prime(reaction)
sys.stdout.write(u'\u0394G\'\u00B0 = %s \u00B1 %s\n' %
(dG0_prime, dG0_uncertainty))
ln_RI = equilibrator.reversibility_index(reaction)
ln_RI = equilibrator.ln_reversibility_index(reaction)
sys.stdout.write(u'ln(Reversibility Index) = %.1f\n' % ln_RI)
else: # treat as a half-reaction
logging.warning('This reaction isn\'t balanced, but can still be treated'
' as a half-reaction')
E0_prime_mV, E0_uncertainty = equilibrator.E0_prime(reaction)
E0_prime_mV, E0_uncertainty = equilibrator.standard_e_prime(reaction)
sys.stdout.write(
u'E\'\u00B0 = %.1f \u00B1 %.1f mV\n' %
u'E\'\u00B0 = %s \u00B1 %s\n' %
(E0_prime_mV, E0_uncertainty))
sys.stdout.flush()
......
......@@ -31,6 +31,7 @@ del get_versions
from equilibrator_cache import Compound, ureg, kJ_per_mol, R
from equilibrator_cache.compound_cache import ccache
from equilibrator_api.bounds import Bounds
from equilibrator_api.pathway import Pathway
from equilibrator_api.phased_reaction import PhasedReaction as Reaction
from equilibrator_api.component_contribution import ComponentContribution
......
......@@ -79,7 +79,7 @@ class ComponentContribution(object):
reactions,
p_h=self.p_h,
ionic_strength=self.ionic_strength,
T=self.temperature,
temperature=self.temperature,
raise_exception=raise_exception
)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment