Nano-drug curative effect biomarker mining system based on multi-omics data

Abstract

The invention relates to the technical field of bioinformatics, in particular to a nano-drug curative effect biomarker mining system based on multi-omics data. The method comprises the following steps: firstly, acquiring multiple omics data of a subject, calculating an initial flux capacity corresponding to each enzymatic reaction, and constructing a directed graph model; further driving the directed graph model to synchronously execute supply-demand bidirectional rate updating at each temperature of the simulated temperature stepping sequence, locking the failure temperature and extracting a speed-limiting enzymatic reaction set limiting the global flux at the failure temperature; finally, mapping the rate-limiting enzymatic reaction set into a sensitive patient screening index combination; according to the method, multiple omics data are obtained to construct an enzymatic reaction directed graph, network supply and demand are deduced through temperature stepping, the failure temperature is locked, the rate-limiting enzymatic reaction is recognized, and the screening index combination is mapped, so that the mining accuracy of the nano-drug curative effect related biomarkers is improved; the technical problem that in the prior art, a sensitive patient screening index combination related to the nano-drug curative effect is difficult to determine is solved.

Description

Technical Field

[0001] This invention belongs to the field of bioinformatics technology, specifically relating to a system for mining biomarkers of nanomedicine efficacy based on multi-omics data. Background Technology

[0002] With the development of multi-omics detection technologies, researchers can obtain a large amount of data related to cellular metabolic state from multiple levels such as gene expression and metabolite concentration, providing a data foundation for analyzing the body's antioxidant metabolic pathways and their impact on drug efficacy. In nanomedicine research, the process of drug-induced reactive oxygen species generation and their clearance by the metabolic network within cells is closely related to the synergistic relationship between various enzymatic reactions and their upstream and downstream metabolites. Therefore, using multi-omics data mining to identify biomarkers that reflect differences in this metabolic process is of great significance for conducting patient stratification analysis.

[0003] However, existing technologies typically only screen and analyze multi-omics data through statistical correlation or single indicators, lacking a holistic characterization of the transmission relationships between multiple enzymatic reactions in antioxidant metabolic pathways, and also making it difficult to assess the limiting effect of different reaction pathways on metabolic capacity at the overall network level. In this case, it is difficult to accurately identify key enzymatic reactions that limit overall responsiveness in complex metabolic networks, thus making it difficult to reliably mine key biomarkers that can reflect differences in the efficacy of nanomedicines from multi-omics data, resulting in low reliability of existing technologies in determining screening indicators for sensitive patients. Summary of the Invention

[0004] The purpose of this invention is to provide a system for mining biomarkers of nanomedicine efficacy based on multi-omics data, which solves the technical problem that existing technologies are unable to accurately identify key enzymatic reactions that limit the overall responsiveness of the antioxidant metabolic network on the basis of multi-omics data, thus making it difficult to determine the combination of sensitive patient screening indicators related to the efficacy of nanomedicine.

[0005] The technical solution adopted in this invention is: a nanomedicine efficacy biomarker mining system based on multi-omics data, comprising: a preprocessing module, a deduction module, and an identification output module. The preprocessing module acquires multi-omics data of the subjects, calculates the initial flux capacity corresponding to each enzymatic reaction, constructs a directed graph model, and inputs it into the deduction module. At each temperature of the simulated temperature step sequence, the deduction module drives the directed graph model to synchronously perform bidirectional supply and demand rate updates, obtains the demand rate for reactive oxygen species neutralization and the global maximum flux at each temperature, locks the failure temperature, and extracts the rate-limiting enzymatic reaction set that limits the global flux at the failure temperature. The identification output module maps the rate-limiting enzymatic reaction set to a combination of sensitive patient screening indicators.

[0006] The invention is further characterized by: Furthermore, the preprocessing module acquires multi-omics data from the subjects and calculates the initial flux capacity corresponding to each enzyme-catalyzed reaction based on a preset standard reference system. A directed graph model of the enzyme-catalyzed reaction is constructed based on the initial flux capacity. Each edge of the directed graph model corresponds to an enzyme-catalyzed reaction and is bound to the initial flux capacity and preset thermal inactivation decay coefficient corresponding to the enzyme-catalyzed reaction.

[0007] Furthermore, methods for obtaining directed graph models include: The metabolic reaction relationships in the antioxidant metabolic pathway are analyzed based on a pre-built biochemical pathway knowledge base. Metabolic products in metabolic reactions are instantiated as nodes in a directed graph model, enzyme-catalyzed reactions are instantiated as directed edges connecting nodes, and the direction of the directed edges is determined by the direction from the substrate to the product. In the nodes, a source node is set to represent the source of reducing power generation; a sink node is set to represent the product pool after reactive oxygen species are neutralized; and an intermediate node is set between the source node and the sink node to represent intermediate metabolites in the metabolic pathway. Construct a directed graph model of enzymatic reactions from source nodes to sink nodes based on nodes and directed edges.

[0008] Furthermore, methods for obtaining the initial flux capacity include: Obtain the enzyme gene expression level and substrate-product concentration ratio corresponding to each enzyme-catalyzed reaction; Based on a pre-defined standard reference system, the fold change in enzyme gene expression level relative to the reference expression level and the fold change in substrate product concentration ratio relative to the reference concentration ratio are calculated. The initial flux capacity of each enzyme-catalyzed reaction is calculated based on the expression fold and change fold.

[0009] Furthermore, the deduction module generates a simulated temperature step sequence. For each temperature, it drives the directed graph model to synchronously perform bidirectional supply and demand rate updates to obtain the demand rate for reactive oxygen species neutralization. Based on the preset thermal inactivation decay coefficient, it deduces the flux limit of each edge to obtain the global maximum flux at each temperature. By comparing the global maximum flux with the demand rate, the failure temperature is locked and the set of rate-limiting enzymatic reactions that limit the global flux at the failure temperature is extracted.

[0010] Furthermore, methods for obtaining the flux cap include: Based on the difference between each temperature and the reference temperature used in the initial flux capacity calculation, and combined with the preset thermal deactivation attenuation coefficient, the attenuation factor of each side is obtained. The initial flux capacity of each side and the corresponding attenuation factor at each temperature are then fused to obtain the flux limit of each side at each temperature. Methods for obtaining demand rate include: The reactive oxygen species (ROS) generation rate of the drug at each temperature is calculated based on the Arrhenius equation. The generation rate is then combined with a preset supply-demand alignment coefficient to obtain the required rate for ROS neutralization.

[0011] Furthermore, methods for obtaining the global maximum flux include: At each temperature, a temporary network instance is constructed based on the flux upper limit of each edge of the directed graph model. The maximum flow algorithm is used to iteratively search for augmenting paths from the source node to the sink node in the temporary network instance. The bottleneck capacity of each augmenting path is accumulated until there are no more augmenting paths. The accumulated result is then used as the global maximum flux at the current temperature.

[0012] Furthermore, methods for obtaining the failure temperature include: The temperature is increased sequentially according to the simulated temperature step sequence. At each temperature, the corresponding global maximum flux and demand rate are obtained. The global maximum flux and demand rate are compared. If the global maximum flux is less than the demand rate, the current temperature is determined to be the failure temperature.

[0013] Furthermore, the identifier output module maps the rate-limiting enzyme reaction set to a combination of sensitive patient screening indicators, as follows: Methods for obtaining the rate-limiting enzyme reaction set include: In the network instance corresponding to the failure temperature, a residual network is constructed based on the calculation results of the maximum flow algorithm, and the corresponding minimum cut edge set is determined based on the maximum flow minimum cut theorem. The enzyme-catalyzed reactions corresponding to the directed edges in the minimum cut edge set are identified as rate-limiting enzyme-catalyzed reactions. By summarizing the rate-limited enzyme reactions, a set of rate-limited enzyme reactions is obtained. Methods for mapping rate-limiting enzyme reaction sets to combinations of sensitive patient screening indicators include: Analyze each edge in the rate-limiting enzyme reaction set and identify the enzyme type corresponding to each edge; By retrospectively analyzing the multi-omics data of the subjects, the expression levels of enzyme genes and the concentration ratios of substrates and products for each enzyme were extracted. The enzyme type, its corresponding enzyme gene expression level, and substrate product concentration ratio are combined to output as a set of sensitive patient screening indicators.

[0014] Furthermore, the identification output module, based on the failure temperature and preset tissue heat dissipation coefficient, combined with the photothermal conversion efficiency of the drug, reverse-engineers the preset power density reference value of the laser treatment device. The method for obtaining the preset power density reference value includes: The target temperature is determined based on the failure temperature and the preset safety redundancy threshold. Calculate the temperature difference between the target temperature and the preset reference temperature, and combine it with the preset tissue heat dissipation coefficient to obtain the heat power requirement required to maintain the target temperature. The thermal power requirement is used as the numerator, and the product of the drug's photothermal conversion efficiency and light absorption coefficient is used as the denominator. The ratio of the fractions is used as the preset power density reference value for the laser treatment equipment.

[0015] The beneficial effects of this invention are: This invention acquires multi-omics data from subjects to establish an analytical foundation. It further calculates the initial flux capacity corresponding to each enzymatic reaction, constructs a directed graph model, and characterizes the basic reaction capacity of different reaction pathways in the antioxidant metabolic network. This provides a structured analytical basis for subsequently extrapolating the overall reaction capacity of the metabolic network under different temperature conditions and identifying key enzymatic reactions that limit overall reaction capacity. Furthermore, at each temperature in the simulated temperature step sequence, the directed graph model is driven to synchronously perform bidirectional supply and demand rate updates. At the network level, the clearance requirements caused by reactive oxygen species generation and the upper limit of the reaction capacity that metabolic pathways can provide under the current temperature conditions are analyzed simultaneously. The failure temperature is locked, and the set of rate-limiting enzymatic reactions that limit global flux at the failure temperature is extracted. Key enzymatic catalytic reactions that constrain overall reaction capacity under specific conditions are identified, providing a basis for subsequent analysis. Finally, the set of rate-limiting enzymatic reactions is mapped to a combination of sensitive patient screening indicators, providing data analysis support for relevant personnel.

[0016] This invention acquires multi-omics data to construct a directed graph of enzyme-catalyzed reactions, uses a temperature step-by-step network to predict supply and demand and lock in failure temperatures, identifies rate-limiting enzyme-catalyzed reactions and maps and screens index combinations, thereby improving the accuracy of discovering biomarkers related to the efficacy of nanomedicines. Attached Figure Description

[0017] To more clearly illustrate the technical solutions and advantages in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0018] Figure 1 This is a system block diagram of a nanomedicine efficacy biomarker mining system based on multi-omics data provided in Embodiment 1 of the present invention; Figure 2 This is a flowchart of a method for synchronously executing supply and demand bidirectional rate updates in a directed graph model, as provided in Embodiment 4 of the present invention.

[0019] In the diagram, 101 is the preprocessing module, 102 is the deduction module, and 103 is the identifier output module. Detailed Implementation

[0020] To further illustrate the technical means and effects of the present invention in achieving its intended purpose, the following detailed description, in conjunction with the accompanying drawings and preferred embodiments, details the specific implementation, structure, features, and effects of a nanomedicine efficacy biomarker mining system based on multi-omics data proposed in accordance with the present invention.

[0021] In the following description, different "one embodiment" or "another embodiment" do not necessarily refer to the same embodiment; furthermore, specific features, structures or characteristics in one or more embodiments may be combined in any suitable form.

[0022] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0023] The following description, in conjunction with the accompanying drawings and embodiments, details a specific scheme for a nanomedicine efficacy biomarker mining system based on multi-omics data provided by the present invention.

[0024] Example 1 Please see Figure 1 The diagram shows a system block diagram of a nanomedicine efficacy biomarker mining system based on multi-omics data provided in this embodiment 1. The system includes: a preprocessing module 101, an inference module 102, and an identification output module 103.

[0025] The preprocessing module 101 is used to acquire multi-omics data of the subjects and calculate the initial flux capacity corresponding to each enzyme-catalyzed reaction based on a preset standard reference system. Based on the initial flux capacity, a directed graph model of the enzyme-catalyzed reaction is constructed. Each edge of the directed graph model corresponds to an enzyme-catalyzed reaction and is bound to the initial flux capacity and preset thermal inactivation decay coefficient corresponding to the enzyme-catalyzed reaction.

[0026] In this embodiment, multi-omics data of the subjects are first acquired to establish the basis for analysis. The subject's omics data collection must be performed under resting and fasting conditions. At the same time, in order to eliminate systematic errors caused by different testing batches and platforms, a preset standard reference system is loaded into the system. Then, based on the preset standard reference system, the initial flux capacity corresponding to each enzyme-catalyzed reaction is calculated. Based on the initial flux capacity, a directed graph model of the enzyme-catalyzed reaction is constructed. In this way, multi-omics data is transformed into flux capacity parameters that can quantify the potential reactivity of each enzyme-catalyzed reaction, thereby characterizing the basic reactivity of different reaction pathways in the antioxidant metabolic network and providing a unified quantitative basis for subsequent network-level reactivity analysis.

[0027] Considering that antioxidant metabolic pathways are essentially formed by multiple enzymatic reactions linked in series or parallel according to the substrate-to-product conversion relationship, each edge in the directed graph model corresponds to an enzymatic reaction, and the initial flux capacity and preset thermal inactivation decay coefficient of the enzymatic reaction are bound together. In this way, the reaction relationship in the metabolic pathway is mapped into a network structure with directional and capacity constraints. The initial flux capacity is used to quantify the maximum reactivity of each reaction pathway under baseline conditions, and the preset thermal inactivation decay coefficient is used to characterize the effect of temperature change on enzyme catalytic ability. Thus, an enzymatic reaction model that can simultaneously reflect the metabolic structure relationship and reactivity limitation is established at the network level, providing a structured analytical basis for subsequent extrapolation of the overall reactivity of the metabolic network under different temperature conditions and identification of key enzymatic reactions that limit the overall reactivity.

[0028] In this embodiment, the multi-omics data includes at least enzyme gene expression levels. substrate-to-product concentration ratio .

[0029] Specifically, the system receives raw data files from high-throughput sequencing platforms and metabolomics detection platforms.

[0030] For transcriptome data, the system reads RNA-Seq records from tumor tissue samples of the subjects and extracts transcript abundance values ​​(FPKM or TPM) of enzymes related to antioxidant metabolic pathways (such as glutathione metabolism and hydrogen peroxide scavenging pathway), which are recorded as enzyme gene expression levels. .

[0031] For metabolomics data, the system reads the mass spectrometry detection records of non-targeted metabolomics, extracts the chromatographic peak areas or relative concentrations of substrates and products in relevant metabolic reactions, and calculates the ratio of substrate concentration to product concentration for each pair of reactions, denoted as the substrate-product concentration ratio. .

[0032] Further loading of a pre-set standard reference system, whose data comes from measurements of standard cell lines (such as HeLa or normal fibroblasts) under a testing procedure completely parallel to that of the subject's sample, including baseline enzyme expression levels. Ratio of the reference metabolite The kinetic constants of each relevant enzyme are read from a pre-set enzyme parameter database (such as the BRENDA database), with a focus on obtaining the enzyme thermal inactivation decay coefficient. , which serves as the preset thermal inactivation decay coefficient; this coefficient describes the rate at which a specific enzyme loses activity when the temperature rises, and is a key parameter for subsequent simulation of thermal effects.

[0033] After data preparation is completed, the system enters the core modeling stage. Preferably, in this embodiment, in order to characterize the overall reactivity of different enzyme-catalyzed reactions in the antioxidant metabolic pathway from multi-omics data and to represent the complex metabolic flow processes within cells in computer memory, a directed graph model of enzyme-catalyzed reactions is constructed based on graph theory. The construction process includes: First, the metabolic reaction relationships in the antioxidant metabolic pathway are analyzed based on the pre-set biochemical pathway knowledge base. For example, the KEGG Pathway Database is selected as the biochemical pathway knowledge base. Then, the metabolites in the metabolic reaction are instantiated as nodes in the graph model, the enzyme-catalyzed reaction is instantiated as the directed edge connecting the nodes, and the direction of the directed edge is determined by the direction from the substrate to the product. Within the nodes, source nodes are set to represent the source of reducing power generation, such as the total intracellular NADPH production pool or the upstream metabolic substrate pool; sink nodes are set to represent the product pool after reactive oxygen species (ROS) are neutralized, i.e., the product pool after the neutralization of ROS generated by nanomedicines (such as water molecules). ); and set intermediate nodes located between the source node and the sink node to represent intermediate metabolites in the metabolic pathway, such as reduced glutathione (GSH), oxidized glutathione (GSSG), and hydrogen peroxide ( )wait; Construct a directed graph model of enzymatic reactions from source nodes to sink nodes based on nodes and directed edges.

[0034] Each edge Each edge uniquely corresponds to a specific enzyme-catalyzed reaction pathway, and the direction of the edge represents the direction of metabolite flow (from substrate to product). For example, a directed edge connecting the source node to the reduced glutathione node represents a reduction reaction catalyzed by glutathione reductase (GR); an edge connecting the reduced glutathione node and the hydrogen peroxide node and pointing to the sink node represents a redox reaction catalyzed by glutathione peroxidase (GPX).

[0035] After constructing the topological skeleton of the directed graph, the system assigns specific numerical weights and physical properties to each edge in the graph.

[0036] Preferably, in this embodiment, the expression level of enzyme genes and the substrate-product concentration ratio corresponding to each enzyme-catalyzed reaction are first obtained, and the basic data are matched with the edges in the directed graph model. Since there are usually individual differences in gene expression levels and metabolite concentration levels among different subjects, it is difficult to reflect the degree of change in responsiveness relative to the normal metabolic state by directly using raw values. Therefore, by introducing a preset standard reference system, the enzyme gene expression level of the subject is compared with the reference expression level to calculate the fold increase, and the substrate product concentration ratio is compared with the reference concentration ratio to calculate the fold change, thereby characterizing the degree of change in enzyme catalytic ability and metabolic reaction driving force relative to the standard metabolic state. Considering that the potential reactivity of enzyme-catalyzed reactions is affected by both the enzyme quantity level and changes in metabolic driving forces, this study characterizes the trend of enzyme catalytic capacity changes by expression folds, and characterizes the impact of substrate supply and product feedback on reaction propulsion by change folds. By combining the combined effects of both on reactivity, the initial flux capacity of each enzyme-catalyzed reaction is calculated, thereby quantifying the basic reactivity of different enzyme-catalyzed reactions in the current subject. This provides a unified capacity parameter for subsequent analysis of the contribution of each reaction pathway to the overall reactivity at the metabolic network level and for identifying key enzyme-catalyzed reactions that limit the overall reactivity.

[0037] The formula for calculating the initial flux capacity is shown in formula (1):

[0038] Where k is the index of the edge in the directed graph. This represents the initial flux capacity of the k-th edge; This represents the expression level of the enzyme gene corresponding to the k-th edge; This represents the substrate-product concentration ratio corresponding to the k-th edge; This represents the baseline enzyme expression level in the standard reference frame corresponding to the k-th edge. ; This represents the baseline metabolite ratio in the standard reference frame corresponding to the k-th edge; Indicates the characteristic adjustment factor, This is a preset positive parameter for division by zero, and its value is taken as... This is used to prevent the denominator from being zero, and can be adjusted according to the dimensions of the denominator. Dimensions; This is a truncation function that restricts the output to a preset range. Internally, to prevent numerical values ​​from jumping around drastically.

[0039] In the formula, The fold increase of the corresponding enzyme gene expression level relative to the reference expression level represents the gene expression contribution. The larger the ratio, the higher the transcriptional abundance of the enzyme in the subject's body, the greater the potential number of enzyme proteins, and thus the greater the reaction throughput. The ratio of substrate-to-product concentration to reference concentration represents the metabolic driving contribution, reflecting the degree of substrate accumulation relative to product in the subject's body. According to the law of mass action, the higher the substrate concentration (or the lower the product concentration), the greater the driving force for the chemical reaction to proceed towards the product.

[0040] This is an empirical index used to adjust the weight of metabolomics data in the overall score. A higher value indicates that the model considers the metabolic environment to have a more significant impact on flux. This was optimized through retrospective fitting analysis based on large-scale historical clinical cohort data and leave-one-out cross-validation. The value is typically between 0.5 and 0.8, and in this embodiment, it is 0.6. Since sequencing or detection errors may lead to extreme outliers (e.g., ratios as high as several thousand times), direct calculation would severely distort the results. Therefore, the system introduces a truncation function to forcibly limit the result of the exponentiation term to a preset range. In this example, take That is, if the calculated result is less than 0.1, then take 0.1; if it is greater than 10, then take 10. This not only prevents the numerical value from going too far, but also conforms to the objective law that there is a saturation effect in enzymatic reactions in organisms.

[0041] The system will calculate The value is assigned to the corresponding edge in the directed graph model as the initial capacity (capacity) of the k-th edge at the baseline temperature, i.e., the initial flux capacity; simultaneously, the system traverses each edge in the directed graph. Based on the enzyme type represented by the edge, the corresponding preset thermal inactivation attenuation coefficient is indexed from the database and bound to the edge. (Preset thermal inactivation attenuation coefficient (enzyme thermal inactivation attenuation coefficient)) The value describes the rate at which enzyme activity decreases with increasing temperature. The higher the value, the more sensitive the enzyme is to heat, and the faster its catalytic ability decreases during heating.

[0042] It should be noted that, since the rate of enzyme-catalyzed reaction and the thermal stability of enzyme proteins are significantly affected by temperature, the default reference temperature in this embodiment is set to 37 degrees Celsius (310.15K) to simulate the metabolic baseline under normal human physiological conditions. At this reference temperature, the system performs topology construction of the directed graph model, standard reference frame calibration, and initial flux capacity calibration. The analysis process for each edge is the same, and only one example is described here without repeating the explanation.

[0043] It should be noted that high-throughput sequencing platforms, metabolomics detection platforms, biochemical pathway knowledge bases (such as KEGG), enzyme parameter databases (such as BRENDA), and methods for acquiring, storing, and retrieving standard cell line data are all conventional techniques well-known to those skilled in the art. The methods for establishing standard reference frames and the specific algorithmic implementations for constructing directed graph models based on graph theory (such as instantiation of nodes and edges, attribute binding, etc.) are also existing mature technologies and will not be elaborated further.

[0044] Example 2 Based on Example 1, since subsequent steps require comparing the absolute rate of reactive oxygen species generated by the drug (unit: mol / s) with the relative flux index of enzyme catalysis (dimensionless) in the computational model, and since the two have different dimensions, they cannot be directly compared; therefore, a preset supply-demand alignment coefficient is obtained in advance here. It serves as a conversion bridge connecting two heterogeneous data spaces.

[0045] The method for setting the supply-demand alignment coefficient is a well-known technique among those skilled in the art. In this embodiment 2, the coefficient is calibrated using an in vitro experimental comparison method based on standard cell lines; the specific implementation is as follows: First, under standard in vitro experimental conditions (37℃), a unit dose of nanomedicine was administered to a standard cell line, and the maximum hydrogen peroxide generation rate that the cells could withstand to maintain redox balance was determined. During the determination process, the baseline cell clearance rate was first measured in a control group without nanomedicine, and then hydrogen peroxide was titrated to the experimental group until the cell survival critical point was reached. The amount of hydrogen peroxide consumed per unit time was calculated. This measured value is a specific physical quantity and is recorded as the absolute value of the measured ROS clearance rate. (Unit: mol / s)

[0046] Secondly, using the directed graph model constructed and the algorithm for calculating the initial flux capacity as described above, a computational model was built based on omics data from standard cell lines, and the initial flux capacity output of this model at 37°C was calculated; this calculated value is a relative value, denoted as the baseline flux index ( ). (dimensionless).

[0047] Finally, the system calculates the ratio of the two to obtain the supply-demand alignment coefficient. : .

[0048] This coefficient Stored as a global constant in the system, it is an engineered relative index (semi-empirical evaluation value) used for internal extrapolation within the computational model. Its purpose is to construct a mapping relationship between gene transcription levels and metabolic microenvironment driving forces within the computational space, rather than being directly equivalent to the absolute reaction rate within the organism. In subsequent calculations targeting specific subjects, the system will utilize the coefficient... The calculated relative flux index is mapped back to the absolute rate in physical space, thus giving the dimensionless calculation results a clear biochemical meaning.

[0049] Example 3 Based on Example 2, the deduction module 102 is used to generate a simulated temperature step sequence. For each temperature, it drives the directed graph model to synchronously perform supply and demand bidirectional rate updates to obtain the demand rate for reactive oxygen species neutralization. Based on the preset thermal inactivation decay coefficient, it deduces the flux upper limit of each side to obtain the global maximum flux at each temperature. By comparing the global maximum flux with the demand rate, the failure temperature is locked and the rate-limiting enzyme reaction set that limits the global flux at the failure temperature is extracted.

[0050] After constructing a graphical model with thermal response properties, the system enters the core dynamic simulation phase.

[0051] Considering that nanomedicines induce reactive oxygen species (ROS) generation during their action in vivo, and that the intracellular antioxidant metabolic network continuously consumes ROS through various enzyme-catalyzed reactions to maintain metabolic homeostasis, and that temperature changes affect enzyme catalytic activity and alter the overall reactivity of the metabolic network, a simulated temperature step sequence is generated. For each temperature, the directed graph model is driven to synchronously perform bidirectional supply and demand rate updates to obtain the required rate for ROS neutralization. Based on a preset thermal inactivation decay coefficient, the flux upper limit of each edge is deduced. Thus, at the network level, the clearance requirement brought about by ROS generation and the upper limit of the reactivity that the metabolic pathway can provide under the current temperature conditions are analyzed simultaneously. Then, the global maximum flux at each temperature is obtained to quantify the maximum overall reactivity that the entire antioxidant metabolic network can achieve under the current temperature conditions. Furthermore, considering that when the demand for neutralizing reactive oxygen species exceeds the overall responsiveness of the metabolic network, some reaction pathways in the network become key links limiting the overall responsiveness, by comparing the global maximum flux with the demand rate, we can identify the critical state in which the metabolic network changes from "able to meet the demand" to "unable to meet the demand," lock the failure temperature, extract the set of rate-limiting enzyme reactions that limit the global flux at the failure temperature, identify the key enzyme catalytic reactions that constrain the overall responsiveness under specific conditions, and provide a clear basis for extracting molecular indicators related to these enzyme reactions from multi-omics data and forming a combination of sensitive patient screening indicators.

[0052] Preferably, in this embodiment, the simulated temperature step sequence is set to a starting temperature of 310.15K and an ending temperature of 321.15K, corresponding to a temperature range of 37°C to 48°C, with a temperature step size of 0.2K. This generates an increasing simulated temperature step sequence, which serves as an external driving variable. Each temperature is sequentially input into the model to trigger subsequent rate updates and flux calculations. The analysis process is the same for each temperature, and only one example is described here without further explanation.

[0053] At the current temperature, the system needs to simultaneously calculate the attack rate on the demand side and the defense rate on the supply side to simulate the rate scissor effect in thermo-chemical synergistic therapy.

[0054] Example 4 Based on Example 3, please refer to Figure 2 This illustrates a flowchart of a method for synchronously updating the supply and demand rates of a directed graph model, as provided in Embodiment 4, specifically including: Step S201: Calculate the reactive oxygen species generation rate of the drug at each temperature based on the Arrhenius equation, combine the generation rate with a preset supply-demand alignment coefficient, and obtain the required rate for reactive oxygen species neutralization.

[0055] This step corresponds to demand-side updates, and the formula for calculating the demand rate is shown in (2):

[0056] in, This represents the rate at which reactive oxygen species are neutralized at the nth temperature; The current temperature is in Kelvin (K). This represents the initial hydrogen peroxide concentration, corresponding to the initial concentration before drug action; Indicates the drug pre-exponential factor; It is an exponential function with the natural constant e as the base; For activation energy, It is the ideal gas constant; This is a preset supply and demand alignment coefficient.

[0057] The formula corresponds to the engineering application of the Arrhenius equation (existing technology), by substituting the current simulation temperature. Kinetic parameters related to drug catalytic reactions ( , , as well as Simulations calculated the nanomedicine at a specific temperature. The theoretical absolute rate of reactive oxygen species (ROS) production under catalytic oxidation. Partially corresponding to the reactive oxygen species generation rate, by multiplying by a coefficient The theoretical drug attack rate derived from the computational model is adjusted to align with the metabolic network defense capability (global maximum flux) calculated based on the directed graph model in terms of both dimensions and physical meaning, thereby obtaining a demand rate for reactive oxygen species neutralization that can be directly used for supply and demand comparison within the unified framework of the model.

[0058] It should be noted that, This represents the baseline oxidative substrate level of the system before drug action, used to characterize the initial hydrogen peroxide concentration available to the nanomedicine in the reaction system, thereby reflecting the substrate supply conditions for potential reactive oxygen species generation reactions; in other embodiments of the invention, the implementer may adjust the dosage according to different experimental systems or pathological microenvironment conditions. Replace it with the measured hydrogen peroxide concentration in the tumor microenvironment, the initial oxidative substrate concentration in the in vitro culture system, or other equivalent oxidative substrate concentration parameters that can reflect the supply level of reactive oxygen species (ROS) generation substrates, such as total peroxide concentration and ROS precursor concentration, to adapt to the reaction kinetic modeling needs of different application scenarios.

[0059] Step S202: Based on the difference between each temperature and the reference temperature used in the initial flux capacity calculation, and combined with the preset thermal deactivation attenuation coefficient, obtain the attenuation factor for each side, fuse the initial flux capacity of each side and the corresponding attenuation factor at each temperature, and obtain the flux limit of each side at each temperature.

[0060] This step corresponds to supply-side updates. Since temperature changes affect the spatial structural stability and catalytic activity of enzyme proteins, and a preset thermal inactivation decay coefficient is used to characterize the sensitivity of different enzyme catalytic reactions to activity decay under increasing temperature conditions, the decay factor of each edge is obtained based on the difference between each temperature and the baseline temperature used to calculate the initial flux capacity, combined with the preset thermal inactivation decay coefficient. This represents the retention ratio of the enzyme catalytic capacity under the current temperature condition relative to the baseline temperature condition. By integrating the initial flux capacity of each edge and the corresponding decay factor at each temperature, the flux upper limit of each edge at each temperature is obtained. This ensures that each enzyme catalytic reaction edge has a reaction capacity constraint that matches the temperature state under the current temperature condition, thereby quantifying the maximum reaction flow that each reaction path in the metabolic network can provide under different temperature conditions, providing a basis for subsequent calculation of the overall maximum flux at the network level.

[0061] As an example, the formula for calculating the flux limit is shown in (3):

[0062] in, This represents the upper limit of flux for the k-th edge under the n-th temperature condition; This represents the initial flux capacity of the k-th edge; This represents the preset thermal deactivation attenuation coefficient of the k-th edge; The current temperature is in Kelvin (K). This represents the lower limit temperature of the simulated temperature step sequence; in this example, it is... ; The attenuation factor corresponding to the k-th edge under the n-th temperature condition.

[0063] In the formula, the flux upper limit of each side at each temperature is obtained by multiplication and fusion.

[0064] At this point, the system has generated a set of side capacities based on the current temperature. A temporary network instance that has undergone specific modifications.

[0065] It should be noted that in other embodiments of the present invention, the implementer may also set the starting temperature preferably to the actual measured temperature of the tumor site at rest in the subject, or 310.15K if the temperature is missing, and the termination temperature and step size may be adjusted according to the simulation accuracy and range.

[0066] Example 5 Based on Example 4, preferably in Example 5, a temporary network instance is constructed at each temperature based on the flux limit of each side. Considering that the antioxidant metabolic pathway is composed of multiple enzyme-catalyzed reaction pathways, and different reaction pathways may be in series or in parallel, when evaluating the reactive oxygen species scavenging capacity from the perspective of the overall network, it is necessary to comprehensively consider the capacity constraints of all reaction pathways. Therefore, the maximum flow algorithm is used to solve the overall reaction capacity of the metabolic network.

[0067] The maximum flow algorithm is used to iteratively search for augmenting paths from the source node to the sink node in a temporary network instance. Each augmenting path represents an effective metabolic reaction pathway that can participate in the reactive oxygen species (ROS) scavenging process. The bottleneck capacity of each augmenting path is accumulated to simulate the overall reaction capacity that different reaction pathways can provide by synergistically participating in the ROS scavenging process. This process continues until no augmenting paths are available, indicating that all metabolic pathways that can participate in the reaction have reached their capacity constraints. The accumulated result is then taken as the global maximum flux at the current temperature. The global maximum flux represents the maximum reactive oxygen species scavenging capacity that the antioxidant metabolic network can achieve under the current temperature conditions and enzyme activity constraints.

[0068] As an example: the execution steps include: ① Initialization: The system sets all edge traffic in the current temporary network instance to 0 and constructs the corresponding residual network (Residual Graph).

[0069] ② Searching for augmenting paths: The system uses breadth-first search (BFS) in the residual network to find a path from the source node. Remittance node The path is such that the remaining capacity of all edges on the path is greater than 0.

[0070] ③ Traffic augmentation: If such a path is found, the system finds the capacity value (bottleneck capacity) of the edge with the smallest remaining capacity on the path, adds this value to the total traffic, and updates the residual network (subtracts the capacity from forward edges and adds the capacity to reverse edges).

[0071] ④ Iterative loop: The system repeats steps ② and ③ until no path from the source node to the sink node can be found in the residual network.

[0072] It should be noted that, in other embodiments of the present invention, in order to avoid the cyclic flow phenomenon caused by reversible reactions in metabolic pathways, the dominant reaction direction can be determined according to the direction of the reaction Gibbs free energy, thereby constraining the direction of the directed edge of the enzyme-catalyzed reaction; the maximum flow algorithm and its iterative solution to the maximum flux allocation process, the method for determining the direction of the reaction Gibbs free energy, etc. are all conventional technical means known to those skilled in the art, and will not be described in detail here.

[0073] Furthermore, the temperature is increased sequentially according to the simulated temperature step sequence. After obtaining the global maximum flux and demand rate at each temperature, the system immediately performs a supply-demand comparison to determine whether the metabolic network has experienced functional failure.

[0074] Compare the global maximum throughput with the demand rate. If the global maximum throughput is greater than or equal to the demand rate, it means that the network still has the capacity to remove ROS generated by the drug. The system then enters the next temperature step point and calculates new parameters. If the global maximum throughput is less than the demand rate, it indicates that the demand pressure has exceeded the network's supply limit. The system determines that the metabolic throughput has failed and identifies the current temperature as the failure temperature. .

[0075] After obtaining the failure temperature, the temperature iteration is terminated, and the sensitive state flag is set to sensitive. Then, in the network instance corresponding to the failure temperature, a residual network is constructed based on the calculation results of the maximum flow algorithm, and the corresponding minimum cut edge set is determined based on the Max-Flow Min-Cut Theorem. The minimum cut edge set constitutes the "structural bottleneck" from the source node to the sink node. Its physical meaning is that the directed edges in this set are all in a flow saturation state at the current temperature (i.e., the actual flux is equal to the upper limit of the flux after temperature correction), and removing this set will cause the source-sink connectivity to be interrupted and the global flux to zero. The enzyme-catalyzed reaction corresponding to the directed edge in the minimum cut edge set is determined as the rate-limiting enzyme-catalyzed reaction; if there are multiple minimum cut edge sets that meet the conditions, the system selects the set with the fewest directed edges according to preset rules; if the number is the same, the set with the smallest sum of edge indices is selected to ensure the determinism of the output result. By summarizing the rate-limited enzyme reactions, we obtain the set of rate-limited enzyme reactions.

[0076] It should be noted that the specific algorithmic implementations of using the maximum flow algorithm to solve for the maximum flow in a network and identifying the minimum cut edge set based on the residual network through breadth-first search (BFS) or depth-first search (DFS) (such as the Edmonds-Karp algorithm, the Dinic algorithm, and their variants) are all well-known and mature common knowledge in the fields of graph theory and computer networks, and will not be elaborated further.

[0077] It should be noted that if the failure temperature is not found after iterating to the upper limit of the simulated temperature step sequence, the sensitive state flag will be set to non-sensitive, and the power calculation and screening index combination process in the subsequent flag output module 103 will be terminated. Only the result of the sensitive state flag will be output to provide reference information for adjusting clinical treatment strategies or formulating combination drug regimens.

[0078] Example 6 Based on Example 5, the identifier output module 103 is used to map the rate-limited enzyme reaction set to a set of sensitive patient screening indicators.

[0079] By mapping the set of rate-limiting enzyme reactions to a combination of sensitive patient screening indicators, and analyzing the changes in the expression of genes or metabolic characteristics related to the corresponding enzyme reactions in multi-omics data, a group of subjects with similar metabolic rate-limiting structures can be identified, thereby providing data analysis basis for relevant personnel and realizing the correlation analysis process from metabolic mechanism analysis to population characteristic identification.

[0080] Preferably, in this embodiment 5, each edge in the rate-limiting enzyme reaction set is analyzed to identify the enzyme type corresponding to each edge (e.g., glucose-6-phosphate dehydrogenase G6PD or glutathione peroxidase GPX4). By retrospectively analyzing the multi-omics data of the subjects, the expression levels of enzyme genes and the concentration ratios of substrates and products for each enzyme were extracted. The enzyme type, its corresponding enzyme gene expression level, and substrate product concentration ratio are combined to output as a set of sensitive patient screening indicators.

[0081] It should be noted that, in another embodiment of the present invention, in order to realize the correlation analysis between the identification results of the vulnerability of the in vivo metabolic network and the external energy input conditions, it is necessary to establish a physical correspondence between the critical temperature at which the metabolic system becomes dysfunctional and the energy input conditions; considering that nanomedicines convert light energy into local heat energy through photothermal conversion under laser irradiation, the actual increase in tissue temperature depends not only on the photothermal conversion efficiency of the drug, but also on the combined influence of tissue heat diffusion and heat dissipation capacity; Therefore, after determining the failure temperature at which the metabolic network's defense capabilities fail, the identification output module 103 also uses the failure temperature and the preset tissue heat dissipation coefficient, combined with the photothermal conversion efficiency of the drug, to reverse-derive the preset power density reference value of the laser treatment device. This allows the acquisition of the energy input level that enables the tissue temperature to reach the failure temperature under model conditions, thereby realizing a quantitative correspondence between the metabolic network vulnerability analysis results and external energy conditions, and providing a reference for the assessment and research of relevant photothermal effects.

[0082] First, based on the failure temperature With preset safety redundancy threshold Determine the target temperature , The safety redundancy threshold can be set as needed. The preset safety redundancy threshold of this invention is set to 0.5K-2.0K. Calculate the temperature difference between the target temperature and the preset reference temperature, and combine it with the preset tissue heat dissipation coefficient. (unit: ), to obtain the thermal power requirement needed to maintain the target temperature; Treating thermal power requirement as a molecule, the photothermal conversion efficiency of the drug... and light absorption coefficient The product of the fractions is used as the denominator, and the ratio of the fractions is used as the preset power density reference value for the laser treatment equipment.

[0083] As an example, the preset power density reference value The calculation formula is shown in (4):

[0084] The preset reference temperature is the reference temperature used for initial flux capacity calculation, corresponding to 310.15K. This is a preset positive parameter for division by zero; The corresponding thermal power requirement to maintain the target temperature; Expressed as the photothermal conversion efficiency of the drug (dimensionless); denoted as the light absorption coefficient of the drug, and denoted as a dimensionless percentage. The photothermal conversion efficiency of the corresponding drug.

[0085] In the formula, the photothermal conversion efficiency is... and light absorption coefficient The two parameters are determined by the physicochemical properties of the drug itself and will not be further specified here. The preset tissue heat dissipation coefficient is set as an empirical constant for engineering estimation, and its value is based on the statistical mean of the average blood perfusion rate of biological tissues under resting conditions. Those skilled in the art can calculate it as needed, and in this example, it is taken as 0.05. .

[0086] In summary, in order to address the technical problem that existing technologies struggle to accurately identify key enzymatic reactions that limit the overall responsiveness of antioxidant metabolic networks based on multi-omics data, thus making it difficult to determine the combination of sensitive patient screening indicators related to the efficacy of nanomedicines, the present invention aims to provide a nanomedicine efficacy biomarker mining system based on multi-omics data.

[0087] This invention first acquires multi-omics data of subjects, calculates the initial flux capacity corresponding to each enzymatic reaction, and constructs a directed graph model; further, at each temperature of the simulated temperature step sequence, it drives the directed graph model to synchronously perform bidirectional supply and demand rate updates, locks the failure temperature, and extracts the rate-limiting enzymatic reaction set that limits the global flux at the failure temperature; finally, it maps the rate-limiting enzymatic reaction set to a combination of sensitive patient screening indicators.

[0088] This invention acquires multi-omics data to construct a directed graph of enzyme-catalyzed reactions, uses a temperature step-by-step network to predict supply and demand and lock in failure temperatures, identifies rate-limiting enzyme-catalyzed reactions and maps and screens index combinations, thereby improving the accuracy of discovering biomarkers related to the efficacy of nanomedicines.

[0089] It should be noted that the order of the above embodiments of the present invention is merely for description and does not represent the superiority or inferiority of the embodiments; the processes depicted in the accompanying drawings do not necessarily require a specific or sequential order to achieve the desired results; in some embodiments, multitasking and parallel processing are also possible or may be advantageous.

[0090] The various embodiments in this specification are described in a progressive manner. The same or similar parts between the various embodiments can be referred to each other. Each embodiment focuses on describing the differences from other embodiments.

Claims

1. A system for mining biomarkers for the therapeutic efficacy of nanomedicines based on multi-omics data, characterized in that, include: The system comprises a preprocessing module, a deduction module, and an identification output module. The preprocessing module acquires multi-omics data from the subjects, calculates the initial flux capacity corresponding to each enzymatic reaction, constructs a directed graph model, and inputs it into the deduction module. The deduction module drives the directed graph model to synchronously perform bidirectional supply and demand rate updates at each temperature in the simulated temperature step sequence, obtains the demand rate for reactive oxygen species neutralization and the global maximum flux at each temperature, locks the failure temperature, and extracts the rate-limiting enzymatic reaction set that restricts the global flux at the failure temperature. The identification output module maps the rate-limiting enzymatic reaction set to a combination of sensitive patient screening indicators.

2. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 1, characterized in that, The preprocessing module acquires multi-omics data from the subjects and calculates the initial flux capacity corresponding to each enzyme-catalyzed reaction based on a preset standard reference system. It then constructs a directed graph model of the enzyme-catalyzed reaction based on the initial flux capacity. Each edge of the directed graph model corresponds to an enzyme-catalyzed reaction and is bound to the initial flux capacity and preset thermal inactivation decay coefficient corresponding to the enzyme-catalyzed reaction.

3. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 2, characterized in that, The method for obtaining the directed graph model includes: The metabolic reaction relationships in the antioxidant metabolic pathway are analyzed based on a pre-built biochemical pathway knowledge base. Metabolic products in metabolic reactions are instantiated as nodes in a directed graph model, enzyme-catalyzed reactions are instantiated as directed edges connecting nodes, and the direction of the directed edges is determined by the direction from the substrate to the product. A source node is set in the node to represent the source of reducing power generation; a sink node is set to represent the product pool after reactive oxygen species are neutralized; and an intermediate node is set between the source node and the sink node to represent intermediate metabolites in the metabolic pathway. Construct an enzymatic reaction directed graph model from the source node to the sink node based on the nodes and the directed edges.

4. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 3, characterized in that, The method for obtaining the initial flux capacity includes: Obtain the enzyme gene expression level and substrate-product concentration ratio corresponding to each enzyme-catalyzed reaction; Based on a preset standard reference system, the fold increase of the enzyme gene expression level relative to the reference expression level and the fold increase of the substrate product concentration ratio relative to the reference concentration ratio are calculated. The initial flux capacity of each enzymatic reaction is calculated based on the expression fold and the change fold.

5. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 4, characterized in that, The inference module generates a simulated temperature step sequence. For each temperature, it drives the directed graph model to synchronously perform bidirectional supply and demand rate updates to obtain the demand rate for reactive oxygen species neutralization. Based on the preset thermal deactivation decay coefficient, it infers the flux upper limit of each side to obtain the global maximum flux at each temperature. By comparing the global maximum flux with the required rate, the failure temperature is determined and the set of rate-limiting enzyme reactions that restrict the global flux at the failure temperature is extracted.

6. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 5, characterized in that, The methods for obtaining the flux upper limit include: Based on the difference between each temperature and the reference temperature used to calculate the initial flux capacity, and combined with the preset thermal deactivation attenuation coefficient, the attenuation factor of each side is obtained. The initial flux capacity of each side and the attenuation factor corresponding to each temperature are then fused together to obtain the flux limit of each side at each temperature. The method for obtaining the demand rate includes: The reactive oxygen species (ROS) generation rate of the drug at each temperature is calculated based on the Arrhenius equation. The generation rate is then combined with a preset supply-demand alignment coefficient to obtain the required rate for ROS neutralization.

7. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 6, characterized in that, The method for obtaining the global maximum flux includes: At each temperature, a temporary network instance is constructed based on the flux upper limit of each edge of the directed graph model. The maximum flow algorithm is used to iteratively search for augmenting paths from the source node to the sink node in the temporary network instance. The bottleneck capacity of each augmenting path is accumulated until there is no augmenting path. The accumulated result is used as the global maximum flux at the current temperature.

8. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 7, characterized in that, The method for obtaining the failure temperature includes: The temperature is increased sequentially according to the simulated temperature step sequence. At each temperature, the corresponding global maximum flux and the required rate are obtained. The global maximum flux and the required rate are compared. If the global maximum flux is less than the required rate, the current temperature is determined to be the failure temperature.

9. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 8, characterized in that, The identifier output module maps the rate-limiting enzyme reaction set to a combination of sensitive patient screening indicators, as follows: The method for obtaining the rate-limiting enzyme reaction set includes: In the network instance corresponding to the failure temperature, a residual network is constructed based on the calculation results of the maximum flow algorithm, and the corresponding minimum cut edge set is determined based on the maximum flow minimum cut theorem. The enzymatic reactions corresponding to the directed edges in the minimum cut edge set are identified as rate-limiting enzymatic reactions. By summarizing the rate-limiting enzyme reactions, the set of rate-limiting enzyme reactions is obtained; Methods for mapping the rate-limiting enzyme reaction set to a combination of sensitive patient screening indicators include: Analyze each edge in the rate-limiting enzyme reaction set and identify the enzyme type corresponding to each edge; By retracing the multi-omics data of the subjects, the expression level of enzyme genes and the substrate-product concentration ratio of each enzyme were extracted; The enzyme types and their corresponding enzyme gene expression levels and substrate-product concentration ratios are combined and output as a set of sensitive patient screening indicators.

10. The nanomedicine efficacy biomarker mining system based on multi-omics data according to claim 9, characterized in that, The identification output module, based on the failure temperature and the preset tissue heat dissipation coefficient, combined with the photothermal conversion efficiency of the drug, reverse-derives the preset power density reference value of the laser treatment device. The method for obtaining the preset power density reference value includes: The target temperature is determined based on the failure temperature and the preset safety redundancy threshold. Calculate the temperature difference between the target temperature and the preset reference temperature, and combine it with the preset tissue heat dissipation coefficient to obtain the heat power requirement required to maintain the target temperature; The thermal power requirement is used as the numerator, and the product of the drug's photothermal conversion efficiency and light absorption coefficient is used as the denominator. The ratio of the fractions is used as the preset power density reference value of the laser treatment device.

Images