Commit 6bb309b8 authored by Elad Noor's avatar Elad Noor

merge

parents 018cd73b 532395b2
Pipeline #81248844 passed with stages
in 7 minutes and 44 seconds
......@@ -66,7 +66,7 @@ class BaseBounds(object):
"""
raise NotImplementedError
def GetLowerBounds(self, compounds: Iterable[Compound]) -> Iterable[float]:
def GetLowerBounds(self, compounds: Iterable[Compound]) -> Iterable[Q_]:
"""Get the bounds for a set of keys in order.
:param compounds: an iterable of Compounds
......@@ -74,7 +74,7 @@ class BaseBounds(object):
"""
return map(self.GetLowerBound, compounds)
def GetUpperBounds(self, compounds: Iterable[Compound]) -> Iterable[float]:
def GetUpperBounds(self, compounds: Iterable[Compound]) -> Iterable[Q_]:
"""Get the bounds for a set of keys in order.
:param compounds: an iterable of Compounds
......@@ -82,7 +82,7 @@ class BaseBounds(object):
"""
return map(self.GetUpperBound, compounds)
def GetBoundTuple(self, compound: Compound) -> Tuple[float, float]:
def GetBoundTuple(self, compound: Compound) -> Tuple[Q_, Q_]:
"""Get both upper and lower bounds for this key.
:param compound: a Compound object
......@@ -92,7 +92,7 @@ class BaseBounds(object):
def GetBounds(
self, compounds: Iterable[Compound]
) -> Tuple[Iterable[float], Iterable[float]]:
) -> Tuple[Iterable[Q_], Iterable[Q_]]:
"""Get the bounds for a set of compounds.
:param compounds: an iterable of Compounds
......@@ -104,13 +104,13 @@ class BaseBounds(object):
@staticmethod
@ureg.check("[concentration]")
def conc2ln_conc(b: float) -> float:
def conc2ln_conc(b: Q_) -> float:
"""Convert a concentration to log-concentration.
:param b: a concentration
:return: the log concentration
"""
return log(b / standard_concentration)
return log(b / standard_concentration).magnitude
def GetLnBounds(
self, compounds: Iterable[Compound]
......@@ -147,7 +147,7 @@ class BaseBounds(object):
return map(self.conc2ln_conc, ubs)
@ureg.check(None, None, "[concentration]", "[concentration]")
def SetBounds(self, compound: Compound, lb: float, ub: float):
def SetBounds(self, compound: Compound, lb: Q_, ub: Q_):
"""Set bounds for a specific key.
Args:
......@@ -170,16 +170,16 @@ class Bounds(BaseBounds):
@ureg.check(None, None, None, "[concentration]", "[concentration]")
def __init__(
self,
lower_bounds: Dict[Compound, float] = None,
upper_bounds: Dict[Compound, float] = None,
default_lb: float = default_conc_lb,
default_ub: float = default_conc_ub,
lower_bounds: Dict[Compound, Q_] = None,
upper_bounds: Dict[Compound, Q_] = None,
default_lb: Q_ = default_conc_lb,
default_ub: Q_ = default_conc_ub,
) -> object:
"""Initialize the bounds object.
Args:
lower_bounds: a dictionary mapping keys to float lower bounds
upper_bounds: a dictionary mapping keys to float upper bounds
lower_bounds: a dictionary mapping keys to lower bounds
upper_bounds: a dictionary mapping keys to upper bounds
default_lb: default lower bound to if no specific one is provided
default_lb: default upper bound to if no specific one is provided
"""
......@@ -198,8 +198,8 @@ class Bounds(BaseBounds):
def from_csv(
cls,
f: TextIO,
default_lb: float = default_conc_lb,
default_ub: float = default_conc_ub,
default_lb: Q_ = default_conc_lb,
default_ub: Q_ = default_conc_ub,
) -> object:
"""Read bounds from a .csv file.
......@@ -212,8 +212,8 @@ class Bounds(BaseBounds):
provided
:return: a Bounds object
"""
lbs: Dict[Compound, float] = dict()
ubs: Dict[Compound, float] = dict()
lbs: Dict[Compound, Q_] = dict()
ubs: Dict[Compound, Q_] = dict()
bounds_df = pd.read_csv(f)
for row in bounds_df.itertuples():
......@@ -247,8 +247,8 @@ class Bounds(BaseBounds):
cls,
df: pd.DataFrame,
bounds_unit: str = None,
default_lb: float = default_conc_lb,
default_ub: float = default_conc_ub,
default_lb: Q_ = default_conc_lb,
default_ub: Q_ = default_conc_ub,
) -> Tuple[object, Dict[Compound, str]]:
"""Read bounds from a Pandas DataFrame.
......@@ -329,7 +329,7 @@ class Bounds(BaseBounds):
new_ub = dict(self.upper_bounds.items())
return Bounds(new_lb, new_ub, self.default_lb, self.default_ub)
def GetLowerBound(self, compound: Compound) -> float:
def GetLowerBound(self, compound: Compound) -> Q_:
"""Get the lower bound for this key.
:param compound: a compound
......
......@@ -33,6 +33,7 @@ from component_contribution import GibbsEnergyPredictor
from . import (
FARADAY,
Q_,
R,
ccache,
default_I,
......@@ -55,10 +56,10 @@ class ComponentContribution(object):
@ureg.check(None, None, None, "[concentration]", "[temperature]")
def __init__(
self,
p_h: float = default_pH,
p_mg: float = default_pMg,
ionic_strength: float = default_I,
temperature: float = default_T,
p_h: Q_ = default_pH,
p_mg: Q_ = default_pMg,
ionic_strength: Q_ = default_I,
temperature: Q_ = default_T,
) -> object:
"""Create a ComponentContribution object.
......@@ -73,7 +74,7 @@ class ComponentContribution(object):
self.temperature = temperature
@property
def RT(self) -> float:
def RT(self) -> Q_:
"""Get the value of RT."""
return R * self.temperature
......@@ -293,3 +294,32 @@ class ComponentContribution(object):
return ComponentContribution.predictor.is_using_group_contribution(
reaction
)
@ureg.check(None, None, None, None, "[length]**2 * [mass] / [time]**2")
def multicompartmental_standard_dg_prime(
self,
reaction: PhasedReaction,
nH: int,
charge: int,
delta_p_h: ureg.Measurement,
delta_phi: ureg.Measurement,
) -> ureg.Measurement:
"""Calculate the transformed reaction energies of a
multi-compartmental reaction.
Based on the equations from
Harandsdottir et al. 2012 (https://doi.org/10.1016/j.bpj.2012.02.032)
:param reaction:
:param nH: the number of hydrogens being transported
:param charge: the net charge being transported
:param delta_p_h: the difference in pH between the two sides
:param delta_phi: the difference in electro-static potential between
the two sides
:return: the transport reaction Gibbs free energy change
"""
dg = self.standard_dg_prime(reaction)
dg -= nH * self.RT * np.log(10.0) * delta_p_h
dg -= FARADAY * charge * delta_phi
return dg
......@@ -104,7 +104,7 @@ class Pathway(object):
if standard_dg_primes is None:
assert dg_sigma is None, (
"If standard_dg_primes are not "
"provided, dg_sigme must also be None"
"provided, dg_sigma must also be None"
)
self.set_aqueous_params(
p_h=p_h,
......@@ -393,7 +393,7 @@ class Pathway(object):
name_to_compound.get, row.ReactionFormula
)
rxn.rid = row.ID
if not rxn.is_balanced(ignore_atoms=("H")):
if not rxn.is_balanced(ignore_atoms=("H",)):
warnings.warn(f"Reaction {rxn.rid} is not balanced")
reactions.append(rxn)
if row.ID in reaction_ids:
......
......@@ -104,7 +104,7 @@ class Condition(object):
I.e. the phase and the abundance.
"""
def __init__(self, phase: str, abundance: ureg.Quantity = None):
def __init__(self, phase: str, abundance: Q_ = None):
"""Create a new condition object."""
assert phase in PHASE_INFO_DICT, f"Unknown phase: {phase}"
self._phase = phase
......@@ -116,17 +116,17 @@ class Condition(object):
return self._phase
@property
def abundance(self) -> ureg.Quantity:
def abundance(self) -> Q_:
"""Return the abundance."""
return self._abundance
@property
def standard_abundance(self) -> ureg.Quantity:
def standard_abundance(self) -> Q_:
"""Return the standard abundance in this phase."""
return PHASE_INFO_DICT[self._phase].standard_abundance
@property
def physiolical_abundance(self) -> ureg.Quantity:
def physiolical_abundance(self) -> Q_:
"""Return the default physiological abundance in this phase."""
return PHASE_INFO_DICT[self._phase].physiolical_abundance
......@@ -148,7 +148,7 @@ class Condition(object):
if self._abundance is None:
return 0.0 # when the phase does not allow relative abundances
else:
return np.log(float(self._abundance / self.standard_abundance))
return np.log(self._abundance / self.standard_abundance).magnitude
@property
def ln_physiological_abundance(self) -> float:
......@@ -158,12 +158,12 @@ class Condition(object):
if self.physiolical_abundance is None:
return 0.0 # when the phase does not allow relative abundances
else:
return np.log(
float(self.physiolical_abundance / self.standard_abundance)
)
return (
np.log(self.physiolical_abundance / self.standard_abundance)
).magnitude
@abundance.setter
def abundance(self, abundance: ureg.Quantity) -> None:
def abundance(self, abundance: Q_) -> None:
"""Change the phase of a specific compound.
:param abundance: the new abundance in the correct units
......@@ -287,12 +287,12 @@ class PhasedCompound(object):
return NON_AQUEOUS_COMPOUND_DICT.get(self.inchi, (AQUEOUS_PHASE_NAME,))
@property
def abundance(self) -> ureg.Quantity:
def abundance(self) -> Q_:
"""Get the abundance."""
return self.condition.abundance
@abundance.setter
def abundance(self, abundance: ureg.Quantity) -> None:
def abundance(self, abundance: Q_) -> None:
"""Change the phase of a specific compound.
:param abundance: the new abundance in the correct units
......@@ -330,11 +330,8 @@ class PhasedCompound(object):
@ureg.check(None, None, "[concentration]", "[temperature]")
def get_stored_standard_dgf_prime(
self,
p_h: ureg.Quantity,
ionic_strength: ureg.Quantity,
temperature: ureg.Quantity,
) -> ureg.Quantity:
self, p_h: Q_, ionic_strength: Q_, temperature: Q_
) -> Q_:
"""Return the stored formation energy of this phased compound.
Only if it exists, otherwise return None (and we will use
......@@ -365,7 +362,7 @@ class PhasedCompound(object):
# for all other compounds, we will use component_contribution
return None
def get_stored_standard_dgf(self) -> ureg.Quantity:
def get_stored_standard_dgf(self) -> Q_:
"""Return the stored formation energy of this phased compound.
Only if it exists, otherwise return None (and we will use
......
......@@ -65,7 +65,12 @@ class PathwayMDFData(object):
:param reaction_prices: shadow prices of reactions
:param compound_prices: shadow prices of compound concentrations
"""
self.mdf = B * R * default_T
if pathway.comp_contrib is not None:
RT = pathway.comp_contrib.RT
else:
RT = R * default_T
self.mdf = B * RT
self.bounds = pathway._bounds.Copy()
concentrations = np.exp(lnC) * standard_concentration
......@@ -78,7 +83,6 @@ class PathwayMDFData(object):
optimized_dg_prime += pathway.dg_sigma @ y
# adjust dG to flux directions and convert to kJ/mol
RT = R * pathway.comp_contrib.temperature
standard_dg_primes = (I_dir @ standard_dg_primes) * RT
physiological_dg_prime = (I_dir @ physiological_dg_prime) * RT
optimized_dg_prime = (I_dir @ optimized_dg_prime) * RT
......
......@@ -77,7 +77,7 @@ def test_water_balancing(comp_contribution, missing_h2o):
assert missing_h2o.is_balanced(ignore_atoms=("H", "O", "e-"))
balanced_rxn = missing_h2o.balance_with_compound(
ccache.get_compound("kegg:C00001"), ignore_atoms=("H")
ccache.get_compound("kegg:C00001"), ignore_atoms=("H",)
)
assert balanced_rxn is not None
......
......@@ -43,13 +43,16 @@ from equilibrator_api import (
(
"kegg:C00002",
"InChI=1S/C10H16N5O13P3/c11-8-5-9(13-2-12-8)15(3-14-5)10-7(17)6("
"16)4(26-10)1-25-30(21,22)28-31(23,24)27-29(18,19)20/h2-4,6-7,10,16-17H,1H2,(H,21,22)(H,23,24)(H2,11,12,13)(H2,18,19,20)/p-4/t4-,6-,7-,10-/m1/s1",
"16)4(26-10)1-25-30(21,22)28-31(23,24)27-29(18,19)20/h2-4,6-7,10,"
"16-17H,1H2,(H,21,22)(H,23,24)(H2,11,12,13)(H2,18,19,20)/p-4/t4-,6-"
",7-,10-/m1/s1",
"ZKHQWZAMYRWXGA-KQYNXXCUSA-J",
),
(
"metanetx.chemical:MNXM7",
"InChI=1S/C10H15N5O10P2/c11-8-5-9(13-2-12-8)15(3-14-5)10-7(17)6("
"16)4(24-10)1-23-27(21,22)25-26(18,19)20/h2-4,6-7,10,16-17H,1H2,(H,21,22)(H2,11,12,13)(H2,18,19,20)/p-3/t4-,6-,7-,10-/m1/s1",
"16)4(24-10)1-23-27(21,22)25-26(18,19)20/h2-4,6-7,10,16-17H,1H2,"
"(H,21,22)(H2,11,12,13)(H2,18,19,20)/p-3/t4-,6-,7-,10-/m1/s1",
"XTWYTFMLZFPYCI-KQYNXXCUSA-K",
),
],
......
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