Molecular substructure attribute contribution degree calculation method and device

Abstract

The invention provides a molecular substructure attribute contribution degree calculation method and device, and relates to the field of molecular attribute prediction. The method comprises the following steps: constructing a substructure set of a target molecule; constructing a layered molecular graph comprising atomic layer nodes, substructure layer nodes and molecular layer nodes, and generating node features and relation edge features; inputting the hierarchical molecular graph into an attribute prediction model to obtain a first prediction value; constructing a covering vector for the target substructure, and performing gating covering on node contribution related to the target substructure through the covering vector; inputting the covering samples into the same attribute prediction model to obtain a second prediction value; and calculating the contribution degree of the substructure to the target attribute based on the difference value between the first predicted value and the second predicted value. Compared with the prior art, the method has the advantages that the prediction accuracy is improved through atom-substructure-molecule multi-granularity unified modeling in combination with covering comparison reasoning.

Description

Technical Field

[0001] This invention relates to the fields of molecular property prediction and machine learning, and particularly to graph representation learning and calculation of the contribution of molecular substructure properties. Specifically, it relates to a method and apparatus for calculating the contribution of molecular substructure properties. Technical Background

[0002] Molecular property prediction is a key task in modern drug development and design, aiming to predict the properties of molecules, such as solubility, toxicity, and bioactivity, based on known structural information. Traditional molecular screening processes rely heavily on biochemical experiments and expert experience, resulting in long cycles and high costs, significantly hindering the efficiency of new drug development. In recent years, with the development of deep learning technology, researchers have begun to use graph neural networks to model molecular structures, treating atoms in the molecule as nodes in a graph and chemical bonds as edges, thereby achieving end-to-end prediction of the overall molecular properties and significantly improving the efficiency of virtual screening and molecular evaluation.

[0003] While overall property prediction plays a crucial role in molecular screening, focusing solely on information at the molecule-wide level remains insufficient for practical applications. Certain key molecular properties are often dominated by specific structural units within the molecule, such as functional groups, heterocycles, and branches. Therefore, fine-grained attribution analysis of the property contributions of specific substructures within a molecule is of great significance in tasks such as molecular generation, property optimization, and model interpretability. By identifying the positive and negative impacts of a particular substructure on target properties, we can more accurately guide molecular optimization, increase the success rate of generated molecules satisfying specific property constraints, and enhance the goal-oriented nature of molecular design.

[0004] While existing research has begun to focus on the importance of substructures in predicting molecular properties and has attempted to achieve finer-grained structure attribution, current methods still have many limitations in their modeling strategies. Some methods attempt to estimate the average contribution of a substructure by averaging the overall property values ​​of all molecules containing that substructure. However, this strategy assumes that the effects of different substructures on molecular properties are independent and linearly additive, ignoring the synergistic effects and contextual dependencies between substructures under specific chemical environments. In fact, molecular properties are often the result of the synergistic effects of multiple structures, and the averaging strategy fails to accurately reflect the functional changes of substructures in different contexts, leading to attribution bias and affecting structural modification decisions during molecule generation or optimization.

[0005] In addition, other methods, based on strategies such as gradient sensitivity and input perturbation, reverse-engineer the inference process after the model has been trained to identify key structural regions. While these techniques improve model interpretability to some extent, they are highly dependent on the gradient paths and structural designs of specific models, exhibiting poor stability across different tasks or model architectures, lacking universality, and being difficult to coordinate with the model training process for optimization, thus lacking good training friendliness. Recent works have proposed dividing graph structures into chemical sub-fragments, achieving some progress in interpretability. However, these methods only apply to the structural division at the input end, failing to link multi-granularity feature representation and hidden layer information fusion at the atomic, substructure, and overall molecular levels, and thus cannot achieve unified modeling and inference of attribute contributions.

[0006] Therefore, there is an urgent need for a unified molecular modeling method that supports multi-granularity fusion, which can clearly express the hierarchical relationship between atoms, substructures and the whole molecule at the structural level, and provide explicit modeling and attribution support for substructure properties at the model mechanism level, so as to achieve accurate identification and interpretable prediction of the contribution of structural properties, and further improve the goal orientation and decision reliability of molecular generation and optimization processes. Summary of the Invention

[0007] The purpose of this invention is to provide a method and apparatus for calculating the contribution of molecular substructure properties. By constructing a hierarchical molecular graph containing atomic-level nodes, substructure-level nodes, and molecular-level nodes, multi-granular structural information and feature representations of atoms, substructures, and molecules are uniformly incorporated into the same graph structure for joint modeling. Furthermore, a masked sample is constructed for the target substructure to mask its structural information expression. Property predictions are performed on both the unmasked and masked samples under the same property prediction model. The difference in prediction results before and after masking is used to determine the contribution score of the target substructure. This solves the problem of bias and instability in contribution assessment caused by the difficulty in simultaneously characterizing the global multi-granular features of atoms, substructures, and molecules in existing technologies during contribution calculation.

[0008] To achieve the above objectives, the present invention adopts the following technical solution.

[0009] A method for calculating the contribution of molecular substructure properties includes:

[0010] S1: Construct a set of substructures for the target molecule;

[0011] Furthermore, the construction of the substructure set includes:

[0012] S11: The target molecule is divided into segments based on the BRICS (Breaking of Retrosynthetically Interesting Chemical Substructures) bond-breaking rules, and the set of atomic indices corresponding to each segment is used as the first substructure set.

[0013] The BRICS-based molecular fragmentation includes:

[0014] The RDKit tool was used to perform BRICS fragment partitioning on the structural representation of the target molecule to obtain a first substructure set consisting of multiple fragments.

[0015] S12: Extract the core skeleton structure from the molecular diagram of the target molecule according to the Bemis–Murcko skeleton extraction rule, and determine the corresponding skeleton atom set in the original molecule as the second substructure set;

[0016] The second substructure set obtained based on the Bemis–Murcko skeleton extraction rule includes:

[0017] The Murcko scaffold core structure of the target molecule is obtained, and the set of skeleton atoms corresponding to the Murcko scaffold core structure is determined in the original molecule. The chemical bonds of the original molecule are traversed to determine the connecting bonds with only one end located in the set of skeleton atoms as the breakpoint connecting bonds between the skeleton and the substituents. Based on the breakpoint connecting bonds, the target molecule is divided into atomic sets to obtain a second substructure set composed of multiple skeleton-related substructures.

[0018] S13: Obtain the third substructure set based on the functional group identification rule, including: identifying functional group substructures in the structural representation of the target molecule through predefined SMARTS matching rules, and using the atomic index set to characterize the substructures as elements in the third substructure set.

[0019] S2: Construct a hierarchical molecular graph and perform feature encoding: Construct a hierarchical molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes; generate node features for atomic type nodes, substructure type nodes, and global type nodes respectively; and perform edge feature encoding for the three types of semantically different relational edges;

[0020] Furthermore, the construction of the hierarchical molecular map also includes:

[0021] S21: Using each atom in the target molecule as an atomic layer node, and generating corresponding atomic feature vectors based on atomic element type, atomic degree, formal charge, hybridization mode, aromaticity, number of hydrogen atoms and chirality information;

[0022] S22: Each substructure in the set of substructures is determined as a substructure layer node, and the features of the substructure node are obtained by aggregating the features of the atomic layer nodes;

[0023] S23: Generate molecular layer node features: Aggregate the features of the substructure layer nodes contained therein to obtain the chemical semantic feature vector of the molecular layer node;

[0024] S24: Atom-to-atom bond edges between atomic layer nodes, representing the chemical bond connections between atoms. The edge characteristics are determined by bond type, conjugation, whether it forms a ring, and stereochemical information.

[0025] S25: Establish substructure-to-atom containment edges between substructure layer nodes and their contained atom layer nodes to characterize the containment relationship between substructures and atoms.

[0026] S26: Substructure-molecule convergence edge between molecular layer nodes and substructure layer nodes, used to represent the information convergence relationship between substructure units and the whole molecule. Substructure-atom inclusion edge and substructure-molecule convergence edge correspond to different preset edge type codes.

[0027] Furthermore, the node feature vectors of the atomic layer nodes, substructure layer nodes, and molecular layer nodes have the same dimension, and are obtained by concatenating the chemical semantic feature vector and the node type identifier vector according to the dimension. The node type identifier vector is used to distinguish whether a node belongs to the atomic layer, substructure layer, or molecular layer.

[0028] S3: Train the graph neural network based on the training data to obtain an attribute prediction model; during prediction, input the feature representation of the molecule into the attribute prediction model to obtain the prediction result of the target attribute value; denoted as the first prediction value;

[0029] Furthermore, the attribute prediction model gates the node weights based on the occlusion vector during the graph-level readout stage to suppress the contribution of the target substructure corresponding nodes to the graph-level representation.

[0030] S4: Construct a masking vector to mask the weight contribution corresponding to the target substructure;

[0031] Furthermore, the masking vector is a masking vector smask corresponding to the constructed substructure. The length of smask is consistent with the number of nodes in the hierarchical molecular graph. The value is 0 for nodes involved in the target substructure and 1 for the rest of the nodes. During the graph-level readout process of the attribute prediction model, smask masks the contribution of the nodes corresponding to the target substructure.

[0032] S5: Input the occluded sample obtained in step S4 into the attribute prediction model. In the graph-level readout stage, use the occlusion vector to gate the weight vector to obtain the second predicted value.

[0033] S6: Perform a difference operation on the first predicted value and the second predicted value to obtain a difference value, and normalize the difference value to obtain the contribution score of the substructure.

[0034] Furthermore, the contribution score of the target substructure is determined based on the difference in prediction output before and after occlusion, including:

[0035] The first predicted value obtained for the unmasked molecule is compared with the second predicted value obtained for the masked sample, and the difference between the two is calculated.

[0036] In the regression task, the first predicted value and the second predicted value are continuous numerical prediction results of the target property;

[0037] In a classification task, the first predicted value and the second predicted value are classification output results, which include any one of the model's class score, class probability, or classification decision value;

[0038] The contribution score is obtained by normalizing the difference between the prediction outputs before and after occlusion.

[0039] A molecular substructure contribution calculation device, used to execute the aforementioned molecular substructure contribution calculation method, comprising:

[0040] The substructure partitioning module is used to obtain the set of substructures of a molecule;

[0041] The layered graph construction and feature extraction module is used to construct a layered molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes, and generate node features and edge features.

[0042] The attribute prediction module is used to reason about the hierarchical molecular graph to obtain a first predicted value, and to reason about the masked molecules to obtain a second predicted value.

[0043] The masking sample construction module is used to construct masking samples for substructures in the substructure set, so that the information of the substructure in the hierarchical molecular diagram is masked.

[0044] The contribution calculation module is used to determine the contribution score of the substructure based on the difference between the first predicted value and the second predicted value.

[0045] Compared with the prior art, the present invention has the following beneficial effects:

[0046] By using hierarchical molecular diagrams to unify the modeling of atoms, substructures, and global molecules, contribution calculation can utilize both local structural information and global contextual information, thereby improving the accuracy of contribution assessment.

[0047] By setting initial features of the same dimension for substructure layer nodes and molecular layer nodes and introducing type distinguishers, multi-granularity nodes are placed in a unified representation space and have distinguishability, which facilitates graph neural networks to perform effective message passing and aggregation on heterogeneous semantic relationship edges.

[0048] By constructing masked samples and comparing the prediction differences before and after masking, the contribution of each substructure to the prediction of the target attribute can be quantified in a feasible and interpretable manner. Furthermore, by repeatedly predicting with multiple random seeds and aggregating the contributions, the impact of the randomness of model inference on the contribution results is reduced, thereby improving the stability and robustness of the results. Attached Figure Description

[0049] To facilitate understanding of the technical concept and implementation of this invention, embodiments are described below with reference to the accompanying drawings. The accompanying figures illustrate key steps in this invention, including molecular structure input, substructure set construction, hierarchical molecular graph modeling, masked sample construction, and contribution calculation based on predicted differences. It should be noted that the figures are for illustrative purposes only and are not drawn to scale; those skilled in the art can make equivalent substitutions or modifications to the structures, steps, or connections without departing from the concept of this invention.

[0050] Figure 1 This is a schematic diagram of the method flow according to an embodiment of the present invention;

[0051] Figure 2 This is a schematic diagram of the layered molecular diagram in an embodiment of the present invention;

[0052] Figure 3 This is a schematic diagram of node feature encoding in an embodiment of the present invention;

[0053] Figure 4 This is a schematic diagram illustrating the calculation of contribution scores for prediction differences before and after occlusion in an embodiment of the present invention.

[0054] Figure 5 This is a schematic diagram of the device structure provided in an embodiment of the present invention. Detailed Implementation

[0055] To enable those skilled in the art to better understand and implement the present invention, the technical solution of the present invention will be further described below in conjunction with the accompanying drawings and specific embodiments. It should be understood that the following embodiments are only used to illustrate the principles and implementation methods of the present invention and are not intended to limit the scope of protection of the present invention; without departing from the concept of the present invention, those skilled in the art can make equivalent substitutions or modifications to the step sequence, parameter settings, module implementation methods, etc., in the embodiments, all of which should fall within the scope of protection of the present invention.

[0056] A method for calculating the contribution of molecular substructure properties includes the following steps:

[0057] S1: The set of substructures that construct molecules.

[0058] In one embodiment of the invention, the structure represents the SMILES (Simplified Molecular Input Line Entry System) string of the target molecule and converts the SMILES into a molecular object.

[0059] When determining the set of substructures based on the molecular object, the molecule is divided into substructures from different chemical perspectives. The set of substructures includes the following three types:

[0060] S11 is a set of BRICS fragment substructures obtained based on the BRICS (Breaking of Retrosynthetically Interesting Chemical Substructures) bond-breaking rules.

[0061] S12 is a set of skeleton-related substructures obtained based on the Bemis–Murcko skeleton extraction rules;

[0062] S13 is a set of functional group substructures obtained based on functional group identification rules.

[0063] The multiple substructure sets mentioned above are deduplicated to obtain the final substructure set used for contribution calculation; each substructure is represented by its corresponding set of atomic indices in the target molecule, so as to align with the subsequent masking sample construction steps.

[0064] S2: Construct a hierarchical molecular graph, generate feature representations for nodes at each level, and generate edge features for the three semantic relationships.

[0065] Based on the target molecule and its substructure set, a hierarchical molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes is constructed:

[0066] Each atom of the target molecule is identified as an atomic layer node;

[0067] Each substructure in the substructure set is identified as a substructure layer node, and a connection edge is established between the substructure layer node and the atomic layer nodes it contains to characterize the substructure-atomic containment relationship.

[0068] Molecular layer nodes are constructed, and connection edges are established between the molecular layer nodes and all substructure layer nodes to characterize the convergence relationship of molecules and substructures.

[0069] Connecting edges are established between atomic layer nodes to characterize the chemical bond connections between atoms.

[0070] The above structure allows atoms, substructures, and the overall molecule to be placed in the same graph structure and to interact and transmit information.

[0071] The node features and edge features are preferably constructed in the following manner:

[0072] S21 atomic layer node features: Composed of atomic element type, atomic degree, formal charge, hybridization mode, aromaticity, hydrogen atom count and chirality information, and formed by one-hot encoding to form a vector representation of the atomic layer node.

[0073] S22 Substructure Layer Nodes: To ensure that nodes at different levels are in a unified feature space and can be distinguished, the initial features of substructure layer nodes are set to vectors with the same dimension as the features of atomic layer nodes, and the type of substructure layer nodes is explicitly identified by a type distinguisher, which is set to (1,0).

[0074] S23 Molecular layer node features: The initial features of the substructure layer nodes are set as vectors with the same dimension as the atomic layer node features, and the type of molecular layer nodes is explicitly identified by a type distinguisher, which is (1,1).

[0075] In the implementation of this invention, the initial features of both substructure layer nodes and molecular layer nodes are set as follows: the first 40 dimensions are obtained by aggregation operation of the atomic-level initial features of the corresponding atom sets, and the last 2 dimensions are type distinguishers.

[0076] Edge features include the bond type and its 3D information encoding of atom-atom bond edges; and the relation type encoding of substructure-atom inclusion edges and molecule-substructure convergence edges, so that the three semantically different relation edges can be distinguished in graph neural network message passing.

[0077] S3: Input the hierarchical molecular diagram into the attribute prediction model to obtain the first predicted value of the target attribute.

[0078] The attribute prediction model can be a graph neural network model, preferably including any one of RGCN (Relational Graph Convolutional Network), GCN (Graph Convolutional Network), or GAT (Graph Attention Network), used to perform message passing, feature aggregation, and output attribute prediction results on the hierarchical molecular graph.

[0079] The attribute prediction model learns parameters based on training data. After training, the model parameters are saved for consistent inference between uncovered and covered samples.

[0080] S4: For the target substructure in the set of substructures, construct a corresponding masking smask list based on the set of atomic indices contained in the target substructure. During the graph-level reading stage, the node weights are gated through the smask list so that the contribution of the masked nodes to the graph-level representation is reduced or becomes zero.

[0081] In one implementation of the present invention, the masking sample is achieved by constructing a node masking vector: a masking vector smask with a length consistent with the number of nodes in the hierarchical molecular graph is constructed, wherein the atomic nodes belonging to the target substructure and their corresponding substructure nodes are assigned a first label for suppression, the remaining nodes are assigned a second label for retention, and the smask is written into the hierarchical molecular graph as a node attribute.

[0082] In the graph-level readout stage of the property prediction model, the node contribution is gated based on the mask, so that the masked nodes do not participate or have reduced participation when the graph-level representation is aggregated, thereby suppressing the influence of the target substructure on the prediction results.

[0083] S5: Input the masked sample into the property prediction model to obtain the second predicted value of the target property.

[0084] The property prediction model is consistent with the attribute prediction model used in step S3, so as to ensure that the difference in prediction before and after masking only comes from the masking of the target substructure information.

[0085] S6: Based on the first predicted value obtained in step S3 and the second predicted value obtained in step S5, determine the degree of influence of the target substructure on the prediction of the target attribute, and obtain the contribution score of the target substructure; and perform the above process on each substructure in the substructure set respectively, and output the contribution result corresponding to each substructure.

[0086] The present invention provides a molecular substructure contribution calculation device for performing the aforementioned molecular substructure contribution calculation method, comprising:

[0087] The substructure partitioning module is used to obtain the set of substructures of a molecule;

[0088] The layered graph construction and feature extraction module is used to construct a layered molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes, and generate node features and edge features.

[0089] The attribute prediction module is used to reason about the hierarchical molecular graph to obtain a first predicted value, and to reason about the masked molecules to obtain a second predicted value.

[0090] The masking sample construction module is used to construct masking samples for substructures in the substructure set, so that the information of the substructure in the hierarchical molecular diagram is masked.

[0091] The contribution calculation module is used to determine the contribution score of the substructure based on the difference between the first predicted value and the second predicted value.

[0092] In an optional implementation of the present invention, in order to improve the robustness of contribution evaluation, steps S3 to S6 can be repeatedly executed on the same target substructure under different random initialization or different random sampling conditions to obtain multiple sets of contribution scores. The multiple sets of contribution scores are then statistically summarized, and the statistical results are used as the final contribution results to characterize the uncertainty of the contribution results and reduce the impact of random fluctuations.

[0093] In this embodiment, the molecular structure uses SMILES (Simplified Molecular Input LineEntry System) strings as input. RDKit is used to convert SMILES into molecular objects, and a hierarchical molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes is constructed based on this. The hierarchical molecular graph realizes multi-granularity information interaction between atoms, substructures, and molecules through node features, edge features, and relation type encoding. The substructure contribution score is calculated by constructing masked samples of the target substructure and comparing the prediction differences before and after masking.

[0094] Examples, such as Figure 1 As shown, proceed with the following steps:

[0095] S1: Obtain the set of substructures of the molecule

[0096] S11: BRICS Substructure Set

[0097] In this embodiment, to construct substructure layer nodes in a hierarchical molecular diagram, the set of substructures must first be determined from the input molecule. This embodiment provides a substructure partitioning method based on BRICS (Breakable Retrosynthetically Interesting Chemical Substructures) rules, the specific operation of which is as follows:

[0098] 11a. For the input string SMILES, call the RDKit's Chem.MolFromSmiles(smiles) function to obtain the molecular object m. If the molecular object is empty, the SMILES is considered invalid and the graphing is skipped.

[0099] 11b. Identify BRICS breakable bonds and form a set of breaking edges. Call the RDKit's BRICS.FindBRICSBonds(m) function to identify chemical bonds in the molecule that conform to the BRICS breaking rules, obtaining a list of breakable bonds, res. Each element contains two parts:

[0100] The atom index pairs (u,v) at both ends of the breakable bond are associated with the corresponding BRICS breakage type markers. In this embodiment, only the atom index pairs (u,v) are used to characterize the breakage location. Each pair of atom indices is converted into a set and summarized into a set of blocking edges: B={{u1,v1},{u2,v2},...}.

[0101] 11c. Extract all breakpoint atom sets, and summarize all atom indices that have appeared in the blocking edge set B into a breakpoint atom set:

[0102]

[0103] 11d. Construct a forbidden adjacency list for breakpoints, treating the atoms at both ends of all BRICS breakable bonds as breakpoint atoms. For each breakpoint atom p, count the set of atoms at the other end of all breakable bonds, forming a blocking mapping table N. p During subsequent traversal of the molecular adjacency structure, if the adjacent atom belongs to N... p If the connection is broken, it is considered a connection across the broken bond and is prohibited, thereby breaking the connectivity at the BRICS bond.

[0104] 11e. Under the condition of prohibiting crossing BRICS broken bonds, partition the molecular atomic graph into connected components: select an unassigned starting atom as a seed, and gradually expand the current atom set by traversing its neighboring atoms; during the expansion process, stop expanding when encountering an adjacency relationship that crosses a BRICS broken bond. After the traversal is completed, a substructure atom index set A is obtained. k Repeat the above process until all atoms have been assigned to a certain substructure set.

[0105] 11f. Obtain the set of atomic indices A for each substructure. k The output is a set of substructures, and the set of breakable BRICS bonds B is recorded as auxiliary information. This ensures that the union of all substructures covers all atoms of the original molecule, and that the substructures do not overlap with each other, thus allowing multiple substructure fragments to combine to form a complete molecular structure.

[0106] 11g. When no BRICS breakable bonds are detected, the entire atomic index of the molecule is output as a unique substructure set A0 to ensure that any input molecule can generate effective substructure layer nodes.

[0107] S12: The MurckoScaffold.GetScaffoldForMol function is called to extract the Murckoscaffold core structure of the target molecule, and the set of skeletal atom indices S in the original molecule is determined by substructure matching. All chemical bonds in the original molecule are traversed, and bonds with only one atom belonging to S are selected, denoted as the breakpoint bond set B, which characterizes the connection positions between the skeleton and substituents. The breakpoint bond set B is removed from the molecular graph, and the remaining graph is divided according to connectivity to obtain several non-overlapping subsets of atom indices as leaf substructures. The union of these leaf substructures covers all atom indices of the original molecule, thus forming a set of substructures based on the Murcko scaffold. If B is empty, the entire molecule's atom set is considered a single leaf substructure.

[0108] S13: Substructure matching of the target molecule is performed using predefined SMARTS (SMILES Arbitrary Target Specification) rules to obtain functional group substructures. A functional group parameter object FragCatParams is constructed using a functional group configuration file based on RDKit, resulting in a set of SMARTS patterns containing functional group name identifiers, denoted as the redundant carbon pattern set; simultaneously, a set of deredundancy patterns is constructed to remove matching redundant carbon atoms.

[0109] For the molecular object `mol` obtained by converting the input `SMILES`, the `GetSubstructMatches()` function is called to perform substructure matching for each type of functional group pattern, obtaining a set of tuples containing the indexes of the hit atoms of that functional group in the molecule. To eliminate redundant carbon atom matching introduced by the `SMARTS` patterns, in one implementation, the hit results of the de-redundancy pattern are used as a set of candidate functional group atoms, and only the matching results whose atom sets are subsets of the hit set of the redundant carbon pattern are retained, thus obtaining a cleaned set of functional group hit atom indices.

[0110] To avoid redundancy in the substructure set due to overlap and nesting between multiple functional group matching results, one implementation sorts all hitting functional groups in descending order of the number of atoms, and then removes smaller matches whose atom index sets are completely contained in other larger functional group sets, resulting in the final set of retained functional group substructures. For each retained functional group substructure, its list of hitting atom indices is recorded.

[0111] S2: Construct a hierarchical molecular graph and generate node and edge features, such as Figure 2 As shown, Figure 2 This is a schematic diagram of the structure of the layered molecular graph in an embodiment of the present invention, including atomic layer nodes, substructure layer nodes, molecular layer nodes, and corresponding relational edges;

[0112] Publicly available molecular property prediction datasets were selected as training and evaluation data sources. These datasets include: ESOL (Estimated SOLubility), Lipop (Lipophilicity), Mutagenicity (mutogenicity prediction), hERG (human Ether-à-go-go-Related Gene, cardiotoxicity risk prediction), and BBBP (Blood-Brain Barrier Penetration prediction).

[0113] Each sample must contain at least: the SMILES (Simplified Molecular Input Line Entry System) string, the corresponding molecular property label field, and the dataset partition field group (training, valid, test). During implementation, the task name is specified through configuration options, and the corresponding data file (ESOL is used as an example below) is read. The SMILES list is used as input, and hierarchical molecular graphs are constructed line by line and saved as graph data files for subsequent loading and use during the graph neural network training phase.

[0114] S21: Atomic Node Feature Extraction and Encoding

[0115] In this embodiment, the RDKit open-source cheminformatics tool is used to extract atomic-level features from molecules. Specifically, each SMILES string is converted into a molecular object Mol using the MolFromSmiles() function of RDKit, and the number N of atomic nodes in the molecule is obtained using GetNumAtoms(). Then, the feature vector of each atomic node in Mol is calculated one by one. The atomic node features are composed of the following concatenations:

[0116] 21a. Atomic element type characteristics: Obtain the element symbol of atom using atom.GetSymbol(), and perform one-hot encoding on the preset element set {'B','C','N','O','F','Si','P','S','Cl','As','Se','Br','Te','I','At','other'}; when an element is not in the above set, it is classified into the other category, thus obtaining a 16-dimensional one-hot vector.

[0117] 21b. Atom Degree Features: Obtain the degree of the atom using atom.GetDegree(), and perform one-hot encoding on the degree set {0,1,2,3,4,5,6} to obtain a 7-dimensional one-hot vector.

[0118] 21c. Charge and Radical Electron Number Characteristics: Obtain the formal charge of atom (atom.GetFormalCharge()) and the number of radical electrons (atom.GetNumRadicalElectrons()) respectively, and add them directly to the feature vector in numerical form to obtain 2D numerical characteristics.

[0119] 21d. Hybrid type features: Obtain the hybrid type of atom using atom.GetHybridization(), and perform one-hot encoding on the hybrid type set {SP, SP2, SP3, SP3D, SP3D2, other} to obtain a 6-dimensional one-hot vector. Hybrid states that do not belong to the above enumerated types are classified into the other category.

[0120] 21e. Aromaticity characteristics: To determine if an atom is an aromatic atom, use atom.GetIsAromatic() and add the boolean value as a 0 / 1 feature to obtain a 1-dimensional feature.

[0121] 21f. Total hydrogen number feature: Obtain the total number of hydrogen atoms of atom using atom.GetTotalNumHs(), and perform one-hot encoding on the set {0,1,2,3,4} to obtain a 5-dimensional one-hot vector.

[0122] 21g. Chirality-related features: Read the CIP chirality tag of the atom using atom.GetProp('_CIPCode') and perform one-hot encoding on the set {'R','S'} to obtain a 2D feature; at the same time, read whether the atom has latent chirality using atom.HasProp('_ChiralityPossible') and add it as a 1D Boolean feature.

[0123] The above method yields the basic feature vectors of the atomic nodes, with a total dimension of 40, which are used as input features for the atomic nodes of subsequent graph neural network models.

[0124] 21h. Type Bit Feature: To achieve unified modeling of multi-granularity nodes and avoid confusion in feature spaces of nodes at different levels, a 2-dimensional node type identifier bit is added to the feature vector of each node. This node type identifier bit is used to distinguish between atomic nodes, substructure nodes, and molecular nodes: the type bit for atomic nodes is (0,0), for substructure nodes it is (1,0), and for molecular nodes it is (1,1). In one implementation, substructure nodes and molecular nodes are virtual nodes. Aggregation operations are performed on the atomic features based on their corresponding atom sets to generate the initial feature vectors for the substructure nodes and molecular nodes. Through this design, hierarchical semantics can be explicitly injected while sharing graph neural network parameters, improving the model's ability to identify and fuse information of different granularities.

[0125] Ultimately, the input feature vector of each node is obtained by concatenating 40-dimensional atomic-level basic features with 2-dimensional node type identifiers, for a total dimension of 42.

[0126] S22: Substructure Node Features

[0127] The feature vectors of the substructure layer nodes are obtained by aggregation operations on the initial atomic features of their corresponding atom sets. This allows the substructure nodes to carry information about their internal atomic composition and local chemical environment, facilitating subsequent information interaction and representation learning at the multi-granularity level of the graph neural network, from the atom to the substructure to the molecule. Specifically, let A be the set of atom indices contained in the k-th substructure. k For any atom i A k The initial eigenvectors of its atomic nodes are denoted as Xi. R dThe X i It is constructed by splicing together element type (one-hot), atomic degree (one-hot), formal charge and free radical electron number numerical characteristics, hybridization type (one-hot), aromaticity Boolean characteristics, total hydrogen number (one-hot), and optional chiral correlation characteristics. Then the initial feature vector m of the k-th substructure node... k R d The feature vector of the substructure is calculated using mean pooling, specifically:

[0128]

[0129] The advantages of using mean pooling are as follows: For Boolean discrete features, the mean pooling can equivalently represent the proportional distribution of each category of features within the substructure, thus reflecting the compositional features of the substructure; for numerical features, such as formal charge and free radical electron number, mean pooling can obtain statistics that are independent of the substructure size, avoiding the bias introduced by the linear growth of feature amplitude with size due to differences in the number of atoms in the factor structure, thereby improving the comparability and stability of the substructure representation.

[0130] S23: Molecular Node Characteristics

[0131] The aforementioned molecular-level sub-features ensure that molecular-level nodes possess global structural and compositional information during the graph construction phase. The initial feature vector of each molecular-level node is obtained by performing mean pooling on the initial features of all atomic nodes. Through this design, molecular-level nodes encode statistical distribution information such as element type, degree, hybridization, and aromaticity at the initial moment, thereby improving the learnability and convergence stability of the global semantics.

[0132] The substructure set can be obtained by any one or more strategies among BRICS splitting, Murcko skeleton splitting, and functional group recognition splitting.

[0133] When constructing a hierarchical molecular graph, the node index uses a unified encoding method of "atomic nodes first, substructure nodes in sequence, and molecular nodes last". Let the number of molecular atomic nodes be N, and the number of substructures be M, then:

[0134] The set of atomic nodes V1=[0,1,...,N-1] has the same index as the RDKit atomic index;

[0135] The set of substructure nodes V2=[N,N+1,...,N+M-1], where the index of the kth substructure node is N+k;

[0136] The number of nodes in the molecular layer is fixed at 1, and its node index is set to g=N+M.

[0137] The above indexing rules allow for unified management of atomic nodes, substructure nodes, and molecular nodes within the same graph object, and provide consistent endpoint numbers for cross-granularity edges.

[0138] S24: Atom-to-atom bond edges and their coding formulas, such as Figure 3 As shown, Figure 3 This is a schematic diagram of node feature encoding in an embodiment of the present invention, wherein the features of substructures and molecular nodes have the same dimensions as the features of atomic nodes, and a type distinguisher is added; its first 40 dimensions can be obtained by aggregating atomic features;

[0139] To enable information exchange between nodes of different granularities, three types of directed edges with distinct semantics are constructed based on molecular chemical bond structures and hierarchical inclusion relationships. Each edge is assigned an integer edge type identifier. These edge type identifiers serve as the relation type input to the relational graph neural network, distinguishing different semantic relationships and driving different parameter transformations or message passing rules.

[0140] For the molecular object Mol obtained by converting the input SMILES using RDKit, iterate through all chemical bonds in Mol. For each bond, read its starting atom index i and ending atom index j, and add two directed edges (i, j) and (j, i) respectively in a "bidirectional undirected bond" manner to support bidirectional message passing in the graph neural network.

[0141] Regarding edge type encoding, attributes such as bond type (single bond, double bond, triple bond, aromatic bond), whether it is conjugated, whether it is cyclic, and stereochemical information are extracted from the RDKit key object, and the above attributes are mapped to integer edge type IDs through preset combination rules.

[0142] In one implementation, let the bond type be mapped to {0,1,2,3} on the set {single bond, double bond, triple bond, aromatic bond}, conjugation be mapped to {0,1}, cyclicity be mapped to {0,1}, and stereochemical information be mapped to {0,1,2,3} on the set {STEREONONE,STEREOANY,STEREOZ,STEREOE}. Then, the edge type id of the atom-to-atom bond edge is defined as:

[0143] t ij =a + 4b + 8c + 16d.

[0144] S25: Atomic-substructure contains edges

[0145] After determining the set of substructures, let A be the set of atomic indices corresponding to the k-th substructure. k And the index of the substructure node is N+k. For any i A k Construct containment edges between atomic nodes and substructure nodes to enable substructure nodes to aggregate their internal atomic information and transmit cross-granularity information back.

[0146] Specifically, in this embodiment, at least one directed edge (i, N+k) is added to each pair (i, N+k), and optionally a reverse edge (N+k, i) is added simultaneously to support bidirectional interaction. An edge type id with a preset constant (e.g., set to 66) is assigned to this type of containing relation edge to explicitly represent the hierarchical semantics of "substructure contains atoms".

[0147] S26: Substructure—Molecular aggregation edge

[0148] To enable molecular layer nodes to globally aggregate all substructure information and realize the feedback propagation of global semantics to lower layer nodes, this embodiment further constructs aggregation relationship edges between substructure nodes and molecular nodes.

[0149] Specifically, for each substructure node N+k V S Add a directed edge (N+k, g) from the substructure to the molecular node, and optionally add a reverse edge (g, N+k) at the same time; and assign another preset constant edge type id to this type of convergence relationship edge, for example, set to 67, to characterize the global organization relationship of "substructure converges to molecule".

[0150] By using the three types of relation edges and their edge type ID encoding, the semantic relationships of "chemical bond connection", "hierarchical containment" and "global convergence" can be explicitly distinguished in the same hierarchical graph object, thereby meeting the needs of hierarchical heterogeneous graph modeling and supporting subsequent relational graph neural networks to perform effective message passing and representation learning at the multi-granularity level of atoms, substructures and molecules.

[0151] S3: Input attribute prediction model to obtain the first predicted value

[0152] In implementation, the hierarchical molecular graphs and their corresponding labels are first read from the graph data file and sample partitioning file, and then divided into training, validation, and test sets according to the group field. Subsequently, a data loader is used to sample in batches, and multiple hierarchical molecular graphs within the same batch are merged into a single batch graph using a batch processing concatenation function, so that they can be input into the model in parallel for training and inference.

[0153] The attribute prediction model employs RGCN (Relational Graph Convolutional Network). During forward propagation, the model uses the node features of the batch graph as node input and the edge type as relation type input, enabling the model to apply differentiated parameter transformations and message passing rules for different semantic relations. Subsequently, the graph-level representation is input into the feedforward network to obtain the attribute prediction output. For regression tasks, the model outputs continuous numerical values ​​as y; for binary classification tasks, the model outputs logits, which are then mapped by a sigmoid function to obtain the probabilistic form of y.

[0154] During the training phase, forward prediction and backpropagation are performed once for each batch of samples to update the model parameters. For regression tasks, the mean squared error loss function MSELoss is used; for classification tasks, the binary cross-entropy loss BCEWithLogitsLoss with class imbalance weights is used, where the weights of positive samples are calculated from the positive and negative proportions of the training set. The model parameters are iteratively updated using the Adam optimizer, and an evaluation metric is calculated on the validation set. The evaluation metric for the regression task is R0. 2 MAE and RMSE are used for classification tasks, while Accuracy is used for other tasks. Early stopping is triggered when validation set metrics fail to improve over several consecutive rounds, and the parameters of the model that performed best on the validation set are rolled back.

[0155] After training, the optimal model saved by early stopping is used for forward inference on the test set, outputting the first predicted value y of the target attribute of the sample. The prediction results of the inference stage are exported and saved separately according to the data partition for use in subsequent contribution calculation and interpretation steps.

[0156] S4: Construct masking samples to obscure the target substructure information, such as Figure 4 As shown:

[0157] In this embodiment, the masked sample construction is implemented using a node-level masking method:

[0158] For the target substructure, obtain its corresponding set of atomic indices, that is, the list of atoms contained in the substructure;

[0159] Construct a masking vector smask_tensor with a length equal to the total number of nodes in the graph, and initialize it with all 1s;

[0160] Set the index position of the atomic layer node belonging to the target substructure to 0, and keep the other nodes at 1; at the same time, set the index position of the substructure layer node corresponding to the target substructure to 0, so as to further suppress the residual information caused by cross-granularity aggregation.

[0161] The smask_tensor is stored in the graph structure as graph node data, thereby suppressing the information of this substructure in the subsequent attribute prediction model;

[0162] In the graph-level readout phase of the attribute prediction model, the node representation or aggregation weights are multiplied and weighted by the smask_tensor node by node, so that nodes with smask_tensor=0 do not participate in message propagation or contribute to the final prediction.

[0163] S5: Input the masked sample to obtain the second predicted value.

[0164] The masked sample obtained in step S4, whose graph node data contains the smask_tensor, is input into the property prediction model to perform one forward inference to obtain the second predicted value y. mask The property prediction model used in step S3 employs the same network structure and model parameters to ensure that the first predicted value y and the second predicted value y are consistent. mask They are comparable.

[0165] The model first updates the node representation h using a multi-layer relational graph convolutional network:

[0166] During the graph-level readout phase, the node importance weight α = sigmoid(Linear(h)) is calculated;

[0167] The effective weights are obtained by element-wise multiplying the importance weights with the masking vector to obtain the effective weights w = α ⊙ smask_tensor;

[0168] The molecular representation h is obtained by weighted aggregation of node representations based on effective weights. G =

[0169] in:

[0170] hi: The representation vector of the i-th node, i.e., the atomic and substructure node.

[0171] wi: The weight of this node

[0172] h G The entire graph, i.e., the graph-level representation vector of the molecules.

[0173] h G Input a fully connected prediction head, output the second predicted value y mask .

[0174] S6: Determine the contribution score using differential analysis and output it.

[0175] For each substructure sample, take the mean prediction y of the complete input of the same molecule and the mean prediction y after covering the substructure. mask , with the difference yy mask This is used as the contribution of the substructure, and the contribution can be further normalized to (-1, 1) using tanh(t).

[0176] Repeat steps S4–S6 for each of the multiple substructures in the substructure set to output the contribution score for each substructure. The contribution scores can then be further sorted to identify key substructures.

Claims

1. A method for calculating the contribution of molecular substructure properties, characterized in that, Includes the following steps: S1: Construct a set of substructures for the target molecule; S2: Construct a hierarchical molecular graph and perform feature encoding: Construct a hierarchical molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes; generate node features for atomic type nodes, substructure type nodes, and global type nodes respectively; and perform edge feature encoding for the three types of semantically different relational edges; S3: Train the graph neural network based on the training data to obtain an attribute prediction model; during prediction, input the feature representation of the molecule into the attribute prediction model to obtain the prediction result of the target attribute value, which is recorded as the first prediction value; S4: Construct a masking vector to mask the weight contribution corresponding to the target substructure; S5: Input the occluded sample obtained in step S4 into the attribute prediction model. In the graph-level readout stage, use the occlusion vector to gate the weight vector to obtain the second predicted value. S6: Perform a difference operation on the first predicted value and the second predicted value to obtain a difference value, and normalize the difference value to obtain the contribution score of the substructure.

2. The method for calculating the contribution of molecular substructure properties according to claim 1, characterized in that, The construction of the substructure set includes: S11: The target molecule is divided into segments based on the BRICS bond-breaking rule, and the set of atomic indices corresponding to each segment is used as the first substructure set. S12: The target molecular skeleton structure is obtained as the second substructure set according to the Bemis–Murcko scaffold molecular core skeleton extraction rule; S13: Obtain the set of third substructures based on the functional group identification rules.

3. The method for calculating the contribution of molecular substructure properties according to claim 2, characterized in that, The segmentation of the target molecule based on the BRICS bond-breaking rule includes: The RDKit tool was used to perform BRICS fragment partitioning on the structural representation of the target molecule to obtain a first substructure set consisting of multiple fragments.

4. The method for calculating the contribution of molecular substructure properties according to claim 2, characterized in that, The obtained second substructure set includes: The Murcko scaffold core structure of the target molecule is obtained, and the set of skeleton atoms corresponding to the Murcko scaffold core structure is determined in the original molecule. The chemical bonds of the original molecule are traversed to determine the connecting bonds with only one end located in the set of skeleton atoms as the breakpoint connecting bonds between the skeleton and the substituents. Based on the breakpoint connecting bonds, the target molecule is divided into atomic sets to obtain a second substructure set composed of multiple skeleton-related substructures.

5. The method according to claim 2, characterized in that, The process of obtaining the third substructure set based on functional group identification rules includes: identifying functional group substructures in the structural representation of the target molecule using predefined SMARTS matching rules, and using an atomic index set to characterize the substructures as elements in the third substructure set.

6. The method according to claim 1, characterized in that, Step S2, which constructs a hierarchical molecular map, includes: S21: Using each atom in the target molecule as an atomic layer node, and generating corresponding atomic feature vectors based on atomic element type, atomic degree, formal charge, hybridization mode, aromaticity, number of hydrogen atoms and chirality information; S22: Each substructure in the set of substructures is determined as a substructure layer node, and the features of the substructure node are obtained by aggregating the features of the atomic layer nodes; S23: Generate molecular layer node features: Aggregate the features of the substructure layer nodes contained therein to obtain the chemical semantic feature vector of the molecular layer node; S24: Atom-to-atom bond edges between atomic layer nodes, representing the chemical bond connections between atoms. The edge characteristics are determined by bond type, conjugation, whether it forms a ring, and stereochemical information. S25: Establish substructure-to-atom containment edges between substructure layer nodes and their contained atom layer nodes to characterize the containment relationship between substructures and atoms. S26: Substructure-molecule convergence edge between molecular layer nodes and substructure layer nodes, used to represent the information convergence relationship between substructure units and the whole molecule. Substructure-atom inclusion edge and substructure-molecule convergence edge correspond to different preset edge type codes.

7. The method for calculating the contribution of molecular substructure properties according to claim 6, characterized in that, The node feature vectors of the atomic layer nodes, substructure layer nodes and molecular layer nodes have the same dimension, and are obtained by concatenating the chemical semantic feature vector and the node type identifier vector according to the dimension. The node type identifier vector is used to distinguish whether a node belongs to the atomic layer, substructure layer or molecular layer.

8. The method for calculating the contribution of molecular substructure properties according to claim 1, characterized in that, The attribute prediction model gates node weights based on the occlusion vector during the graph-level readout stage to suppress the contribution of the target substructure corresponding nodes to the graph-level representation.

9. The method for calculating the contribution of molecular substructure properties according to claim 1, characterized in that, The constructed masking vector is the masking vector smask corresponding to the constructed substructure. The length of smask is consistent with the number of nodes in the hierarchical molecular graph. The value is 0 for nodes involved in the target substructure and 1 for the rest of the nodes. During the graph-level readout process of the attribute prediction model, smask masks the contribution of the nodes corresponding to the target substructure.

10. The method for calculating the contribution of molecular substructure properties according to claim 1, characterized in that, The contribution score of the substructure is determined based on the difference in prediction output before and after occlusion, including: The first predicted value obtained for the unmasked molecule is compared with the second predicted value obtained for the masked sample, and the difference between the two is calculated. In the regression task, the first predicted value and the second predicted value are continuous numerical prediction results of the target property; In a classification task, the first predicted value and the second predicted value are classification output results, which include any one of the model's class score, class probability, or classification decision value; The contribution score is obtained by normalizing the difference between the prediction outputs before and after occlusion.

11. A device for calculating the contribution of molecular substructure, characterized in that, The apparatus is used to perform the method according to any one of claims 1-10, the apparatus comprising: The substructure partitioning module is used to obtain the set of substructures of a molecule; The layered graph construction and feature extraction module is used to construct a layered molecular graph containing atomic layer nodes, substructure layer nodes, and molecular layer nodes, and generate node features and edge features. The attribute prediction module is used to reason about the hierarchical molecular graph to obtain a first predicted value, and to reason about the masked molecules to obtain a second predicted value. The masking sample construction module is used to construct masking samples for substructures in the substructure set, so that the information of the substructure in the hierarchical molecular diagram is masked. The contribution calculation module is used to determine the contribution score of the substructure based on the difference between the first predicted value and the second predicted value.

Images