Essays, Posts & Presentations

Analyzing Multiway Chemical Reaction Networks Using Chemical Space Networks

Nicholas Frieler

Author

Nicholas Frieler

Title

Analyzing Multiway Chemical Reaction Networks Using Chemical Space Networks

Description

Introduction

Visualization is one of the most important tools when it comes to data analysis. It allows researchers to see clear trends in the data they are working with, and different visualization techniques allow for different features of the data to be analyzed. In the past decade, chemical space networks have emerged as a powerful visualization tool for chemists and biologists. A chemical space network is a collection of molecules, each represented as a node, with an edge between two molecules if they are sufficiently similar. Chemical space networks are very versatile as they can be used on systems with a very large number of molecules, and they can be used with any similarity metric, allowing for different features of molecules to be analyzed (1).

In the world of computation, exhaustive searching is an essential tool for fully exploring a given system. In the Wolfram Language, this can be achieved with multiway graphs, which take in a list of states and iteratively apply a set of rules for a specified number of steps. The generality of multiway computation allows it to be applied in a wide array of scenarios. Previously, Ariana-Dalia Vlad completed a project developing a multiway model of chemical reactions. This allows for a great number of related molecules to be generated based on only a few reactions and a small number of starting molecules (2). In a multiway chemical system, all molecules present will be structurally related in some way, with the level of similarity depending on the reactions used. This makes the molecules from any multiway chemical reaction network a perfect candidate for analyzing using chemical space networks.

This project is heavily inspired by the work done by Robert Nachbar, especially what he presented in his 2022 Wolfram Technology Conference talk “Modeling Chemical Reaction Networks with Token-Event Graphs.” In fact, his talk presented multiway analysis on monoterpene biosynthesis, which is the main chemical system explored in this post (3).

MoleculePlotOutputForm

This is the code that automatically displays Molecule objects using MoleculePlot rather than the typical descriptor box. See the output generated below this cell for information on how to disable this temporarily.

In[43]:=

moleculeImageSize[mol_Molecule,scale_:0.5]:= Module[{atoms,bonds,coords,dist,minmax}, atoms=mol["AtomList"]/.Atom[a_,___]:>a// AssociationThread[Range@Length@#,#]&; bonds=Select[mol["BondList"]/.Bond[ij_,___]:>ij, FreeQ[atoms/@#,"H"]&]; coords=QuantityMagnitude@mol["AtomDiagramCoordinates"]; dist=EuclideanDistance@@coords[[#]]&/@bonds// Replace[{}->{0}]//Mean; minmax=coords[[#]]&/@bonds//Catenate//Map[First]//MinMax; Max[36,Min[672,(Subtract@@Reverse@minmax)/Max[1,dist]scale72]] ]$RBNMoleculeFormat=True;Print["Set $RBNMoleculeFormat = False to restore built-in behavior."];Chemistry`Formatting`moleculeBox[mol_,fmt_]/;TrueQ[$RBNMoleculeFormat]:= If[Chemistry`Common`has3DCoordinates[mol], With[{graphics=MoleculePlot[mol]}, {plot=ToBoxes[graphics,fmt]}, Replace[plot,{HoldPattern[Graphics3DBox[x_,y___]]:> Graphics3DBox[TagBox[x,BoxForm`ChemistryTag[HoldComplete[mol]]],y, ContentSelectable->False, Selectable->False, DefaultBaseStyle->{FrontEnd`GraphicsHighlightColor-> RGBColor[0.269281,0.535948,0.63268]}], _:>$Failed} ] ], With[{graphics=MoleculePlot[mol,ImageSize->moleculeImageSize[mol,0.375]]}, {plot=ToBoxes[graphics,fmt]}, Replace[plot,{HoldPattern[GraphicsBox[x_,y___]]:> GraphicsBox[TagBox[x,BoxForm`ChemistryTag[HoldComplete[mol]]],y, ContentSelectable->False, Selectable->False, DefaultBaseStyle->{FrontEnd`GraphicsHighlightColor-> RGBColor[0.269281,0.535948,0.63268]}], _:>$Failed} ] ] ]

Set $RBNMoleculeFormat = False to restore built-in behavior.

Pattern Reaction Paradigms

The

PatternReaction

function utilizes graph rewriting to model chemical reactions. There are several options for specifying a chemical reaction using PatternReaction, but the most powerful of these, by far, is a SMARTS string. SMARTS is an acronym that stands for SMILES (simplified molecular-input line-entry system) arbitrary target specification. The power of this method comes from allowing the user to be as specific or general as they wish to be.

For example, a condensation reaction involves the formation of a bond between two hydroxyl groups, forming an ether and giving off water as a byproduct. Utilizing SMARTS allows the user to ignore all other atoms in each molecule (or the singular molecule if it is an intramolecular condensation reaction) and only specify the atoms directly relevant to the reaction.

In[47]:=

generalCondensationReaction=PatternReaction["[*:1][OH1:2].[OH1:3][*:4]>>[*:1][O:2][*:4].[OH2:3]"]

Out[47]=

SMARTS also provides the ability to specify characteristics of the atoms that are not directly affected by the reaction. For example, we could specify that one of the hydroxyls is bonded to a nitrogen atom that is present in a 6-membered ring and that the other hydroxyl is bonded to a carbon atom that is present in exactly 2 rings and that these two must be in two separate molecules before the reaction takes place.

In[48]:=

specificCondensationReaction=PatternReaction["([N;r6:1][OH1:2]).([OH1:3][C;R2:4])>>([N:1][O:2][C:4]).([OH2:3])"]

Out[48]=

Although

PatternReaction

does not show that we have specified the ring properties for the carbon and nitrogen atoms, we can test it it on a series of molecules using the

ApplyReaction

function. Here we have 4 sets of molecules, all of which fit the pattern for

generalCondensationReaction

and only the 4th one fitting the pattern for

specificCondensationReaction

In[49]:=

(patternReactionMoleculeSets=Map[Molecule,{{"NO","CO"},{"C1CCCCN1O","CO"},{"NO","C1CCCC2CCCCC21O"},{"C1CCCCN1O","C1CCCC2CCCCC21O"}},{2}])//Column

Out[49]=

















By applying

generalCondensationReaction

to all 4 sets of molecules, we expect the reactions to be carried out and receive 4 sets of reactants.

In[50]:=

Row[{Plus@@#,"⟶",Plus@@ApplyReaction[generalCondensationReaction,#]},Spacer[9]]&/@patternReactionMoleculeSets//Column

Out[50]=

⟶

By applying

specificCondensationReaction

to all 4 sets of molecules, we expect only the 4th set to undergo the reaction, with the first 3 reactions each returning an empty list.

In[51]:=

Row[{Plus@@#,"⟶",Plus@@ApplyReaction[specificCondensationReaction,#]/.(0"∅")},Spacer[9]]&/@patternReactionMoleculeSets//Column

Out[51]=

⟶

∅

⟶

∅

⟶

∅

⟶

As we can see above,

specificCondensationReaction

was applied to all four sets of reactants because each reactant contained at least one hydroxyl group. However, only the fourth set of reactants specified the requirements of specificCondensationReaction and thus was the only reaction carried out.

Another benefit of

PatternReaction

is that it does not require conservation of mass or atomic identity. That is, we are able to turn lead into gold!

In[52]:=

PatternReaction["[Pb]>>[Au]"]

Out[52]=

While this might seem silly at first, this allows for simpler specification of some reactions and can speed up applying reactions to large sets of molecules as checking if a molecule satisfies the requirements of a SMARTS string is computationally costly. In the following sections, we explore a few paradigms of

PatternReaction

that utilize slightly different implementations.

Construction

The construction paradigm involves taking a starter molecule and adding onto it without specifying what is added as a reactant. This is particularly useful for building polymers using

ApplyReaction

, especially if we wish to restrict the binding sites for polymerization. For example, if we take a glucose molecule (note that we have not specified the stereochemistry here although we are able to do so), we can use the construction paradigm to specify that every time an additional glucose molecule is added, its 1-carbon must attach to the 6-carbon of one of the already-present moieties.

In[53]:=

glucose=Molecule["O1C(O)C(O)C(O)C(O)C1CO"];glucosePolymerization=PatternReaction["[C;!R:1][OH1:2]>>[C:1][O:2][C]1[C]([O])[C]([O])[C]([O])[C]([C]([O]))[O]1"];glucoseTetramer=ApplyReaction[glucosePolymerization,ApplyReaction[glucosePolymerization,ApplyReaction[glucosePolymerization,glucose]]];Column@{glucose,glucosePolymerization,glucoseTetramer}

Out[56]=

Notice that we have conveniently lost track of the resultant water molecule from this condensation reaction. This allows for less processing between polymerization steps, which would help simplify code carrying out polymerization in this way.

Self-Interaction

The next

PatternReaction

paradigm we will discuss is self-interaction. As the name suggests, self-interaction takes two parts of the same molecule and reacts them in some way. For example, if we take the tetramer from the previous section, there is now enough flexibility to allow for a condensation reaction to occur between two of the hydroxyls of the tetramer. A condensation reaction of this type will produce a larger ring structure. There are many ways for this reaction to occur.

In[57]:=

intramolecularCondensation=PatternReaction["([*:1][OH:2].[OH][*:3])>>[*:1][O:2][*:3]"];ringClosedGlucoseTetramer=ApplyReaction[intramolecularCondensation,glucoseTetramer,All][[34,1]];Column@{intramolecularCondensation,{ringClosedGlucoseTetramer,ringClosedGlucoseTetramer//MoleculePlot3D}}

Out[59]=





One problem that occurs with this implementation is a reaction between parts of the molecule that are too close to each other. While this is no problem for

ApplyReaction

, some of the resulting molecules do not make sense chemically. To remedy this problem, we can repeatedly use the wildcard character “ * ” in SMARTS to represent a series of unknown atoms. There are two problems with this method, however. The number of wildcards used must exactly match the number of atoms, i.e. there is no way to represent a range of unknown numbers of atoms without explicitly writing out each scenario. The other problem is that the number of wildcards corresponds to a possible path between two atoms, not necessarily the shortest. So between two atoms in a different molecule, there could be paths of 12, 15, and 18 atoms. Because we care about atoms being too close to each other we may want to exclude shortest paths of length 12 or less, but these two atoms could still react as there are longer possible paths between the two atoms.

Combination

The combination paradigm is the main way in which

PatternReaction

was originally intended to be used. In this method, we react two (or more) separate molecules to get a list of products, while maintaining atomic count for each element present. One example is that of the Diels-Alder reaction, in which a diene and a dienophile react together to form a ring, pictured below.

	+		⟶
diene		dienophile		product

This reaction can be represented as follows:

In[60]:=

PatternReaction["([C:1]=[C:2][C:3]=[C:4]).([C:5]=[C:6])>>[C:1]1[C:2]=[C:3][C:4][C:5][C:6]1"]

Out[60]=

This method is especially useful in contexts of repeated reactions, which will be discussed in more detail in the following section on “Multiway Chemical Reaction Networks.” One of the benefits of this method is that now every reactant and product atom can be given a label and kept track of during the reaction. This helps reduce the possibility of losing or adding atoms unintentionally or putting atoms in the wrong place in the products. However, as alluded to earlier, specifying every atom in these reactions can be computationally costly for reactions involving much larger reactants.

Multiway Chemical Reaction Network

Explicit reaction pathways are of great importance in chemistry for building molecules from other more accessible or more easily manipulable molecules. Some molecules of interest require a long reaction pathway and can often be reached via multiple pathways. This feature of some reaction pathways lends itself to multiway computational analysis. Multicomputation can be used with any of the previously discussed PatternReaction paradigms, but we shall focus on the self-interaction method.

For chemical reactions, we utilize the

TokenEventGraph

resource function as this allows us to visualize the reactants, products, and the reactions themselves. The reactants and products are represented with tokens and each reaction is an event in the graph (a small circle). Pictured below is a multiway graph featuring monoterpenes and several rearrangement reactions taken two layers deep.

This is all wrapped in the function

MultiwayChemicalReactionNetwork

, which we use both to generate the set of molecules that will later be analyzed and to display the structural similarities of those molecules.

In[61]:=

MultiwayChemicalReactionNetwork[reactions_,startingMolecules_,numSteps_,OptionsPattern[]]:=Module[{},ResourceFunction["TokenEventGraph"][{{mol1_}:>Catenate@(Catenate[ApplyReaction[#,{mol1},All]]&/@reactions)},{startingMolecules},numSteps,"EventVertices"->OptionValue["EventVertices"],"EventDeduplication"->OptionValue["EventDeduplication"],"TokenDecompositionFunction"->(DeleteDuplicates[#,MoleculeMatchQ]&),"TokenDeduplication"->OptionValue["TokenDeduplication"],"TokenEquivalenceFunction"->OptionValue["TokenEquivalenceFunction"],"TokenRenderingFunction"OptionValue["TokenRenderingFunction"],"TokenLabeling"->OptionValue["TokenLabeling"],"MultiOutputEvents"->OptionValue["MultiOutputEvents"]]]

Note that the above function currently only works with reaction sets involving self-interaction reactions only, so more work needs to be done to generalize this to reaction sets with varying numbers of reactants.

In order to work with the molecules generated by

MultiwayChemicalReactionNetwork

, we extract the molecules from the multiway graph and access the

SMILES

attribute, which returns the SMILES string for each molecule.

In[62]:=

GetMultiwayChemicalReactionNetworkSMILES[multiwayGraph_]:=Module[{reactedMolecules,reactedMoleculeSMILES},#["SMILES"]&/@VertexList@multiwayGraph]

This also provides a convenient and data-efficient format for storing all of the molecules generated by the multiway graph.

Chemical Space Network

A chemical space network is a graph where each node corresponds to an individual molecule and two nodes are connected by an edge if they have a similarity score (ranging between 0 and 1 inclusive) above a specified threshold. Depending upon the context, there are many possible choices for a similarity metric. All of the analysis in this post utilizes the TopologicalFingerprint function from WolframChemistry and compares these fingerprints using the JaccardDissimilarity function.

The code below is the first step in creating a chemical space network. Taking the list of SMILES strings from the multiway chemical reaction network as input, this function returns an association between pairs of molecules and their similarity scores. The pairs are stored as undirected edges to make the creation of graphs later on easier.

In[63]:=

ChemicalSpaceNetwork[moleculeSMILES_]:=Module[{moleculeFingerprints},moleculeFingerprints=(#->TopologicalFingerprint[Molecule[#]]&/@moleculeSMILES)//Association;ParallelMap[UndirectedEdge[First[#],Last[#]]->N@(1-(JaccardDissimilarity@@{moleculeFingerprints[First[#]],moleculeFingerprints[Last[#]]}))&,Subsets[Normal[moleculeSMILES],{2}]]//Association]

This function is the only one in this analysis that makes use of the

TopologicalFingerprint

and the

JaccardDissimilarity

functions, making it easy to implement a different similarity measure for future analyses.

Statistical Analysis

There are several statistical features of the chemical space network that are of interest to graph. The three statistics I chose to analyze were:

◼

The number of remaining non-singleton molecules

◼

i.e. the number of molecules similar to at least one other molecule at a given similarity threshold

◼

The number of non-singleton connected components

◼

The size of the largest connected component

The first and third statistics are monotonically decreasing as the similarity threshold is increased from 0 to 1. In general, we would like for the number of non-singletons to remain relatively high over the whole range of similarity thresholds as if too many molecules are singletons, we cannot draw as many meaningful conclusions about the similarities and differences between the molecules we are analyzing.

Typically, the second statistic increases to a maximum at a similarity threshold slightly below 1 and then decreases a varying amount. However, there is a lot of variability in this statistic depending on the system of reactions and reactants chosen. For example, there is the possibility for several local maxima, and this seems to occur more often in smaller systems.

The function below takes in the association we initially created as well as a

similarityIncrement

that controls how fine each of the plots are and returns a list of three

ListLinePlot

objects labeled with the appropriate information.

In[64]:=

ChemicalSpaceNetworkAnalysis[chemicalSpaceNetworkIn_,similarityIncrement_]:=Module[{chemicalSpaceNetwork=chemicalSpaceNetworkIn,graph,iVals={},numMoleculesRemaining={},numConnectedComponents={},largestConnectedComponentSize={},connectedComponents},Do[chemicalSpaceNetwork=Select[chemicalSpaceNetwork,#>=i&];graph=chemicalSpaceNetwork//Keys//Graph[#&/@Union@@List@@@#,#]&;iVals=iVals~Join~{i};numMoleculesRemaining=numMoleculesRemaining~Join~{graph//VertexList//Length};connectedComponents=ConnectedComponents@graph;numConnectedComponents=numConnectedComponents~Join~{connectedComponents//Length};largestConnectedComponentSize=largestConnectedComponentSize~Join~{Max[{0,Length/@connectedComponents}//Flatten]};,{i,0,1,similarityIncrement}];{ListLinePlot[MapThread[List,{iVals,numMoleculesRemaining}],PlotLabel->"Non-Singletons Remaining",PlotRange->Full],ListLinePlot[MapThread[List,{iVals,numConnectedComponents}],PlotLabel->"Number of Connected Components",PlotRange->Full],ListLinePlot[MapThread[List,{iVals,largestConnectedComponentSize}],PlotLabel->"Size of Largest Connected Component",PlotRange->Full]}]

Here are the plots generated for the same monoterpene system mentioned in the “Multiway Chemical Reaction Network” section, this time carried out four layers deep.

We can see that at a similarity threshold of 1, i.e. molecules are only similar if they are topologically identical (according to their topological fingerprints), approximately one third of the molecules are singletons. The majority of the chemical reaction systems I investigated had between 10-30% of their molecules as singletons at a similarity threshold of 1.

Here we can see that the number of connected components does not drop off as the similarity threshold approaches 1. We can also see that the chemical space network remains connected until we impose a similarity threshold slightly over 0.5 and then steadily grows to around 40 clusters.

This graph is very similar to a step function as we can see several thresholds where the largest connected component fractionates into two or more clusters. Another behavior that can be seen with the size of the largest connected component is a steep dropoff at a high similarity threshold. based on my investigation, this seems to be more common in polymer-building systems.

Graphs

We would also like to display the chemical space network at various similarity thresholds to aid in our analysis of these systems. The function below creates a list of chemical space network graphs, also including the ability to click on any of the nodes to get the SMILES string of its associated molecule.

In[65]:=

ChemicalSpaceNetworkGraphs[moleculeSMILES_,chemicalSpaceNetworkEdges_,{minSimilarityThreshold_,maxSimilarityThreshold_,thresholdIncrement_}]:=Module[{graphs={},edges=chemicalSpaceNetworkEdges},Do[edges=Select[edges,#>=i&];graphs=graphs~Join~{Graph[Button[#,CopyToClipboard["\""<>#<>"\""]]&/@Normal[moleculeSMILES],Keys@edges,GraphLayout->"GravityEmbedding",PlotLabel->i]};,{i,minSimilarityThreshold,maxSimilarityThreshold,thresholdIncrement}];graphs]

Pictured below are the chemical space networks for the four layer deep monoterpene system described earlier at six different similarity thresholds.

Combining the Two

Colored Clusters on the Multiway Graph

One way we can combine the multiway chemical reaction network and chemical space network is to display the clusters (connected components) from the chemical space network on the multiway graph using a unique color for each cluster.

Although the only difference between this graph and the original multiway graph is the coloring, we create it separately as it is much easier to modify this way. The

ClusterColoredMuktiwayChemicalReactionNetwork

function also allows us to only color clusters of a given size or larger. Any clusters too small to be colored are given a gray background

Shown below is the three layer deep monoterpene system analyzed at a similarity threshold of 0.9.

We can see that, for the most part, molecules of the same cluster (those with the same colors) are genealogically close. There are some notable exceptions, such as light green, which has only five molecules in its cluster, but those molecules are scattered across the graph and not very genealogically close to each other.

Genealogical Distance on Chemical Space Networks

One other way we can combine the information from the multiway graph and the chemical space network is using the distance between the nodes in the multiway graph, which will be referred to as the genealogical distance. For each cluster in the chemical, one node is chosen as the central node, and the genealogical distances from that central node to all other nodes in the cluster are calculated. Each node is then given a color corresponding to its genealogical distance from its central node, with central nodes all colored red. The redder a node is, the closer it is genealogically to its cluster’s central node.

Shown below is the 4 layer deep monoterpene system at a similarity threshold of 0.9 with the genealogical distance indicated by the colors.

As we can see the image, most of the nodes within smaller clusters are genealogically very close to each other with only a few of the smaller clusters having colors besides red and orange. For the larger cluster, we see more yellow and green, which shows that the molecules within the clusters are more genealogically spread out. It is possible that these clusters would appear more red if we allowed for genealogical distance to go down trough the multiway graph, i.e. calculate the distance through “descendant” molecules, but for most analytical purposes, the definition of genealogical distance used in this analysis seems more applicable.

Branchial Analysis

Because we are investigating the genealogical distance between molecules in these reaction networks, branchial graphs are of much interest. The

BranchialGraphs

resource function provides a list of graphs where two nodes are connected if they share a “parent” node, i.e. a node such that there exists a directed edge from the “parent” node to the two nodes of interest (called “child” nodes).

Shown below are the branchial graphs for a two layer deep Diels-Alder reaction network discussed in the “Pattern Reaction Paradigms” section.

Even going only two layers deep, we can see the beginnings of a cluster separation in the genealogy of the molecules. Each side of the branchial graph on the right has many more internal connections than connections with the other side.

Unfortunately, the

BranchialGraphs

function only works on acyclic graphs, so any multiway chemical reaction network that has reversible or cyclic reactions cannot be analyzed using branchial graphs.

Full Analysis

Example: Geranyl Cation/Monoterpene Biosynthesis Reactions (3 steps)

Statistical analysis graphs:

Chemical space network at a similarity threshold of 0.9:

Multiway chemical reaction network colored using the clusters from above:

Chemical space network with nodes colored according to genealogical distance from cluster centers:

Principal component analysis on genealogical distance matrix from multiway graph:

Code

This section is a potential workflow to generate the analysis for a multiway chemical reaction system of your choice.

The variables in this input cell are for:

◼

fullAnalysisMolecules

: list of molecules to start your multiway chemical reaction system

◼

fullAnalysisReactions

: list of reactions for your system

◼

fullAnalysisNumSteps

: the number of reaction steps your system will undergo

◼

fullAnalysisSMILES

: the number of reaction steps your system will undergo

In[261]:=

fullAnalysisMolecules={geranylCation};fullAnalysisReactions=monoTerpeneBiosynthesisReactions;fullAnalysisNumSteps=4;fullAnalysisSMILES=GetMultiwayChemicalReactionNetworkSMILES[MultiwayChemicalReactionNetwork[fullAnalysisReactions,fullAnalysisMolecules,fullAnalysisNumSteps]];

The variable in this input cell is for:

◼

fullAnalysisChemicalSpaceNetwork

: association of molecule pairs to their similarity scores

In[607]:=

fullAnalysisChemicalSpaceNetwork=ChemicalSpaceNetwork[fullAnalysisSMILES];ChemicalSpaceNetworkAnalysis[fullAnalysisChemicalSpaceNetwork,0.01]

The variable in this input cell is for:

◼

fullAnalysisChemicalSpaceNetworkGraphs

: chemical space network graphs at specified similarity thresholds

In[609]:=

fullAnalysisChemicalSpaceNetworkGraphs=ChemicalSpaceNetworkGraphs[fullAnalysisSMILES,fullAnalysisChemicalSpaceNetwork,{0.5,1,0.05}]

The variables in this input cell are for:

◼

fullAnalysisMinimumClusterSize

: smallest cluster size that will be given a color on the multiway graph

◼

fullAnalysisChemicalSpaceNetworkGraph

: chemical space network chosen from fullAnalysisChemicalSpaceNetworkGraphs to take clusters from

◼

fullAnalysisMultiwayGraph

: multiway graph colored according to clusters of fullAnalysisChemicalSpaceNetworkGraph

In[610]:=

fullAnalysisMinimumClusterSize=3;fullAnalysisChemicalSpaceNetworkGraph=fullAnalysisChemicalSpaceNetworkGraphs[[9]];fullAnalysisMultiwayGraph=ClusterColoredMultiwayChemicalReactionNetwork[fullAnalysisReactions,fullAnalysisMolecules,fullAnalysisNumSteps,fullAnalysisChemicalSpaceNetworkGraph,fullAnalysisMinimumClusterSize,"EventVertices"->False];

This input cell returns a chemical space network with the nodes colored according to their genealogical distance from the central node of their respective connected components. Note the “0.9” must match the similarity threshold of

fullAnalysisChemicalSpaceNetworkGraph

In[614]:=

ChemicalSpaceNetworkAndGenealogicalDistanceGraph[fullAnalysisSMILES,fullAnalysisChemicalSpaceNetwork,fullAnalysisMultiwayGraph,0.9]

This input cell returns a scatterplot of the principal components analysis on the genealogical distance matrix generated from the multiway chemical reaction network.

In[613]:=

ColoredDimensionReducedClusters[fullAnalysisSMILES,fullAnalysisChemicalSpaceNetworkGraph,fullAnalysisMinimumClusterSize,2]

Future Work

Branchial Analysis

As mentioned previously, branchial analysis, as it currently exists in Wolfram Language, is incompatible with most reaction networks due to the existence of directed cycles. It would be an interesting endeavor to create a version of BranchialGraph that is able to work with cyclic graphs. I suspect this would have applications far beyond the analysis done for this project, so I look forward to seeing someone tackle this problem.

Multiway Chemical Reaction Network

Number of Reactants

As of now, MultiwayChemicalReactionNetwork can only handle reaction sets that have the same number of reactants for all reactions, and even then, the function has to be modified manually to accommodate different numbers of reactants. A logical next step would be to rewrite the rules given to the token event graph function in such a way as to allow for a set of reactions with different numbers of reactants.

Pruning Function

Many reaction systems of interest quickly explode with the number of possible molecules, slowing any analysis beyond the point of usefulness. To remedy this, a function that removes molecules at each step in the token event graph based on some input criteria would allow for larger depth investigations of some of these systems.

There could also be a Random option for the pruning function that only allows the multiway graph to keep a certain number of molecules at each time step, choosing these at random.

Probabilistic Multiway Computation

One of the downsides of using multiway computation to model chemical reaction networks is that the edges in the generated graphs all have the same weights, which can be interpreted as the reactions being equiprobable. This is not very realistic, and it would be of great use to incorporate a probabilistic aspect into the multiway computation. For example, if each reaction was assigned a probability of occurring, each molecule could be paired with a list of probabilities associated with how likely it is to reach that molecule from a given state after a given number of steps. The greatest challenge with this would likely be visualization, so a great deal of creativity would be needed to incorporate this into the Wolfram Language.

Chemical Space Network

Other Similarity Metrics

The analysis with the current implementation can only be done using the topological fingerprints and the JaccardDissimilarity function for similarity comparisons. To expand the analytical abilities of this method, a more general form for the ChemicalSpaceNetwork function allowing for different metrics to be used should be implemented

GUI

A WSS ’23 project by Joshua Khorsandi explored chemical space networks in the context of comparing bio activity to chemical structure. As part of his project, he developed a GUI that allowed for individual cluster investigation and molecular inspection among other capabilities (4). His work should be modelled after and expanded upon to incorporate the multiway chemical reaction network and analyses described in this post.

Concluding Remarks

The versatility of multiway computation has great potential in chemical analysis, especially when paired with more traditional analysis methods, such as chemical space networks. From the results presented, it is clear that multiway chemical reaction networks and chemical space networks provide very similar, yet distinct, information about the chemical systems that they are analyzing. This makes pairing them together a very powerful tool: verifying the validity of each other while also providing more information than either method could on its own.

Acknowledgments

I would like to acknowledge Robert Nachbar, who served as my mentor throughout this project. He consistently helped me focus my project towards goals that were both significant and achievable in the time frame of the Wolfram Summer School and swiftly resolve any coding challenges I ran into during my project. I would also like to thank Stephen Wolfram and Mads Bahrami for directed me towards this project and providing their insights.

References

G. Maggiora and J. Bajorath (2014), “Chemical space networks: a powerful new paradigm for the description of chemical space,” J Comput Aided Mol Des, https://pubmed.ncbi.nlm.nih.gov/24925682/

[FSWS23] A Multiway Model of Chemical Reactions, by Ariana-Dalia Vlad, Wolfram Community, January 13, 2023, https://community.wolfram.com/groups/-/m/t/2778062

Modeling Chemical Reaction Networks with Token-Event Graphs, by Robert Nachbar, Wolfram Technology Conference, October 21, 2022, https://www.wolfram.com/broadcast/video.php?c=104&p=6&disp=list&v=3812

Chemical space networks from molecule data, by Joshua Khorsandi, Wolfram Community, July 12, 2023, https://community.wolfram.com/groups/-/m/t/2958091

Cite This Notebook

Analyzing Multiway Chemical Reaction Networks Using Chemical Space Networks
by Nicholas Frieler
https://community.wolfram.com/groups/-/m/t/3210073

Initialization


Cite this as: Nicholas Frieler, "Analyzing Multiway Chemical Reaction Networks Using Chemical Space Networks" from the Notebook Archive (2024), https://notebookarchive.org/2024-07-488otbg

Introduction

MoleculePlotOutputForm

Pattern Reaction Paradigms

Construction

Self-Interaction

Combination

Multiway Chemical Reaction Network

Chemical Space Network

Statistical Analysis

Graphs

Combining the Two

Colored Clusters on the Multiway Graph

Genealogical Distance on Chemical Space Networks

Branchial Analysis

Full Analysis

Example: Geranyl Cation/Monoterpene Biosynthesis Reactions (3 steps)

Code

Future Work

Branchial Analysis

Multiway Chemical Reaction Network

Number of Reactants

Pruning Function

Probabilistic Multiway Computation

Chemical Space Network

Other Similarity Metrics

GUI

Concluding Remarks

Acknowledgments

References

Cite This Notebook

Initialization

Initialization
