Drug-drug interaction prediction model and training method thereof
Patent Information
- Application Number
- CN202610762891.7
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-05-29
- Publication Date
- 2026-08-18
AI Technical Summary
第一,模型普遍缺乏因果推理能力
Smart Images

Figure CN122598831A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the interdisciplinary field of pharmaceutical informatics and artificial intelligence, specifically to drug-drug interaction prediction models and their training methods. Background Technology
[0002] Drug-drug interaction (DDI) refers to the phenomenon where the pharmacological effect of one drug is altered due to the presence of another drug when two or more drugs are used simultaneously or sequentially. This alteration may lead to reduced efficacy, increased toxicity, or even life-threatening serious adverse events. Statistics show that drug interactions are one of the main causes of drug withdrawal from the market and clinical medication accidents. Therefore, accurately predicting and assessing potential drug interactions in the early stages of drug development and in clinical medication decisions has significant scientific and clinical value.
[0003] In recent years, with the rapid development of deep learning technology, data-driven drug interaction prediction methods have made significant progress. Early methods mainly relied on traditional features such as drug molecule fingerprints and structural similarity, combined with machine learning models such as support vector machines and random forests for prediction. Subsequently, graph neural networks (GNNs) have been widely used in DDI prediction tasks, which can directly take drug molecule graphs as input and automatically learn the topological features related to interactions in the molecular structure through message passing mechanisms. Furthermore, multimodal fusion methods have begun to integrate multi-source heterogeneous information such as drug chemical structures, target proteins, gene pathways, and electronic health records to improve the accuracy and robustness of predictions.
[0004] However, existing technologies still have the following significant shortcomings: First, the models generally lack causal reasoning ability. Most existing methods rely on statistical correlation for prediction, failing to distinguish between "causal structures leading to interactions" and "related structures merely co-occurring with interactions." This results in insufficient generalization ability when faced with novel drug combinations not adequately covered in the training data, and makes it difficult to provide pharmacologically meaningful mechanistic explanations. Second, interpretability is severely lacking. Most high-performance DDI prediction models are "black box" models, only outputting interaction categories or probabilities, unable to explicitly identify which structural segments in the drug molecule are key drivers of interactions. This deficiency severely restricts the practical application value of the models in drug structure optimization and clinical decision support. Third, the training process lacks adaptability. Existing methods typically employ fixed training strategies; key hyperparameters such as the generation intensity of counterfactual samples and the range of feature perturbations cannot be dynamically adjusted according to the model's learning state, leading to low training efficiency or insufficient causal signal extraction. Fourth, cross-sample knowledge sharing mechanisms are lacking. The causal discovery process of each training sample is independent, making it impossible to utilize learned global causal patterns to accelerate inference and causal localization of new samples, limiting the model's performance in data-sparse scenarios.
[0005] In summary, there is an urgent need for a drug interaction prediction method and system that can achieve both high-precision prediction and provide molecular-level interpretable causal reasoning, in order to make up for the shortcomings of existing technologies in causal discovery, interpretability, adaptive training, and cross-sample knowledge reuse. Summary of the Invention
[0006] To address the problems of existing technologies, this invention provides a drug-drug interaction prediction model and its training method.
[0007] To solve the above-mentioned technical problems, the present invention is achieved through the following technical solution: Firstly, a training method for a drug-drug interaction prediction model, comprising: Acquire sample data for multiple drug pairs, each drug pair comprising a first drug molecule map, a second drug molecule map, and interaction labels; An initial model is constructed, which includes a graph neural network encoder, a causal mask generation network, an interaction prediction network, a dynamic perturbation strength controller, and a global causal prior knowledge base; The node features of the first drug molecule graph and the second drug molecule graph are extracted using the graph neural network encoder. The causal mask generation network generates a current causal mask vector for each drug molecule graph based on the node features. The causal mask vector is used to indicate the probability that each node or edge in the drug molecule graph belongs to a causal substructure. The drug molecule diagram is input into the interaction prediction network to obtain the first prediction result; The first prediction result is used as a feedback signal, fused with the node features, and then re-input into the causal mask generation network to generate an updated causal mask vector, thereby realizing the dynamic guidance of the prediction result for the discovery of causal substructures. Based on the updated causal mask vector, the structure of the drug molecule graph is perturbed to generate a counterfactual molecule graph; The counterfactual molecular graph is input into the interaction prediction network to obtain a second prediction result; The prediction loss is calculated based on the first prediction result and the label, and the causal calibration loss is calculated based on the difference between the first prediction result and the second prediction result. The dynamic perturbation intensity controller automatically adjusts the degree of structural perturbation of the drug molecule map in the next training round based on the magnitude of the causal calibration loss, thereby achieving adaptive linkage between perturbation intensity and model prediction capability. With the goal of minimizing the joint loss, the parameters of the initial model are optimized, and the joint loss includes at least the prediction loss and the causal calibration loss, until the convergence condition is met, thus obtaining a trained drug interaction prediction model.
[0008] In one specific implementation of the first aspect, the step of using the first prediction result as a feedback signal, fusing it with the node features, and then re-inputting it into the causal mask generation network includes: The first prediction result is encoded as a prediction state vector, which contains features of the interaction category probability distribution or the causal calibration loss. The predicted state vector is concatenated with or the node features are fused with attention weighting to obtain node features with enhanced feedback. The updated causal mask vector is generated based on the node features enhanced by the feedback.
[0009] In one specific embodiment of the first aspect, the automatic adjustment of the disturbance level by the dynamic disturbance intensity controller includes: A threshold is set, and when the causal calibration loss is lower than a first threshold, the perturbation level of the next round of training is increased, including increasing the proportion of masked nodes, expanding from masked single atoms to masked functional groups, or increasing the variance of the substitution noise. When the causal calibration loss is higher than the second threshold, the perturbation level of the next round of training is reduced to avoid excessive loss of causal substructures leading to training instability.
[0010] In one specific embodiment of the first aspect, a substructure antagonistic substitution step across the drug pair is further included: From samples of other drug pairs, high-confidence causal substructures identified by corresponding causal mask vectors are extracted and associated with known pharmacological function identifiers and stored in the global causal prior knowledge base; When generating the counterfactual molecular graph in this round, the causal substructures masked in the drug molecular graph are replaced with antagonistic substructures that are topologically similar but have different pharmacological function identifiers retrieved from the global causal prior knowledge base, so as to generate the antagonistic counterfactual molecular graph. The adversarial counterfactual molecular graph is input into the interaction prediction network to obtain a third prediction result; Based on the difference between the third prediction result and the label, the adversarial calibration loss is calculated and added to the joint loss.
[0011] In one specific implementation of the first aspect, the global causal prior knowledge base is continuously updated during the training process, storing the high-confidence causal substructures of each drug pair confirmed after training and their corresponding causal effect strengths. When generating causal mask vectors for new drug pairs, the causal mask generation network retrieves topologically similar known causal substructures from the global causal prior knowledge base and integrates their prior features into the discovery process of the current causal substructure through an attention mechanism, forming a synergistic enhancement loop between global knowledge and local discovery.
[0012] In one specific implementation of the first aspect, the causal mask vector is obtained by discretizing the causal importance scores of nodes or edges using the Gumbel-Softmax reparameterization technique, so that the mask value approaches 0 or 1.
[0013] In one specific implementation of the first aspect, the causal calibration loss is constructed based on maximizing the distributional difference between the first prediction result and the second prediction result, using maximum mean difference (MMD) or Kullback-Leibler divergence (KL divergence), such that when the causal substructure is removed, the model's prediction result is substantially reversed.
[0014] In one specific implementation of the first aspect, when the trained drug interaction prediction model is performing inference, the causal mask vector is output synchronously to visually highlight the causal substructures driving the interaction in each drug molecule.
[0015] Secondly, a method for predicting drug interactions includes: Obtain the first and second drug molecule maps of the drug pair to be predicted; Input the first drug molecule map and the second drug molecule map into the drug interaction prediction model trained using the training method of the drug-drug interaction prediction model; The model outputs the interaction prediction results of the drug pair, as well as the corresponding causal substructure identifiers.
[0016] Thirdly, a training system for a drug interaction prediction model based on multi-closed-loop feedback includes: The data acquisition module is used to acquire sample data from multiple drug pairs; The encoding module is used to extract node features of the drug molecule graph using a graph neural network encoder; The causal mask generation module is used to generate a causal mask vector based on the node features, and receive the first prediction result fed back by the interaction prediction module, and generate an updated causal mask vector after fusing it with the node features. The counterfact generation module is used to perturb the original drug molecule graph based on the updated causal mask vector to generate a counterfact molecular graph, and can retrieve antagonistic substructures from the global causal prior knowledge base for replacement. The prediction module is used to predict the interaction between the original drug molecule map and the counterfactual molecule map, and outputs the first prediction result and the second prediction result. A dynamic disturbance intensity control module is used to automatically adjust the degree of structural disturbance performed by the counterfactual generation module based on the magnitude of the causal calibration loss. The joint loss calculation module is used to calculate the predicted loss and the causal calibration loss; A global causal prior knowledge base is used to store high-confidence causal substructures and their effect features discovered during training, and to provide the causal mask generation module with query and fusion capabilities. The parameter optimization module is used to update the model parameters based on the joint loss.
[0017] The beneficial effects of this invention are as follows: 1. This invention upgrades drug interaction prediction from unidirectional static inference to a learning paradigm with self-driven causal reasoning capabilities by constructing a four-pronged closed loop: prediction feedback fusion, dynamic perturbation intensity control, cross-sample adversarial substructure replacement, and global causal knowledge collaboration. Specifically, the prediction feedback mechanism inversely encodes the output state of the interaction prediction network and injects it into the causal mask generation network, enabling the search for causal substructures to be dynamically optimized based on the model's current prediction confidence, significantly suppressing the interference of irrelevant structures on causal discovery. The dynamic perturbation intensity controller continuously adjusts the generation difficulty of counterfactual samples between conservative masking and functional group-level adversarial replacement based on the adaptive change of the causal calibration loss, ensuring that the model always learns within the optimal difficulty range, thereby fully exploring the deep causal relationships in molecular structures while maintaining training stability. Cross-sample adversarial substructure replacement utilizes topologically similar but pharmacologically different substructures in the global knowledge base for cross-validation, forcing the model to learn the specific mapping between structural features and pharmacological effects, effectively eliminating the misleading influence of surface correlation on causal inference; the global causal prior knowledge base continuously accumulates high-confidence causal substructures and their effect patterns, injecting cross-sample causal priors into the mask generation process of new samples in an attention-weighted manner, achieving synergistic enhancement of local causal discovery and global knowledge evolution, and significantly improving the model's generalization ability to rare drug pairs and novel substructures.
[0018] 2. The organic integration of the technical means in this invention enables it to identify concise causal substructures driving interactions while outputting high-precision interaction prediction results, fundamentally overcoming the shortcomings of traditional black-box models that cannot provide decision-making basis. The joint constraint of causal calibration loss and sparsity regularization ensures that the discovered substructures possess both real causal effects and high chemical interpretability, directly providing molecular-level mechanistic explanations for drug structure optimization and clinical combination drug risk assessment. The multi-level interconnected self-driven training architecture eliminates dependence on external intervention for parameter tuning, and the model's causal reasoning ability evolves autonomously with the training process, significantly improving training efficiency and model robustness. The continuous accumulation and reuse mechanism of the global knowledge base enables the model's causal discovery ability to continuously enhance with the increase of training data, exhibiting good scalability and continuous learning characteristics, providing reliable and transparent basic technical support for large-scale, multi-source drug interaction knowledge discovery. Attached Figure Description
[0019] Figure 1 This is a schematic diagram of the overall system architecture of the present invention.
[0020] Figure 2 This is a schematic diagram of the overall system architecture of the present invention. Detailed Implementation
[0021] The technical solutions of the present invention will be clearly and completely described below with reference to the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.
[0022] like Figures 1 to 2 The drug-drug interaction prediction model and its training method are shown.
[0023] I. This invention proposes a training method and system for a drug interaction prediction model based on multi-closed-loop feedback causal reasoning. Its core innovation lies in constructing four intertwined, dynamically linked closed loops, upgrading the traditional unidirectional "feature extraction-prediction" process into a self-driven learning system with "diagnosis-verification-adaptive evolution" capabilities.
[0024] This invention uses the prediction results of the interaction prediction network as a feedback signal, which is then injected back into the causal mask generation network. This allows the network to dynamically adjust its search strategy for causal substructures based on the current prediction state, forming a closed loop of "diagnosis-feedback-optimization". Simultaneously, this invention introduces a dynamic perturbation intensity controller, which automatically adjusts the degree of structural perturbation on the drug molecule map based on the magnitude of the causal calibration loss, ensuring that training remains at an appropriate difficulty level, forming a closed loop of "validation-difficulty adjustment". Furthermore, this invention constructs a global causal prior knowledge base. When generating counterfactual samples, the identified causal substructures can be replaced with topologically similar but pharmacologically different antagonistic substructures from other samples, forcing the model to learn the specificity of the substructures rather than their correlation, forming a closed loop of "cross-validation-knowledge identification". Furthermore, during the training process, this invention continuously stores high-confidence causal substructures into the global knowledge base and retrieves prior information from the knowledge base when processing new samples to accelerate causal discovery. At the same time, new local discoveries continuously enrich the global knowledge base, forming a positive reinforcement loop of "local discovery, global accumulation, and feedback to the local".
[0025] II. Core Module Technical Solution.
[0026] 2.1 Graph Neural Network Encoder.
[0027] Functional Description: The graph neural network encoder is used to map drug molecule graphs into low-dimensional dense vector representations while preserving node-level and full-graph-level feature information, providing basic features for subsequent causal mask generation and interaction prediction.
[0028] Input data: Drug molecule graph G=(V,E). Where V is the set of atomic nodes, each node v_i has an initial feature vector x_i; E is the set of chemical bond edges, each edge has features such as bond type.
[0029] The initial feature dimensions for each node are 78. Specifically, these include: atom type (44 dimensions, covering common atoms such as C, N, O, S, P, F, Cl, Br, and I), atomic chirality (4 dimensions), formal charge (1 dimension), hybridization state (5 dimensions: sp, sp2, sp3, sp3d, sp3d2), aromaticity (1 dimension), number of hydrogen atoms (5 dimensions, 0-4), ring properties (1 dimension), rotatable bond properties (1 dimension), and 16 extended features for compatibility with different molecular characterization formats. These feature dimensions can be adjusted as needed, but it is recommended to keep them within the above range to ensure model performance.
[0030] Edge feature dimensions: 12 dimensions. These include chemical bond types (single, double, triple, and aromatic bonds, 4 dimensions), conjugation properties (1 dimension), intracyclic properties (1 dimension), and stereochemical information (6 dimensions).
[0031] Network Structure: A Graph Isomorphism Network (GIN) is used as the backbone network. The GIN has 3 layers, with each layer having a hidden dimension of 128. The node update formula for the l-th layer is: h_i^(l)=MLP^(l)((1+ε^(l))·h_i^(l-1)+Σ_{j∈N(i)}h_j^(l-1)) Where N(i) is the set of neighboring nodes of node i; ε^(l) is the learnable parameter of the l-th layer, initially set to 0 and automatically adjusted during training; MLP^(l) is the l-th layer multilayer perceptron, consisting of two fully connected layers with a structure of "128 to 128", and the activation function is the rectified linear unit (ReLU function). The role of this MLP is to perform a nonlinear transformation between the aggregated neighborhood information and its own information to extract higher-level semantic features.
[0032] The experimental basis for parameter settings: The 128-dimensional hidden dimension represents a balance between model expressive power and computational efficiency. Ablation experiments verified that increasing the hidden dimension from 64 to 128 dimensions significantly improved model performance (Macro-F1 improved by approximately 5.2%); while increasing it from 128 to 256 dimensions resulted in a more gradual performance improvement (only about 1.1%), but the computational cost increased by approximately 1.8 times. Therefore, 128 dimensions represent the optimal cost-effectiveness choice. The GIN layer count was set to 3 layers because drug molecule maps typically have a diameter of 3 to 5 hops, and a 3-layer GIN can cover most functional group-level substructure perception ranges. Deeper layers would lead to oversmoothing issues, which would reduce the discriminative power of node features.
[0033] Output data: After global pooling, a global representation vector h_G is obtained for each drug molecule graph, with a dimension of 256. This 256-dimensional vector is concatenated by global summation pooling (128 dimensions) and global max pooling (128 dimensions) to capture the overall semantic features and local salient features of the molecule graph. Simultaneously, the embedding matrix H of all nodes is retained, with a dimension of |V|×128, where |V| is the number of atomic nodes in the molecule graph. The node embedding matrix will be used for subsequent causal mask generation.
[0034] 2.2 Causal Mask Generation Network.
[0035] Functional Description: The causal mask generation network generates a causal importance score for each node in the drug molecule graph and obtains a causal mask vector through discretization sampling. This mask vector indicates the probability that each node in the molecule graph belongs to a causal substructure and is a core component for achieving interpretability.
[0036] Input data: Initial stage: Node embedding matrix H output by the graph neural network encoder, with dimensions |V|×128.
[0037] Feedback phase: The node feature matrix H' is enhanced by feedback, with a dimension of |V|×192.
[0038] Network Structure: The causal mask generation network consists of two fully connected layers. The first layer maps the input features from 128 (or 192) dimensions to 64 dimensions, followed by a batch normalization layer, a GELU activation function, and a dropout layer with a dropout rate of 0.2. The second layer maps the 64-dimensional features to 1 dimension and outputs the raw causal importance score s_i for each node. The aforementioned dropout layer is used to prevent overfitting, and its dropout rate of 0.2 is the optimal value determined through grid search on the validation set.
[0039] Discretization Sampling Mechanism: To achieve discrete selection of causal substructures, this invention employs the Gumbel-Softmax reparameterization technique. Specifically, for the original causal importance score s_i of each node, it is first converted into a probability value π_i = sigmoid(s_i) using the Sigmoid function, and then Gumbel distribution noise is introduced and divided by the temperature parameter τ for discretization sampling. m_i=σ((log(π_i)+g_i) / τ) Where g_i is noise independently sampled from the Gumbel(0,1) distribution, τ is the temperature parameter, and σ is the Sigmoid function. This sampling method makes the mask value approach a continuous value (about 0.5) in the early stages of training, which is beneficial for gradient flow and model exploration; as the temperature decreases, the mask value gradually approaches a discrete value of 0 or 1, achieving accurate substructure selection.
[0040] Annealing strategy for temperature parameter τ: The initial value of temperature parameter τ is set to 1.0. During training, τ gradually decreases according to an exponential decay strategy: every 10 training epochs, τ is multiplied by a decay factor of 0.9. The lower limit of τ is set to 0.1. When τ is large (close to 1.0), the mask value m_i tends to a continuous value near 0.5, and the gradient flow is smooth, which is beneficial for the early exploration of the model; as τ gradually decreases, the mask value m_i tends to a discrete value of 0 or 1, achieving accurate substructure selection. This strategy effectively balances the exploration capability in the early stage of training and the decision accuracy in the later stage. In the 200 training epoch setting of this invention, the temperature decays from 1.0 to about 0.12, so that the mask maintains continuity and differentiability in the first half of training and gradually discretizes in the second half, matching the overall convergence rhythm of the model.
[0041] Output data: Causal mask vector M = {m_1, m_2, ..., m_|V|}, where m_i ranges from [0, 1]. The dimension of the mask vector is equal to the number of nodes in the molecular graph, |V|.
[0042] Masking threshold settings: During the training phase, nodes with m_i values greater than 0.5 are included in the causal substructure category for counterfactual sample generation. During the inference phase, nodes with m_i values greater than 0.7 are identified as high-confidence causal substructures for interpretability output and visualization. The reason for using different thresholds in the two phases is that the training phase requires a relatively lenient threshold to cover enough candidate structures for causal exploration; the inference phase requires a stricter threshold to ensure the credibility and conciseness of the output results.
[0043] 2.3 Prediction and Feedback Fusion Module.
[0044] Functional Description: The prediction feedback fusion module is used to implement the first layer of linkage closed loop. This module encodes the prediction results of the interaction prediction network and injects them back into the input features of the causal mask generation network, so that the mask generation process is no longer static, but can be dynamically adjusted based on the current prediction state of the model, realizing a "diagnosis-feedback-optimization" closed loop.
[0045] Input data: First prediction result P1: a vector of length C (C is the total number of interaction categories), where each element p_c represents the predicted probability value of the drug pair belonging to category c, satisfying Σp_c=1.
[0046] The node embedding matrix H has a dimension of |V|×128 and comes from the output of the graph neural network encoder.
[0047] Processing flow: First, the first prediction result P1 is encoded through a fully connected layer. This fully connected layer has an input dimension of C, an output dimension of 64, and uses ReLU as the activation function. After passing through this layer, a prediction state vector z_pred is obtained, with a dimension of 64. This vector encodes the model's current prediction state, including the distribution characteristics of the predicted class, prediction confidence information, and classification boundary information.
[0048] The second step involves copying and expanding the predicted state vector z_pred along the node dimension. Specifically, z_pred is copied |V| times (where |V| is the number of atomic nodes in the drug molecule graph), resulting in the expanded predicted state matrix Z_pred with a dimension of |V|×64. This operation ensures that each node in the molecule graph receives the same global predicted state information.
[0049] The third step involves concatenating the expanded prediction state matrix Z_pred with the original node embedding matrix H along the feature dimension to generate a feedback-enhanced node feature matrix H' with dimensions |V|×192 (composed of 128-dimensional original node features and 64-dimensional prediction state features). This matrix contains both local topological information of the molecular structure and integrates the global prediction state context.
[0050] The fourth step involves inputting the feedback-enhanced node feature matrix H' into the causal mask generation network to generate an updated causal mask vector M1. Since the input features now incorporate prediction information, the mask network can adjust the search direction and accuracy of the causal substructure based on the current prediction state of the model.
[0051] Analysis of the linkage effect: The key innovation of this module lies in realizing the reverse guidance of predictive information for causal discovery. When the predicted probability P1 is close to the classification decision boundary (e.g., the maximum predicted probability value is between 0.4 and 0.6), it indicates that the model's prediction of the sample is highly uncertain. At this time, the information contained in the feedback signal is the largest, which can effectively drive the masking network to more actively search for causal substructures to find the root cause of uncertainty. Conversely, when the predicted probability P1 is close to 0 or 1, it indicates that the model has sufficient confidence in its prediction, and the strength of the feedback signal naturally weakens. The masking network tends to be conservative, avoiding excessive perturbation of the determined results. This adaptive mechanism significantly improves the causal search efficiency of the model. Ablation experiments show that simply adding this linkage can improve the causal substructure accuracy from 0.312 to 0.401.
[0052] 2.4 Dynamic disturbance intensity controller.
[0053] Function Description: The dynamic perturbation intensity controller is used to implement the second-level linkage closed loop. Based on the magnitude of the causal calibration loss L_causal, this controller automatically adjusts the perturbation level generated by counterfactual samples, ensuring the model is always trained within the "zone of proximal development." This avoids insufficient causal signals due to excessively weak perturbations, or unstable training due to excessively strong perturbations. This module constructs a negative feedback adjustment loop of "validation-difficulty adjustment."
[0054] Input data: the causal calibration loss value L_causal of the current training batch (see Section 2.8 for the calculation method of this loss), and the exponential moving average L_avg of the historical loss values.
[0055] Exponential moving average calculation: The update formula for L_avg is L_avg_new = α·L_causal_current + (1-α)·L_avg_old, where the smoothing coefficient α is set to 0.9. This value allows the historical mean to quickly respond to changes in loss trends while maintaining a certain smoothness, avoiding frequent switching of disturbance levels caused by abnormal fluctuations in a single batch.
[0056] Disturbance Level Definition: This invention defines four disturbance levels, and the correspondence between each level and the loss range, as well as the specific disturbance strategies, are shown in the table below: ; The experimental basis for the classification: The aforementioned loss thresholds (0.1, 0.3, 0.8) were determined by observing the distribution of causal calibration loss during training. On the DrugBank dataset, the mean L_causal loss was approximately 0.85 in the early stages of training, decreasing to approximately 0.15 in the later stages. Losses above 0.8 were defined as "too high" (insufficient perturbation), 0.3 to 0.8 were defined as "normal range," and losses below 0.3 were defined as "too low" (potentially insufficient perturbation). Losses below 0.1 for multiple rounds were considered to indicate that the model had adapted to the current perturbation, requiring the introduction of adversarial enhancement.
[0057] Update mechanism: After each training epoch, the mean of the causal calibration loss L_causal for all batches within that epoch is calculated, and the exponential moving average L_avg is updated accordingly. Based on the updated L_avg, the default perturbation level for all samples in the next training epoch is determined.
[0058] Furthermore, this invention also establishes a sample-level fine-tuning mechanism: for a single sample, if its individual causal calibration loss value deviates from the current batch mean by more than 2 standard deviations, the perturbation level of that sample is adjusted individually. Specifically, if the individual loss is significantly lower than the batch mean (deviation < -2σ), it indicates that the sample is not sensitive enough to the current perturbation, and its perturbation level is increased by one level; if the individual loss is significantly higher than the batch mean (deviation > +2σ), it indicates that the sample has been over-perturbed, and its perturbation level is decreased by one level. This fine-tuning mechanism ensures that the training difficulty of an individual sample is precisely matched to its current learning state.
[0059] Analysis of the linkage effect: The core principle of this controller is to construct a negative feedback adjustment loop. The causal calibration loss L_causal essentially measures the "difference in model prediction results before and after removing the causal substructure." When the loss is too low, it indicates that the removed causal substructure has little impact on the prediction results, suggesting that the currently found substructure may not be a true causal structure, requiring increased perturbation intensity to explore more core regions. When the loss is too high, it indicates that the perturbation has been excessive, potentially leading to severe damage to the molecular graph structure, requiring a reduction in intensity to maintain the stability and rationality of training. This adaptive mechanism allows the model to autonomously maintain itself within the optimal training difficulty range without human intervention. Ablation experiments show that adding linkage two to linkage one further improves the causal substructure accuracy from 0.401 to 0.463.
[0060] 2.5 Counterfactual Sample Generation Module.
[0061] Function Description: The counterfactual sample generation module is used to perform precise structural intervention on the original drug molecule map based on the updated causal mask vector, generating a counterfactual molecular map. This module is also responsible for implementing the antagonistic substructure replacement function in Linkage 3, and is the core execution component for realizing causal effect verification and specificity learning.
[0062] Input data: original drug molecule graph G, updated causal mask vector M1, current perturbation level Lv, and global causal prior knowledge base (for adversarial replacement at level 4).
[0063] Conventional counterfactual generation strategies (perturbation levels 1 to 3): For disturbance levels 1 to 3, the following two basic disturbance methods or a combination thereof are used: (1) Zero-value feature masking: For a node v_i in the causal mask vector whose value exceeds the threshold corresponding to the current level, all 78 dimensions of its atomic feature vector x_i are set to 0, i.e., x_i' = 0 vector. This operation simulates the effect of "removing" the atom and all the chemical information it carries from the molecule, and is used to verify the causal necessity of the node.
[0064] (2) Gaussian noise replacement: For the selected node v_i, its atomic feature vector x_i is replaced with a random vector independently sampled from a multivariate Gaussian distribution with a mean of 0 and a variance of 0.1, i.e., x_i'~N(0,0.1·I), where I is a 78×78 identity matrix. This operation simulates the effect of replacing key atoms with "unknown" or "interference" states. The variance of 0.1 is set based on the fact that the normalized values of the original node features are usually between -1 and 1, and a variance of 0.1 can introduce sufficient perturbation without causing the feature values to overflow a reasonable range.
[0065] Adversarial counterfactual generation strategy (perturbation level 4): When perturbation level 4 is triggered, in addition to performing the regular perturbation of level 3, the following adversarial substructure replacement procedure is also executed: The first step is candidate substructure extraction: extracting high-confidence causal substructures from the causal mask vector M1 of the current drug. Specifically, nodes with m_i values greater than 0.7 and their edges are extracted as one or more connected subgraphs, and each connected subgraph is a candidate causal substructure.
[0066] The second step is fingerprint calculation and retrieval: For each candidate causal substructure, its extended connectivity fingerprint (ECFP4) is calculated. The fingerprint parameters are set as follows: radius = 2 (i.e., considering neighborhood information within a 2-hop range around the atom), length = 2048 bits. Using this fingerprint as the query key, a similarity search is performed in the global causal prior knowledge base.
[0067] The third step is the screening of antagonistic substructures: Entries meeting the following criteria are selected from the search results as antagonistic substructures: Tanimoto similarity greater than 0.6 (ensuring topological similarity), and the pharmacological function label associated with the entry differs from the interaction label of the current drug pair (ensuring functional difference). If multiple entries meet the criteria, the one with the strongest causal effect is selected.
[0068] The fourth step is substructure replacement: All causal substructures identified in the original molecular diagram are deleted, and the selected antagonistic substructures are replaced in their original positions using chemical bonds. During replacement, ensure that the bond type and number of connection points of the replaced substructure are consistent with the rest of the molecule. If a perfect match cannot be achieved, select the closest matching scheme and optimize the bond lengths and angles to a reasonable range according to molecular mechanics.
[0069] Output data: Counterfactual molecular graph G' (standard counterfactual version), and optional adversarial counterfactual molecular graph G_adv (output only at perturbation level 4). The graph topology and node characteristics of these counterfactual molecular graphs have been modified according to the above strategy and will be used for subsequent causal effect verification and adversarial calibration.
[0070] Loss calculation for Linkage 3: After the adversarial substructure replacement is triggered, the generated adversarial counterfactual molecular graph G_adv is input into the prediction network to obtain the third prediction result P3. The adversarial calibration loss is calculated as: L_adv = CrossEntropy(P3, y_label). This loss-penalized model, even after replacing topologically similar but pharmacologically different substructures, still predicts the behavior of the original label, forcing the model to learn the specificity of the substructure rather than the surface structural correlation. Ablation experiments show that adding Linkage 3 to Linkages 1 and 2 improves the causal substructure accuracy from 0.463 to 0.521.
[0071] 2.6 Global causal prior knowledge base.
[0072] Function Description: The Global Causal Prior Knowledge Base is a dynamic storage and retrieval system used to accumulate, manage, and reuse discovered causal knowledge during training. This module is a core component of Linkage 3 and Linkage 4, enabling cross-sample knowledge sharing and co-evolution.
[0073] Data structure: The knowledge base uses a key-value pair storage structure.
[0074] Key: Extended connection fingerprint of the substructure (ECFP4). The fingerprint parameters are set as follows: radius = 2 (considering neighborhood information within a 2-hop range around the atom), length = 2048 bits. This fingerprint is used for fast similarity retrieval, and the 2048-bit length achieves a good balance between retrieval accuracy and storage efficiency.
[0075] Value: A tuple containing the following fields: The substructure specification SMILES stands for: a unique identifier and derivation of substructures, ensuring the accurate preservation of chemical structure information.
[0076] Causal effect strength vector: The dimension is the same as the number of interaction categories C, recording the average effect strength of the substructure under various interaction types. The effect strength is defined as ΔP = |P1[this category] - P2[this category]|, that is, the change in predicted probability before and after removing the substructure.
[0077] Associated pharmacological function categories: string tags, such as "CYP2C9 inhibitor structure", "P-gp substrate structure", "aromatic hydroxylation site", etc., are used to identify the known or inferred pharmacological functions of the substructure.
[0078] Cumulative frequency of occurrence: The cumulative number of times this substructure is discovered during training, used to measure its prevalence and importance.
[0079] Confidence score: The calculation formula is frequency / (total number of entries + 1). At the initial entry, the confidence score is 1 / (current total number of entries + 1), and it gradually increases as the frequency increases.
[0080] Knowledge base update strategy: After each training cycle, the system iterates through all samples whose causal calibration loss L_causal exceeds 0.5 within that cycle. For each eligible sample, its high-confidence causal substructure (connected components consisting of nodes with m_i greater than 0.7) is extracted. For each extracted substructure: Calculate its ECFP4 fingerprint. If the fingerprint has a similarity of 0.9 or greater with the Tanimoto entry in the knowledge base, it is considered to be the same substructure and merged and updated: the effect intensity vector is weighted by frequency of occurrence, the frequency is incremented by 1, and the confidence level is recalculated.
[0081] If the similarity between the fingerprint and all entries in the knowledge base is less than 0.9, it will be directly added to the knowledge base as a new entry, with an initial frequency of 1 and an initial confidence level of 1 / (current total number of entries added to the knowledge base + 1).
[0082] Knowledge base capacity management: To avoid the knowledge base growing indefinitely and causing a decrease in retrieval efficiency, a maximum capacity of 100,000 entries is set. When the capacity limit is reached, entries are sorted from low to high according to the product of "confidence level × frequency", and the bottom 10% of entries are removed, retaining the most valuable causal knowledge.
[0083] Knowledge base query and prior fusion: During the generation and updating of the causal mask (step S6), the system uses the ECFP4 fingerprint of the candidate causal substructure identified in the initial mask M0 as the query key to perform a similarity search in the knowledge base, returning all entries with a Tanimoto similarity greater than 0.8. If the number of returned entries exceeds the preset Top-K value (K is set to 5), they are sorted in descending order by the product of "similarity × confidence", and the top 5 entries are selected.
[0084] The retrieved causal effect strength vectors of the prior substructures are mapped to prior feature vectors through a fully connected layer (input dimension C, output dimension 64). Then, the attention fusion weights α_k are calculated based on similarity and confidence. α_k = softmax(sim_k·conf_k), where k = 1, 2, ..., K Finally, the original node embedding matrix H is fused with the attention-weighted sum of prior knowledge: H_fused=H+Σ_{k=1}^{K}α_k·W_prior·effect_k W_prior is a 64×128 learnable projection matrix used to map 64-dimensional prior features to a 128-dimensional node feature space, achieving dimensional alignment and semantic fusion between prior knowledge and original features.
[0085] Analysis of the linkage effect: This module implements a knowledge enhancement loop of "local discovery → global accumulation → feedback to local". In the early stage of training, the knowledge base is sparse and has little auxiliary effect on the model. As training progresses, the knowledge base is continuously enriched, and the causal discovery of new samples can draw on the learned global causal patterns, and the convergence speed of causal substructure discovery gradually accelerates. According to the experimental results of this invention, in the later stage of training (after the 100th cycle), the model using prior fusion improves the causal discovery speed of new samples by about 23% compared with the control model without prior fusion. Ablation experiments show that adding linkage four to the first three linkages improves the causal substructure accuracy from 0.521 to 0.587 and the Macro-F1 from 0.763 to 0.782.
[0086] 2.7 Interaction Prediction Network.
[0087] Functional Description: The interaction prediction network predicts whether an interaction exists between two drug molecules and, if so, what type of interaction it is, based on their global representation vectors. The network receives the global representation from the encoder output and captures the complex feature relationships between drug pairs through bilinear interactions.
[0088] Input data: global representation vector h_G1 (256-dimensional) of the first drug molecule graph, and global representation vector h_G2 (256-dimensional) of the second drug molecule graph. Both vectors are from the global pooling output of the graph neural network encoder.
[0089] Network architecture design: Step 1, bilinear interaction layer. h_G1 and h_G2 interact through a learnable bilinear transformation: h_interact=h_G1^T·W_bilinear·h_G2 Where W_bilinear is a 256×256 learnable parameter matrix. Bilinear interaction can capture the second-order feature relationship between two drug molecule representations. Compared with simple vector concatenation or dot product, it has a stronger representational ability and can model complex patterns of "the combined effect of a certain type of feature of drug A and a certain type of feature of drug B".
[0090] The second step is the feature concatenation layer. The result of the bilinear interaction is concatenated with the original global representation vectors of the two drugs along the feature dimension: h_combined=Concat(h_G1,h_G2,h_interact) The concatenated vector h_combined has a dimension of 513 (composed of 256-dimensional h_G1, 256-dimensional h_G2, and 1-dimensional h_interact). This concatenation operation preserves the independent and interactive features of the drugs, providing a complete feature foundation for subsequent classification.
[0091] The third step is the classification layer. A two-layer fully connected network maps the concatenated features to the classification output. The first layer maps 513 dimensions to 128 dimensions, followed by a GELU activation function and a Dropout layer with a dropout rate of 0.3. The second layer maps the 128 dimensions to C dimensions (C being the total number of interacting classes), outputting the raw scores (logits) for each category.
[0092] The fourth step is probability normalization. The output C-dimensional score vector is normalized by applying the Softmax function to obtain the final probability distribution P∈R^C, where each element p_c satisfies 0≤p_c≤1 and Σp_c=1.
[0093] The dropout rate setting is based on the following: the dropout rate of the classification layer is set to 0.3, higher than the 0.2 in the causal mask generation network. This is because the prediction network directly processes the global representation, has a larger number of parameters, and is more prone to overfitting. Validation using a grid search (search ranges of 0.1, 0.2, 0.3, 0.4, and 0.5) shows that a dropout rate of 0.3 achieves optimal performance on the validation set, striking the best balance between preventing overfitting and maintaining the model's expressive power.
[0094] Output data: The normalized probability distribution P∈R^C represents the predicted probability of the drug pair belonging to each type of interaction. For binary classification tasks, C=2, and the output is the probability of the presence or absence of an interaction; for multi-class classification tasks, the output is the probability of belonging to each specific interaction category.
[0095] 2.8 Joint Loss Function.
[0096] Function Description: The joint loss function integrates four optimization objectives: prediction accuracy, causal effect significance, substructure sparsity, and adversarial robustness. It guides the update direction of model parameters through weighted summation, achieving synergistic optimization of multiple objectives.
[0097] Overall formula: L_total=λ1·L_pred+λ2·L_causal+λ3·L_adv+λ4·L_sparse Detailed definitions of each loss item: (1) Predicted loss L_pred (weight λ1=1.0): L_pred=CrossEntropy(P1,y)=-Σ_{c=1}^{C}y_c·log(p_c) Where P1 is the first prediction result, y is the one-hot encoded vector of the true label, and C is the total number of categories. The loss function uses the standard cross-entropy function to ensure the model's prediction accuracy. The weight 1.0 is the base reference value, and all other weights are adjusted and calibrated based on this.
[0098] (2) Causal calibration loss L_causal (weight λ2=0.5): L_causal=MMD(P1,P2) Wherein, MMD stands for Maximum Mean Discrepancy. This invention uses a multi-scale RBF (Radial Basis Function) kernel to calculate MMD, and the kernel function formula is as follows: k(p,q)=Σ_{σ}w_σ·exp(-||pq|| 2 / (2σ 2 )) The scale parameter σ takes values of {0.1, 1.0, 10.0}, with corresponding weights w_σ of {0.2, 0.5, 0.3}. The multi-scale kernel design enables MMD to simultaneously capture the differences in prediction distribution across fine, medium, and coarse scales. The MMD calculation formula is: MMD(P1,P2)=||E_{p~P1}[φ(p)]-E_{q~P2}[φ(q)]|| 2 _H Where φ is the latent eigenmap of the RBF kernel, and H is the Hilbert space of the regenerated kernel. In actual computation, it is estimated using the quadratic form of the kernel matrix. The goal of L_causal is to maximize the difference in the predicted distribution before and after removing a causal substructure, thereby ensuring that the substructures identified by the masked network have genuine causal effects. If the predicted distribution remains almost unchanged after removing a substructure, it indicates that the substructure is not a true causal driver.
[0099] The causal calibration loss weight λ² = 0.5 is set based on the following: Through grid search (search range 0.1, 0.3, 0.5, 0.7, 1.0), if λ² is too large (e.g., greater than 0.8), the model will overemphasize creating predictive discrepancies at the expense of predictive accuracy; if λ² is too small (e.g., less than 0.2), the discriminative power of causal substructures is insufficient, failing to effectively drive causal discovery. The value of 0.5 has been shown in experiments to achieve the optimal balance between predictive accuracy and the significance of causal effects.
[0100] (3) Adversarial calibration loss L_adv (weight λ3=0.3): L_adv=CrossEntropy(P3,y) P3 represents the prediction result of the adversarial counterfactual molecular graph. This loss is only enabled when perturbation level 4 is triggered. Its design intent is that when the model replaces causal substructures in the original molecular graph with topologically similar but pharmacologically different antagonistic substructures, the model should predict results different from the original label. This loss penalizes the model for still predicting the behavior of the original label after substructure replacement, forcing the model to learn the specificity of the substructure rather than surface correlation.
[0101] The rationale for setting the adversarial calibration loss weight λ3=0.3 is that this weight is less than the main prediction loss and the causal calibration loss because adversarial training is an auxiliary regularization method and should not dominate the training direction. Through grid search (search range 0.1, 0.2, 0.3, 0.5), the weight setting of 0.3 ensures that adversarial training enhances model robustness without compromising the main prediction and causal discovery tasks.
[0102] (4) Sparsity regularization term L_sparse (weight λ4=0.01): L_sparse=||M1||1+||M2||1=Σ_i|m_1i|+Σ_j|m_2j| Here, M1 and M2 are the causal mask vectors for the two drugs, respectively. The regularization term uses the L1 norm, encouraging mask values to approach 0, thus making the identified causal substructure as concise as possible. A concise causal substructure has better interpretability, avoiding the identification of too many irrelevant atoms as "causal," which would weaken the accuracy of the interpretation.
[0103] The sparsity regularization weight λ4 = 0.01 is set based on the fact that this is a small value, intended to impose a moderate sparsity constraint without significantly affecting the performance of the main task. Experiments (search range: 0.001, 0.005, 0.01, 0.05, 0.1) show that when λ4 = 0.01, the number of nodes in the causal substructure is reduced by approximately 35% compared to the model without sparsity constraints, while the prediction performance does not significantly decrease. An excessively large λ4 can lead to an overly sparse causal substructure, missing some true causal nodes.
[0104] Summary of weight determination method: The four weight values (λ1=1.0, λ2=0.5, λ3=0.3, λ4=0.01) were determined through a systematic grid search on the validation set. The search spaces are: λ2∈{0.1,0.3,0.5,0.7,1.0}, λ3∈{0.1,0.2,0.3,0.5}, λ4∈{0.001,0.005,0.01,0.05,0.1}. Finally, the weight combination that maximizes the Macro-F1 score on the validation set is selected as the optimal configuration for this invention.
[0105] III. Detailed description of the training process.
[0106] 3.1 Training hyperparameter configuration.
[0107] The key hyperparameters involved in the training process of this invention and their configurations are shown in the table below: ; Experimental basis for key hyperparameters: A batch size of 64 is chosen by considering both GPU memory utilization and gradient estimation stability. On an NVIDIA A100 40GB GPU, a batch size of 64 results in a training time of approximately 0.32 seconds per batch and a memory usage of approximately 18GB, which fully utilizes hardware resources while leaving reasonable margins.
[0108] A 1:1 ratio of positive to negative samples is used to mitigate class imbalance. The original DDI dataset typically has a 1:3 to 1:5 ratio of positive to negative samples, which can cause the model to be biased towards predicting the negative class if trained directly. This invention performs dynamic equalization sampling within each batch to ensure that positive and negative samples each account for 50%.
[0109] The Gumbel temperature decay strategy is designed based on the following considerations: it takes approximately 220 cycles for the temperature to decay from 1.0 to 0.1 (22 decays in total, each cycle multiplied by 0.9). In the 200-cycle training setting of this invention, the temperature is approximately 0.12 at the end of training, which ensures that the mask maintains continuity and differentiability in the first half and gradually becomes discretized in the second half, matching the overall convergence rhythm of the model.
[0110] 3.2 Single-batch training process.
[0111] In a training batch of this invention, the complete process of the model processing a batch of drug pair samples is performed according to the following steps. Suppose the current batch contains B drug pair samples (B=64), each sample consists of a first drug molecule map G1_i, a second drug molecule map G2_i, and an interaction label y_i, where i=1,2,…,B.
[0112] Step 1: Feature Extraction. Input the B first-stage drug molecule graphs and B second-stage drug molecule graphs from the batch into a graph neural network encoder. The encoder performs message passing and graph convolution operations on each drug molecule graph, outputting node embedding matrices H1_i (dimension |V1_i|×128) and H2_i (dimension |V2_i|×128), and global representation vectors h_G1_i (dimension 256) and h_G2_i (dimension 256) concatenated by global sum pooling and global max pooling. The computational complexity of this step is linearly related to the number of nodes and edges in the molecule graph, making it one of the most computationally intensive steps in the training process.
[0113] Step 2: Initial Causal Mask Generation. The node embedding matrices H1_i and H2_i are input into the causal mask generation network. The network calculates the causal importance score for each node and performs Gumbel-Softmax discretization sampling, outputting initial causal mask vectors M1_0i and M2_0i. Each element in these mask vectors takes a value between 0 and 1, initially indicating the probability that each node belongs to a causal substructure.
[0114] Step 3: First round of interaction prediction. The global representation vectors h_G1_i and h_G2_i are input into the interaction prediction network. After processing through bilinear interaction, feature concatenation, and a fully connected classification layer, the first prediction result P1_i is obtained, which is a C-dimensional probability vector representing the predicted probability distribution of the drug pair belonging to various interaction types.
[0115] Step Four: Implementing Linkage One – Prediction Feedback Fusion. First, the first prediction result P1_i is mapped to a prediction state vector z_pred_i through a fully connected layer (C-dimensional input, 64-dimensional output, ReLU activation). This vector encodes the model's current prediction state. Second, z_pred_i is copied and expanded in the node dimension to the same dimension as H1_i and H2_i. Then, the expanded prediction state matrix is concatenated with the original node embedding matrices H1_i and H2_i in the feature dimension to generate feedback-enhanced node feature matrices H1'_i (dimension |V1_i|×192) and H2'_i (dimension |V2_i|×192).
[0116] Step 5: Update the causal mask generation. The feedback-enhanced node feature matrices H1'_i and H2'_i generated in Step 4 are input again into the causal mask generation network. Since the input features now incorporate prediction information, the causal mask vectors M1_i and M2_i generated by the network are more accurate and can specifically focus on the substructures that truly drive the interactions. This step completes the "diagnosis-feedback-optimization" closed loop of Linkage 1.
[0117] Step Six: Implement Linkage Four – Global Knowledge Base Query and Prior Fusion. Using nodes with values greater than 0.5 in the update masks M1_i and M2_i, extract candidate causal substructures (taking connected components). Calculate the ECFP4 fingerprint of each candidate substructure and perform similarity retrieval in the global causal prior knowledge base to obtain relevant prior knowledge features. This prior information is fused with node features within the causal mask generation network through an attention mechanism to assist in subsequent causal substructure refinement and counterfactual sample generation.
[0118] Step 7: Determine the disturbance level. The dynamic disturbance intensity controller reads the currently maintained causal calibration loss exponential moving average L_avg and determines the default disturbance level for this batch of samples according to the level determination rules described in Section 2.4. For anomalous samples where the causal calibration loss deviates from the batch mean by more than 2 standard deviations, the disturbance level is fine-tuned individually.
[0119] Step 8: Generate counterfactual samples. The counterfactual sample generation module performs structural intervention on the original drug molecule diagrams G1_i and G2_i based on the perturbation level determined in Step 7 and the updated causal mask vectors M1_i and M2_i output in Step 5.
[0120] For perturbation levels 1 to 3, standard counterfactual molecular graphs G1_cf_i and G2_cf_i are generated by using zero-value masking and / or Gaussian noise replacement.
[0121] For perturbation level 4, in addition to the standard perturbation, an adversarial substructure replacement is performed: high-confidence causal substructures are extracted from the original molecular graph, topologically similar but pharmacologically different adversarial substructures are retrieved from the global knowledge base, and the substructure replacement operation is performed to generate adversarial counterfactual molecular graphs G1_adv_i and G2_adv_i.
[0122] Step Nine: Second and Third Rounds of Prediction. The standard counterfactual molecular graphs G1_cf_i and G2_cf_i are input again into the graph neural network encoder and interaction prediction network to obtain the second prediction result P2_i.
[0123] If this batch triggers adversarial enhancement (level 4), then the adversarial counterfactual graphs G1_adv_i and G2_adv_i are also input into the prediction network to obtain the third prediction result P3_i.
[0124] Step 10: Joint Loss Calculation. Calculate the loss value for each item according to the joint loss function defined in Section 2.8: Based on the first prediction result P1_i and the true label y_i, the prediction loss L_pred is calculated using the standard cross-entropy function.
[0125] Based on the first prediction result P1_i and the second prediction result P2_i, the causal calibration loss L_causal is calculated using the maximum mean difference of the multi-scale RBF kernel. This loss aims to maximize the difference in the predicted distribution before and after removing the causal substructure.
[0126] Based on the causal mask vectors M1_i and M2_i, the sum of the absolute values of all mask elements is calculated to obtain the sparsity regularization term L_sparse.
[0127] If a third prediction result P3_i is generated in this batch, the adversarial calibration loss L_adv is calculated based on P3_i and the true label y_i.
[0128] The total loss L_total for this batch is obtained by weighting and summing the losses according to the preset weights (λ1=1.0, λ2=0.5, λ4=0.01, and optional λ3=0.3).
[0129] Step 11: Backpropagation and Parameter Update. Calculate the gradient of the total loss L_total with respect to all learnable parameters in the model. Using the AdamW optimizer, update the model parameters based on the gradient, the current learning rate, and the weight decay coefficients. The update scope covers all learnable parameters of the graph neural network encoder, causal mask generation network, prediction feedback fusion fully connected layer, interaction prediction network, and prior knowledge fusion projection matrix.
[0130] Step 12: Implement state updates for Linkage 2 and Linkage 4. Post-processing steps at the end of this batch of training: (1) Update the global causal prior knowledge base: Traverse the samples in this batch with a causal calibration loss L_causal_i greater than 0.5, extract their high-confidence causal substructures (connected components with m_i greater than 0.7), and store or update the substructure information to the global knowledge base according to the update strategy described in Section 2.6.
[0131] (2) Update the dynamic disturbance intensity controller: Feed back the mean of the causal calibration loss of this batch to the controller and update the exponential moving average L_avg for the next batch to determine the disturbance level.
[0132] Step Thirteen: Iterate through the loop. Repeat steps one through twelve, processing all data in the entire training set in batches to complete one training cycle. Evaluate the current model performance on the validation set. If the validation set performance does not improve for 30 consecutive training cycles, trigger the early stopping mechanism, terminate training, and save the model parameters with the best performance on the validation set.
[0133] 3.3 Reasoning process.
[0134] After model training is complete, the workflow during the inference (testing or real-world application) phase is as follows: Step 1: Obtain Input. Obtain the chemical structure information of the drug pair to be predicted (drug A and drug B), usually provided in the form of SMILES strings. Use cheminformatics tools such as RDKit to convert the SMILES into a molecular graph representation, constructing a graph data structure containing atomic and edge features.
[0135] Step 2: Feature Extraction and Prediction. Input the molecular graphs of drug A and drug B into the trained model. The model executes the following sequentially: a graph neural network encoder extracts node features and global representations; a causal mask generation network (with an inference threshold set to 0.7) generates causal mask vectors; and an interaction prediction network outputs the predicted probability distribution. Counterfactual generation and loss calculation are not performed during the inference phase; only forward propagation is performed.
[0136] Step 3: Output Results. The model outputs two pieces of information: Interaction prediction results: include the predicted interaction categories and their probability values. For multi-class classification tasks, the top K most likely categories can be output in descending order of probability.
[0137] Interpretability Labeling: The causal structures driving the interaction are identified on each drug molecule diagram in the form of highlights or heatmaps (corresponding to atomic nodes with a mask value m_i > 0.7 and their connected regions). This labeling information allows medicinal chemists or clinical pharmacists to intuitively understand the basis of the model's predictions, thereby increasing confidence in the prediction results and assisting in subsequent drug structure optimization decisions.
[0138] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.
Claims
1. A training method for a drug-drug interaction prediction model, characterized in that, include: Acquire sample data for multiple drug pairs, each drug pair comprising a first drug molecule map, a second drug molecule map, and interaction labels; An initial model is constructed, which includes a graph neural network encoder, a causal mask generation network, an interaction prediction network, a dynamic perturbation strength controller, and a global causal prior knowledge base; The node features of the first drug molecule graph and the second drug molecule graph are extracted using the graph neural network encoder. The causal mask generation network generates a current causal mask vector for each drug molecule graph based on the node features. The causal mask vector is used to indicate the probability that each node or edge in the drug molecule graph belongs to a causal substructure. The drug molecule diagram is input into the interaction prediction network to obtain the first prediction result; The first prediction result is used as a feedback signal, fused with the node features, and then re-input into the causal mask generation network to generate an updated causal mask vector, thereby realizing the dynamic guidance of the prediction result for the discovery of causal substructures. Based on the updated causal mask vector, the structure of the drug molecule graph is perturbed to generate a counterfactual molecule graph; The counterfactual molecular graph is input into the interaction prediction network to obtain a second prediction result; The prediction loss is calculated based on the first prediction result and the label, and the causal calibration loss is calculated based on the difference between the first prediction result and the second prediction result. The dynamic perturbation intensity controller automatically adjusts the degree of structural perturbation of the drug molecule map in the next training round based on the magnitude of the causal calibration loss, thereby achieving adaptive linkage between perturbation intensity and model prediction capability. With the goal of minimizing the joint loss, the parameters of the initial model are optimized, and the joint loss includes at least the prediction loss and the causal calibration loss, until the convergence condition is met, thus obtaining a trained drug interaction prediction model.
2. The method according to claim 1, characterized in that, The step of using the first prediction result as a feedback signal, fusing it with the node features, and then re-inputting it into the causal mask generation network includes: The first prediction result is encoded as a prediction state vector, which contains features of the interaction category probability distribution or the causal calibration loss. The predicted state vector is concatenated with or the node features are fused with attention weighting to obtain node features with enhanced feedback. The updated causal mask vector is generated based on the node features enhanced by the feedback.
3. The method according to claim 1, characterized in that, The automatic adjustment of the disturbance level through the dynamic disturbance intensity controller includes: A threshold is set, and when the causal calibration loss is lower than a first threshold, the perturbation level of the next round of training is increased, including increasing the proportion of masked nodes, expanding from masked single atoms to masked functional groups, or increasing the variance of the substitution noise. When the causal calibration loss is higher than the second threshold, the perturbation level of the next round of training is reduced to avoid excessive loss of causal substructures leading to training instability.
4. The method according to claim 1, characterized in that, It also includes antagonistic substitution steps of substructures across drug pairs: From samples of other drug pairs, high-confidence causal substructures identified by corresponding causal mask vectors are extracted and associated with known pharmacological function identifiers and stored in the global causal prior knowledge base; When generating the counterfactual molecular graph in this round, the causal substructures masked in the drug molecular graph are replaced with antagonistic substructures that are topologically similar but have different pharmacological function identifiers retrieved from the global causal prior knowledge base, so as to generate the antagonistic counterfactual molecular graph. The adversarial counterfactual molecular graph is input into the interaction prediction network to obtain a third prediction result; Based on the difference between the third prediction result and the label, the adversarial calibration loss is calculated and added to the joint loss.
5. The method according to claim 4, characterized in that, The global causal prior knowledge base is continuously updated during the training process, storing the high-confidence causal substructures of each drug pair confirmed after training and their corresponding causal effect strengths. When generating causal mask vectors for new drug pairs, the causal mask generation network retrieves topologically similar known causal substructures from the global causal prior knowledge base and integrates their prior features into the discovery process of the current causal substructure through an attention mechanism, forming a synergistic enhancement loop between global knowledge and local discovery.
6. The method according to claim 1, characterized in that, The causal mask vector is obtained by discretizing the causal importance scores of nodes or edges using the Gumbel-Softmax reparameterization technique, so that the mask value approaches 0 or 1.
7. The method according to claim 1, characterized in that, The causal calibration loss is constructed based on maximizing the distribution difference between the first and second prediction results, using maximum mean difference (MMD) or Kullback-Leibler divergence (KL divergence) to ensure that the model's prediction results are substantially reversed when the causal substructure is removed.
8. The method according to any one of claims 1 to 7, characterized in that, When the trained drug interaction prediction model performs inference, the causal mask vector is output synchronously to visually highlight the causal substructures driving the interaction in each drug molecule.
9. A method for predicting drug interactions, characterized in that, include: Obtain the first and second drug molecule maps of the drug pair to be predicted; Input the first drug molecule map and the second drug molecule map into the drug interaction prediction model trained using the method described in any one of claims 1-8; The model outputs the interaction prediction results of the drug pair, as well as the corresponding causal substructure identifiers.
10. A training system for a drug interaction prediction model based on multi-closed-loop feedback, characterized in that, include: The data acquisition module is used to acquire sample data from multiple drug pairs; The encoding module is used to extract node features of the drug molecule graph using a graph neural network encoder; The causal mask generation module is used to generate a causal mask vector based on the node features, and receive the first prediction result fed back by the interaction prediction module, and generate an updated causal mask vector after fusing it with the node features. The counterfact generation module is used to perturb the original drug molecule graph based on the updated causal mask vector to generate a counterfact molecular graph, and can retrieve antagonistic substructures from the global causal prior knowledge base for replacement. The prediction module is used to predict the interaction between the original drug molecule map and the counterfactual molecule map, and outputs the first prediction result and the second prediction result. A dynamic disturbance intensity control module is used to automatically adjust the degree of structural disturbance performed by the counterfactual generation module based on the magnitude of the causal calibration loss. The joint loss calculation module is used to calculate the predicted loss and the causal calibration loss; A global causal prior knowledge base is used to store high-confidence causal substructures and their effect features discovered during training, and to provide the causal mask generation module with query and fusion capabilities. The parameter optimization module is used to update the model parameters based on the joint loss.