nskinetics.ReactionSystem#

class nskinetics.ReactionSystem(ID, reactions, species_system)[source]#

Bases: object

Represents a system of one or more chemical reactions involving dynamic changes in species concentrations over time.

This class allows the construction, simulation, and analysis of reaction systems governed by mass action kinetics or other user-defined rate laws. It supports features such as species spikes, event-triggered integration, parameter fitting, and nesting of reaction systems.

Parameters:
  • ID (str) – Unique identifier for the reaction system.

  • reactions (list of str, Reaction, or ReactionSystem) – A list defining the reactions in the system. Each item can be: - A Reaction object - Another ReactionSystem object (for hierarchical systems) - A string describing a chemical equation with kinetic information

  • species_system (SpeciesSystem) – The species system containing all species involved in the reaction system.

__init__(ID, reactions, species_system)[source]#

Methods

C_at_t(t[, species])

Get interpolated species concentrations at a given time.

GUI([close_immediately])

__init__(ID, reactions, species_system)

add_reaction(reaction)

Add a reaction to the system.

change_reaction(index[, ...])

Modify a reaction at a given index.

fit_reaction_kinetic_parameters_to_data(...)

Fit reaction kinetic parameters to experimental time-series data.

from_drawing([ID])

Launch an interactive reaction-network drawing GUI to construct a ReactionSystem visually.

get_dconcs_dt()

Compute the total rate of change of concentrations for all species in the reaction system by summing the contributions from all reactions.

plot_solution([fig, ax, show_events, ...])

Plot the concentrations of selected species over time.

save_solution_df(filename[, save_events_df])

Save the latest simulation solution to an Excel file and cache DataFrames.

set_reaction_kinetic_params(param_vector)

Set kinetic parameters (kf, kb) for all reactions from a flat vector.

solve(t_span[, sp_conc_for_events, events, ...])

Simulate the reaction system and get concentration vs.

Attributes

reaction_kinetic_param_keys

Keys corresponding to the kinetic parameters.

reaction_kinetic_params

Current kinetic parameter values for all reactions.

reactions

List of reactions and nested reaction systems in this reaction system.

reactions_flattened

All reactions in the system, flattened to a single list.

species_system

The SpeciesSystem associated with this reaction system.

C_at_t(t, species=None)[source]#

Get interpolated species concentrations at a given time.

Parameters:
  • t (float or list of floats) – Time or times at which to evaluate concentrations.

  • species (Species, str or None) – Species object or ID of species to evaluate. If None, returns all concentrations.

Returns:

Concentration(s) at time t.

Return type:

float or np.ndarray

GUI(close_immediately=False)[source]#
__repr__()[source]#

Return the formal string representation of the object.

Returns:

Same as __str__.

Return type:

str

__str__()[source]#

Return a string representation of the reaction system.

Returns:

String listing the reaction system and its reactions.

Return type:

str

add_reaction(reaction)[source]#

Add a reaction to the system.

Parameters:

reaction (Reaction, ReactionSystem, or str) – The reaction object or string to add.

change_reaction(index, new_equation_string=None, new_kf=None, new_kb=None)[source]#

Modify a reaction at a given index.

Parameters:
  • index (int) – Index of the reaction to modify.

  • new_equation_string (str, optional) – New chemical equation.

  • new_kf (float, optional) – New forward rate constant.

  • new_kb (float, optional) – New backward rate constant.

fit_reaction_kinetic_parameters_to_data(data, bounds, solve_method='BDF', de_workers=1, de_popsize=200, de_maxiter=10, do_basinhop=True, bh_niter=10, bh_stepsize=0.25, local_method='L-BFGS-B', polish_top_k=7, seed=42, p0=None, all_species_tracked=False, use_only=None, normalize=True, show_output=True, show_progress=False, plot_during_fit=False, call_before_each_solve=None, timeout_solve_ivp=0.5, **kwargs)[source]#

Fit reaction kinetic parameters to experimental time-series data. (i.e., inverse modeling).

This method optimizes the reaction kinetic parameters of the system to best fit provided concentration data over time using nonlinear least squares optimization. The function assumes a shared set of parameters governing all species.

Parameters:
  • data (pandas.DataFrame, dict, str, or list) – A pandas DataFrame, a dictionary, a path to a .csv or .xlsx file, or a list containing any combination of those items. Each item represents data containing time (‘t’) and concentrations of one or more species as columns. Time should be in the column labeled ‘t’; all other columns are interpreted as species concentrations. If data is a list of data items, each data item must have the same format (i.e., the same columns and order of columns).

  • p0 (array_like, optional) – Initial guess for the kinetic parameters. If not provided, the current reaction_kinetic_params will be used if they are finite and not NaN.

  • all_species_tracked (bool, default False) – If True, assumes all species are tracked and uses the full concentration vector when computing predictions. Otherwise, uses only individual species.

  • show_output (bool, default True) – If True, prints fit summary including R² score and success status.

  • use_only (list of str or Species, optional) – List of Species object or species IDs to fit against. If not provided, all species in the input data will be used. It is recommended to include at least one species involved in each reaction in reactions.

  • method (str, default 'Powell') – Optimization method to use for parameter fitting. Passed to scipy.optimize.minimize.

  • n_minimize_runs (int, optional) – Number of local optimization (minimization) runs to perform, each with a different starting point (first uses p0, others use DE output).

  • n_de_runs (int, optional) – Number of differential evolution runs per minimization run to generate good starting guesses; the best DE result is used to initialize the local minimization.

  • minimize_kwargs (dict) – Additional keyword arguments passed to scipy.optimize.minimize.

  • differential_evolution_kwargs (dict) – Additional keyword arguments passed to scipy.optimize.differential_evolution_kwargs.

  • timeout_solve_ivp (float, int, or None, optional) – Enforce timeout of scipy.integrate.solve_ivp when it exceeds this value. Creates and passes a timeout_function as an event to scipy.integrate.solve_ivp. Defaults to 0.5 (seconds).

Returns:

Updates the reaction kinetic parameters of the system in place and stores the fit result in self._fitsol as a tuple: (best-fit parameters, R² score, success flag).

Return type:

None

Notes

  • The fit is performed by normalizing each species’ concentration to its maximum observed value to ensure scale invariance.

  • The objective function maximizes the mean R² score across selected species.

  • This method relies on solve() and C_at_t() methods to simulate the system under trial parameters and look up interpolated results, and on fit_multiple_dependent_variables() to perform optimization.

classmethod from_drawing(ID='rxn_sys')[source]#

Launch an interactive reaction-network drawing GUI to construct a ReactionSystem visually.

Parameters:

ID (str) – Identifier for the new reaction system.

Returns:

A new ReactionSystem object defined by the drawn network.

Return type:

ReactionSystem

get_dconcs_dt()[source]#

Compute the total rate of change of concentrations for all species in the reaction system by summing the contributions from all reactions.

Returns:

Array representing the net rate of change of concentrations for each species.

Return type:

numpy.ndarray

plot_solution(fig=None, ax=None, show_events=True, sps_to_include=None, x_ticks=None, y_ticks=None, auto_ticks=False, show=True)[source]#

Plot the concentrations of selected species over time.

Parameters:
  • show_events (bool, default True) – Whether to annotate event times on the plot.

  • sps_to_include (list of str or Species, optional) – Species to include in the plot. If None, includes all.

  • x_ticks (list or array-like, optional) – Tick values for the x-axis. If None, selected automatically.

  • y_ticks (list or array-like, optional) – Tick values for the y-axis. If None, selected automatically.

property reaction_kinetic_param_keys#

Keys corresponding to the kinetic parameters. Excludes frozen parameters if self._exclude_frozen_params is True.

Returns:

Parameter keys in the format ‘reaction_ID.kf’ or ‘reaction_ID.kb’. Excludes frozen parameters if self._exclude_frozen_params is True.

Return type:

list of str

property reaction_kinetic_params#

Current kinetic parameter values for all reactions. Excludes frozen parameters if self._exclude_frozen_params is True.

Returns:

Array of kinetic parameters (kf and kb) in the system. Excludes frozen parameters if self._exclude_frozen_params is True.

Return type:

numpy.ndarray

property reactions#

List of reactions and nested reaction systems in this reaction system.

Returns:

List of Reaction or ReactionSystem objects.

Return type:

list

property reactions_flattened#

All reactions in the system, flattened to a single list.

Returns:

All individual reactions, including those from nested ReactionSystem objects.

Return type:

list of Reaction

save_solution_df(filename, save_events_df=True)[source]#

Save the latest simulation solution to an Excel file and cache DataFrames.

Parameters:
  • filename (str) – Path to the output Excel file. If None, data is not saved to disk.

  • save_events_df (bool, default True) – Whether to include event data in the output.

set_reaction_kinetic_params(param_vector)[source]#

Set kinetic parameters (kf, kb) for all reactions from a flat vector.

Parameters:

param_vector (array_like) – Flat list or array of parameters.

solve(t_span, sp_conc_for_events=None, events=None, spikes=None, method='BDF', t_eval=None, atol=1e-06, rtol=0.0001, dense_output=False, y0=None, dt_spike=1e-06, remove_negative_concs=True, filename=None, save_events_df=True)[source]#

Simulate the reaction system and get concentration vs. time data.

Parameters:
  • t_span (list or tuple) – Two-item iterable consisting of the desired start time and end time.

  • sp_conc_for_events (dict, optional) – Keys are Species objects or Species ID strings. Values are concentrations. When a Species concentrations achieves the corresponding value, the time at which this event occurs is stored in ReactionSystem._solution[‘t_events’] and the concentrations array at that time is stored in ReactionSystem._solution[‘y_events’]. Event indices are stored in ReactionSystem._solution[‘events’].

  • events (callable, list or tuple, optional) – Function or list of functions that accepts arguments t (time) and concentrations (ndarray of species concentrations ordered the same as ReactionSystem.species_system.all_sps). Each function represents an event (triggered when the function returns zero), When the function returns zero, the time at which this event event occurs is stored in ReactionSystem._solution[‘t_events’] and the concentrations array at that time is stored in ReactionSystem._solution[‘y_events’]. Event indices are stored in ReactionSystem._solution[‘events’].

  • spikes (dict, optional) –

    A dictionary with time stamps (t) as keys. Values for spikes can assume one of the following forms:

    1. string or list of strings each of the form f’Change; {species_ID}; {conc}’ (e.g., ‘Change; Substrate; 0.001’), where conc is the desired change (spike) in the concentration of the corresponding Species at time t;

    2. string or list of strings each of the form f’Target; {species_ID}; {conc}’ (e.g., ‘Target; Substrate; 0.001’), where conc is the concentration of the corresponding Species to be targeted via a change (spike) at time t;

    3. string or list of strings each of the form f’Target (negative change allowed); {species_ID}; {conc}’ (e.g., ‘Target (negative change allowed); Substrate; 0.001’), where conc is the concentration of the corresponding Species to be targeted via a change (spike) at time t, with negative values for the change (i.e., dips) permitted;

    4. array of floats equal to the change (spike) in concentrations of all Species in the species_system (0 if none) at time t; or

    5. function that accepts the time and current concentrations array as arguments, and returns the desired change (spike) in concentrations of all Species in the species_system (0 if none) at time t.

  • filename (str, optional) –

    Filename or path including filename to save results. Results will be saved to an Excel file in the following format: Sheet: Main

    t [Species ID 1] [Species ID 2] … - ———— ————– . . . … . . . … . . . …

    Sheet: Events

    event t_event [Species ID 1] at t_event [Species ID 2] at t_event … —– ——- ————————- ————————- . . . . … . . . . … . . . . …

Examples

>>> import nskinetics as nsk
>>> sp_sys = nsk.SpeciesSystem('sp_sys',
>>>                        ['E', 'S', 'ES', 'P'], # enzyme, substrate, enzyme-substrate complex, product
>>>                        concentrations=[1e-4, 1e-4, 0., 0.])
>>> reactions = [
>>>             'E + S <-> ES; kf = 12.0, kb = 10.0', # kf = kon, kb = koff
>>>             'ES -> E + P; kf = 32.0' # kf = kcat (enzyme turnover number)
>>>             ]
>>> rxn_sys = nsk.ReactionSystem(ID='rxn_sys',
>>>                                  reactions=reactions,
>>>                                  species_system=sp_sys)
>>> rxn_sys.solve(t_span=[0, 2*24*3600], # I want to simulate the system over 2 days
>>>                  sp_conc_for_events={'S':1e-6}, # In addition to a full simulation,
>>>                  )                              # I want to know the time at which [S] drops to 1e-6
>>> np.allclose(rxn_sys._solution['t_events'], np.array([42219.44616989]), rtol=1e-5, atol=1e-8)
True
>>> np.allclose(rxn_sys._solution['y_events'],
>>>        np.array([[9.99998909e-05],
>>>                 [1.00000000e-06],
>>>                 [1.09091871e-10],
>>>                 [9.89998909e-05]]),
>>>                 rtol=1e-5, atol=1e-8)
True
property species_system#

The SpeciesSystem associated with this reaction system.

Returns:

Object containing all species in the system.

Return type:

SpeciesSystem