Intelligent Analysis Method and System for Endothelial Cell Heterogeneity in Atherosclerosis
Through single-cell sequencing technology and intelligent analysis methods, the subpopulation differences of endothelial cells in atherosclerosis are accurately identified, which solves the problem of insufficient understanding of endothelial cell heterogeneity in traditional methods, and achieves in-depth research on the molecular mechanism of the disease and the discovery of therapeutic targets.
Patent Information
- Application Number
- CN202510145547.9
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-02-10
- Publication Date
- 2025-07-04
- Estimated Expiration
- 2045-02-10
AI Technical Summary
The prior art is difficult to accurately identify and characterize the heterogeneity of endothelial cells in atherosclerosis, resulting in insufficient understanding of the molecular mechanism of the disease and limiting the discovery of potential molecular markers and therapeutic targets.
Gene expression information of endothelial cells was obtained through single-cell sequencing technology, and subpopulations were divided using Gaussian nuclear similarity and density peak clustering methods. Combined with support vector machine regression model and gene set enrichment analysis, significantly differentially expressed genes were identified and potential drug targets were verified.
Accurate analysis of endothelial cell heterogeneity is achieved, functional changes and potential molecular markers at the single-cell level, comprehensively understand the molecular mechanisms of the disease and discover new therapeutic targets.
Smart Images

Figure CN119580839B_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical field of data analysis, and particularly to an intelligent analysis of endothelial cell heterogeneity in atherosclerosis. Background Art
[0002] In recent years, in the field of atherosclerosis research, traditional analysis methods mainly rely on bulk RNA sequencing (bulkRNA-seq) technology. By comprehensively analyzing the gene expression of whole tissue samples, the overall gene expression changes under disease states are revealed. However, this method is based on the average expression level after mixing multiple cells, making it difficult to capture the differences between individual cells and unable to effectively reflect the heterogeneity characteristics within the cell population, especially within the endothelial cell population. In addition, traditional analysis means usually focus on the detection of known markers or preset target genes, which greatly limits the discovery of new subpopulations and potential biomarkers.
[0003] Atherosclerosis is a highly complex vascular disease, and its pathological mechanism involves multiple cell types and their interactions. Among them, endothelial cells play an important role as key regulators in the occurrence and development of the disease. However, gene expression analysis at the tissue level is difficult to clearly reveal the functional changes and molecular characteristic differences of specific endothelial cell subpopulations during the disease process, resulting in insufficient understanding of the pathological mechanism. Existing studies have begun to attempt to use single-cell sequencing technology to analyze cell diversity in complex tissues, providing a new means for revealing endothelial cell heterogeneity in atherosclerosis. However, the current technology has significant limitations: the processing and analysis process of single-cell data is complex, and an efficient and systematic method has not been formed, making it difficult to accurately identify and characterize endothelial cell subpopulations and their specific roles in the dynamic changes of the disease.
[0004] In particular, for single-cell sequencing studies on atherosclerosis-related endothelial cells, existing technologies mostly focus on preliminary exploration and lack a highly targeted and high-resolution analysis strategy. The significant heterogeneity within the endothelial cell population, such as the differences in function and molecular level among different subpopulations, is often ignored. This deficiency directly limits the comprehensive understanding of the role of endothelial cells in the atherosclerotic process and hinders the discovery of potential molecular markers and new therapeutic targets.
[0005] Therefore, the existing technical solutions are mainly limited to gene expression analysis at the tissue level or preliminary single-cell sequencing attempts, lacking the ability to accurately analyze the heterogeneity of atherosclerotic endothelial cells. This deficiency not only restricts the in-depth study of the molecular mechanism of the disease but also makes it difficult to support the development of precision treatment strategies based on cell subpopulations. There is an urgent need for an innovative method to fill the above technical gaps. Summary of the Invention
[0006] An embodiment of the present invention provides a method and system for intelligent analysis of endothelial cell heterogeneity in atherosclerosis, which can accurately analyze the endothelial cell heterogeneity in atherosclerosis.
[0007] An embodiment of the present invention provides a method for intelligent analysis of endothelial cell heterogeneity in atherosclerosis, including:
[0008] Obtaining sample data of endothelial cells in atherosclerosis of a patient; the sample data includes the original sequencing data of the gene expression information of single endothelial cells, and the original sequencing data is obtained by pre - analyzing the transcriptome level of endothelial cells using single - cell sequencing technology;
[0009] Performing feature extraction on the sample data of the endothelial cells to obtain the gene expression features of the endothelial cells;
[0010] Dividing the endothelial cells into different sub - populations according to the gene expression features of the endothelial cells and the gene expression similarity, determining marker genes for distinguishing different sub - populations, and defining and characterizing the features of each sub - population by analyzing the expression patterns of the marker genes;
[0011] Comparing the gene expression differences between different endothelial cell sub - populations, and identifying differentially expressed genes that are significantly up - regulated or down - regulated in the process of atherosclerosis in the relevant endothelial cell sub - populations, and the differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells.
[0012] As an improvement of the above solution, the step of dividing the endothelial cells into different sub - populations according to the gene expression features of the endothelial cells and the gene expression similarity, determining marker genes for distinguishing different sub - populations, and defining and characterizing the features of each sub - population by analyzing the expression patterns of the marker genes includes:
[0013] Similarity measurement and probability calculation steps: determining the vector of cells in the high - dimensional gene expression feature space with the obtained gene expression features, calculating the Gaussian kernel similarity between cells, and calculating the conditional probability and joint probability based on this similarity;
[0014] Low-dimensional vector layout optimization steps: In a low-dimensional space of a preset dimension, set low-dimensional vectors corresponding to high-dimensional cell vectors, and determine the low-dimensional vector layout by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix, and joint probability based on gene expression feature data, then randomly initialize the low-dimensional vectors so that they are evenly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. Calculate the gradient of the objective function with respect to the low-dimensional vectors at each iteration and update the low-dimensional vectors accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule based on the change of the objective function during the iteration. The number of iterations is set within a predetermined range and convergence is judged based on the comparison of the change amount of the objective function with a preset convergence threshold;
[0015] Density peak-based clustering steps: After obtaining the optimized low-dimensional vectors, use a density peak-based clustering algorithm in the low-dimensional space to calculate the local density of cells and the minimum distance to cells with higher density. The local density calculation is based on the Euclidean distance between cells and combines statistical parameters of the distance matrix to determine the calculation scale. For cases where there are no cells with higher density, set the minimum distance to a predetermined large value. Construct a decision graph based on the local density and the minimum distance and select cells with predetermined characteristics as cluster centers, and assign the remaining cells to the clusters to which the nearest cluster centers belong, thereby realizing the division of endothelial cell subsets;
[0016] Marker gene determination steps: For each subset after division, calculate the change multiple of the average expression level ratio of genes within the subset to genes in other subsets. Set a threshold based on this ratio and the average expression level of genes within the subset to screen marker genes, and use the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
[0017] As an improvement to the above solution, comparing the gene expression differences between different endothelial cell subsets, identifying the differentially expressed genes that are significantly up-regulated or down-regulated in relevant endothelial cell subsets during the process of atherosclerosis, and using the differentially expressed genes as key molecular markers to understand the functional changes and heterogeneity of endothelial cells, including:
[0018] Subset gene expression level calculation sub-step: Based on the characteristics of the marker genes of different subsets, for two selected endothelial cell subsets, count the number of cells within each subset respectively, and then calculate the expression level of genes in the corresponding subsets, which is used as the basis data preparation for subsequent model construction and analysis;
[0019] Model construction and parameter estimation operator steps based on negative binomial distribution: Construct a negative binomial distribution model suitable for describing the gene expression level distribution. The model expression contains terms associated with cell characteristic variables, model parameters to be determined, and discrete parameters unique to the negative binomial distribution; Use the maximum likelihood estimation method to carry out model parameter estimation work to provide a basis for subsequent gene expression level estimation;
[0020] Support vector machine regression model construction sub-steps: Set the training data to consist of cell characteristic vectors and gene expression levels, and construct a support vector machine regression model. The objective function of the support vector machine regression model comprehensively considers the relevant norms of the weight vector and terms related to slack variables, and at the same time sets a constraint condition system constructed based on cell characteristic vectors, weight vectors, bias terms, slack variables, and insensitive loss parameters; By solving the corresponding quadratic programming problem, determine the parameters of the support vector machine regression model;
[0021] Differential expression analysis and screening sub-steps: Use the gene expression level data estimated by the negative binomial distribution model as the input information of the support vector machine regression model, and use the support vector machine regression model to deeply learn and analyze the gene expression difference relationship between cell subpopulations; Calculate the differential expression fold change of genes between different subpopulations based on the adjusted expression level data output by the support vector machine regression model;
[0022] Testing sub-steps: Use a testing method based on normal distribution approximation to perform a significance test on the calculated differential expression fold change, and screen out genes with significant differential expression according to the preset significance level threshold.
[0023] As an improvement of the above solution, after comparing the gene expression differences between different endothelial cell subpopulations and identifying the differentially expressed genes that are significantly upregulated or downregulated in the process of atherosclerosis in the relevant endothelial cell subpopulations, and using the differentially expressed genes as key molecular markers to understand the functional changes and heterogeneity of endothelial cells, the method further includes:
[0024] Select the protein products corresponding to the significantly changed differentially expressed genes as potential drug targets, and perform simulated functional verification on the selected potential drug targets.
[0025] As an improvement of the above solution, the step of selecting the protein products corresponding to the significantly changed differentially expressed genes as potential drug targets and performing simulated functional verification on the selected potential drug targets includes:
[0026] Sub-step for screening and ranking differentially expressed genes: For the obtained set of differentially expressed genes, calculate the fold change of differential expression and statistical test values of each gene through statistical test methods, and construct a comprehensive scoring function in combination with the gene function importance index. This function importance index is obtained based on a multi-factor comprehensive evaluation of the degree of gene participation in key biological pathways and the degree of association with known disease-related genes. Rank the differentially expressed genes in descending order according to the comprehensive scoring function;
[0027] Sub-step for protein product analysis and preliminary target screening: Use bioinformatics tools to obtain the information of the protein product set corresponding to the ranked differentially expressed genes. Based on multiple screening factors including the participation in cell signal transduction, the interaction with lipid metabolism and inflammatory factors, and the localization characteristics of membrane proteins in subcellular compartments of the protein's functional characteristics, construct a screening function that includes the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor. Screen out the candidate set of potential drug targets through this screening function in combination with a set threshold. Among them, the weight of each screening factor is set according to its relative importance in drug target screening, and the scoring function determines the value according to the matching degree between the protein characteristics and the screening criteria, so as to achieve the preliminary screening of protein products to determine potential drug target candidates;
[0028] Sub-step for calculating binding affinity based on molecular docking: For each potential drug target, select candidate ligand molecules from a specific molecular library. First, perform structural preprocessing on the target protein and ligand molecules, then perform conformational search and docking in the active pocket region of the target protein. Use the molecular docking algorithm formula that includes the number of atoms of the ligand and target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the control angle range parameters to calculate the binding energy. Determine the partial charges and van der Waals interaction parameters according to the atomic chemical element type, set the electrostatic interaction adjustment parameter according to experience, and set the angle range parameter according to the characteristics of the protein active pocket. Screen out the ligand-target combinations with strong binding affinity by comparing the binding energy with a preset threshold to preliminarily verify the druggability of the target;
[0029] Sub-step for target stability analysis based on molecular dynamics simulation: For the ligand-target combinations with good binding energy, construct a simulation system that includes the ligand, target protein, and solvent molecules. Determine the atomic mass parameters according to the atomic chemical composition, set the atomic potential energy function and interatomic interaction forces according to the molecular mechanics force field parameters, determine the external forces according to the simulation environment temperature and pressure conditions, use the molecular dynamics simulation algorithm formula for simulation, set the simulation time and step size, monitor the root mean square deviation of the target protein structure stability index, and judge the structural stability of the target protein after binding to the ligand based on the comparison of the root mean square deviation with a preset range to further verify the effectiveness of the target;
[0030] Sub - steps for verifying the target regulation relationship based on network analysis: Construct a gene regulation network including differentially expressed genes, potential drug targets, and other related genes. The network nodes are genes and targets, and the edges represent regulation relationships with weight assignments. For each potential drug target, determine the gene set with regulation relationships with the target and the shortest path length between the gene and the target according to the network topology structure. Use the network analysis algorithm formula to calculate the importance index of the target in the gene regulation network. Judge the key degree of the target in the gene regulation network based on the comparison between this index and a preset threshold, so as to verify the rationality of the target from the gene regulation level.
[0031] Another embodiment of the present invention correspondingly provides an intelligent analysis system for endothelial cell heterogeneity in atherosclerosis, including:
[0032] An acquisition module, configured to acquire sample data of endothelial cells in the atherosclerosis of a patient; the sample data includes the original sequencing data of the gene expression information of a single endothelial cell, and the original sequencing data is obtained by pre - analyzing the transcriptome level of endothelial cells using single - cell sequencing technology;
[0033] A feature extraction module, configured to extract features from the sample data of the endothelial cells to obtain the gene expression features of the endothelial cells;
[0034] A sub - population division module, configured to divide the endothelial cells into different sub - populations according to the gene expression features of the endothelial cells and according to gene expression similarity, to determine marker genes for distinguishing different sub - populations, and to define and characterize the features of each sub - population by analyzing the expression patterns of the marker genes;
[0035] A differential analysis module, configured to compare the gene expression differences between different endothelial cell sub - populations, and identify differentially expressed genes that are significantly up - regulated or down - regulated in the process of atherosclerosis in the relevant endothelial cell sub - populations, and the differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells.
[0036] As an improvement of the above - mentioned solution, the sub - population division module is specifically configured to:
[0037] Determine the vector of the cell in the high - dimensional gene expression feature space with the obtained gene expression features, calculate the Gaussian kernel similarity between cells, and calculate the conditional probability and joint probability based on this similarity;
[0038] In a low-dimensional space of a preset dimension, a low-dimensional vector is set corresponding to a high-dimensional cell vector, and the layout of the low-dimensional vector is determined by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix, and joint probability according to the gene expression feature data, then randomly initialize the low-dimensional vector so that it is uniformly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. At each iteration, calculate the gradient of the objective function with respect to the low-dimensional vector and update the low-dimensional vector accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule based on the change of the objective function during the iteration. The number of iterations is set within a predetermined range and it is judged whether to converge according to the comparison between the change amount of the objective function and the preset convergence threshold;
[0039] After obtaining the optimized low-dimensional vector, use the density peak-based clustering algorithm in the low-dimensional space to calculate the local density of the cells and the minimum distance to the cells with higher density. The local density calculation is based on the Euclidean distance between cells and combines the statistical parameters of the distance matrix to determine the calculation scale. For the case where there are no cells with higher density, set the minimum distance to a predetermined large value. Construct a decision graph based on the local density and the minimum distance and select the cells with predetermined characteristics as the clustering centers, and assign the remaining cells to the clusters belonging to the nearest clustering centers, so as to achieve the division of endothelial cell subsets;
[0040] For each subgroup after division, calculate the change multiple of the average expression level ratio of the genes within the subgroup to the genes of other subgroups. Set a threshold to screen the marker genes based on this ratio and the average expression level of the genes within the subgroup, and use the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
[0041] As an improvement of the above solution, the differential analysis module is specifically used for:
[0042] Based on the characteristics of the marker genes of different subgroups, for the selected two endothelial cell subgroups, respectively count the number of cells within each subgroup, and then calculate the expression level of the genes in the corresponding subgroups, which is used as the basis data preparation for subsequent model construction and analysis;
[0043] Construct a negative binomial distribution model suitable for describing the gene expression level distribution. The model expression contains terms associated with cell characteristic variables, as well as model parameters to be determined and discrete parameters unique to the negative binomial distribution; use the maximum likelihood estimation method to carry out model parameter estimation work to provide a basis for subsequent gene expression level estimation;
[0044] The set training data consists of cell feature vectors and gene expression levels. A support vector machine regression model is constructed. The objective function of the support vector machine regression model comprehensively considers the relevant norms of the weight vector and the terms related to the slack variables. At the same time, a constraint condition system based on the cell feature vector, weight vector, bias term, slack variables, and insensitive loss parameter is set up; by solving the corresponding quadratic programming problem, the parameters of the support vector machine regression model are determined;
[0045] The gene expression level data estimated by the negative binomial distribution model is used as the input information of the support vector machine regression model. The support vector machine regression model is used to deeply learn and analyze the gene expression difference relationship between cell subsets; based on the adjusted expression level data output by the support vector machine regression model, the differential expression fold of the gene between different subsets is calculated;
[0046] A significance test is performed on the calculated differential expression fold using a test method based on the normal distribution approximation. Genes with significant differential expression are screened out according to the preset significance level threshold.
[0047] As an improvement to the above solution, the system further includes:
[0048] A drug target verification module, which is used to select the protein products corresponding to the significantly differentially expressed genes as potential drug targets and perform simulated functional verification on the selected potential drug targets.
[0049] As an improvement to the above solution, for the obtained set of differentially expressed genes, the differential expression fold change and statistical test values of each gene are calculated through statistical test methods, and a comprehensive scoring function is constructed in combination with the gene function importance index. This function importance index is obtained based on a multi-factor comprehensive evaluation of the degree of participation of the gene in biological pathways and the degree of association with known disease-related genes. The differentially expressed genes are sorted in descending order according to the comprehensive scoring function;
[0050] Using bioinformatics tools to obtain the information of the protein product set corresponding to the differentially expressed genes after sorting. Based on the functional characteristics of the protein, including its participation in cell signal transduction, interaction with lipid metabolism and inflammatory factors, and subcellular localization characteristics such as membrane protein localization, a screening function including the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor is constructed. Through this screening function and a set threshold, a candidate set of potential drug targets is screened out. The weight of each screening factor is set according to its relative importance in drug target screening, and the scoring function determines the value according to the matching degree between the protein characteristics and the screening criteria, so as to realize the preliminary screening of protein products to determine potential drug target candidates;
[0051] For each potential drug target, candidate ligand molecules are selected from a specific molecular library. First, the structures of the target protein and ligand molecules are preprocessed, and then conformational search and docking are performed in the active pocket region of the target protein. The binding energy is calculated using a molecular docking algorithm formula that includes parameters such as the number of atoms of the ligand and target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the parameter for controlling the angle range. The partial charges and van der Waals interaction parameters are determined based on the atomic chemical element types, the electrostatic interaction adjustment parameter is set according to experience, and the angle range parameter is set based on the characteristics of the protein active pocket. Ligand-target combinations with strong binding affinity are screened out by comparing the binding energy with a preset threshold to preliminarily verify the druggability of the target;
[0052] For ligand-target combinations with good binding energy, a simulation system containing the ligand, target protein, and solvent molecules is constructed. The atomic mass parameters are determined based on the atomic chemical composition, the atomic potential energy function and the interatomic interaction forces are calculated according to the molecular mechanics force field parameters, the external forces are determined according to the simulation environment temperature and pressure conditions, and the simulation is carried out using the molecular dynamics simulation algorithm formula. The simulation time and step size are set, and the root mean square deviation (RMSD), an indicator of the structural stability of the target protein, is monitored. The structural stability of the target protein after binding to the ligand is judged based on the comparison of the RMSD with a preset range to further verify the effectiveness of the target;
[0053] A gene regulatory network containing differentially expressed genes, potential drug targets, and other related genes is constructed. The network nodes are genes and targets, and the edges represent regulatory relationships with weight assignments. For each potential drug target, the gene set with regulatory relationships with the target and the shortest path length between the gene and the target are determined based on the network topology structure. The importance index of the target in the gene regulatory network is calculated using a network analysis algorithm formula, and the key degree of the target in the gene regulatory network is judged based on the comparison of this index with a preset threshold, so as to verify the rationality of the target from the gene regulatory level.
[0054] Compared with the prior art, the embodiments of the present invention have the following beneficial effects:
[0055] Based on single-cell sequencing technology, gene expression information of endothelial cells from patients with atherosclerosis is obtained. By extracting and analyzing the features of these original sequencing data, endothelial cells are classified into different subpopulations according to gene expression similarity, and the marker genes and their expression patterns of each subpopulation are determined. Subsequently, by comparing the gene expression differences between different subpopulations, differentially expressed genes that are significantly upregulated or downregulated during the process of atherosclerosis are identified as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells. Compared with the prior art that cannot accurately reveal the complex heterogeneity and dynamic changes within endothelial cells during the process of atherosclerosis, the embodiments of the present invention can accurately analyze endothelial cell heterogeneity through in-depth analysis at single-cell resolution, thereby revealing functional changes and potential molecular markers at the single-cell level, and ultimately achieving a comprehensive understanding of the molecular mechanism of the disease and discovering new therapeutic targets. BRIEF DESCRIPTION OF THE DRAWINGS
[0056] Figure 1 is a schematic flowchart of a method for intelligent analysis of endothelial cell heterogeneity in atherosclerosis provided by an embodiment of the present invention;
[0057] Figure 2 is a schematic structural diagram of a system for intelligent analysis of endothelial cell heterogeneity in atherosclerosis provided by an embodiment of the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS
[0058] The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all of the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.
[0059] See Figure 1 , which is a schematic flowchart of a method for intelligent analysis of endothelial cell heterogeneity in atherosclerosis provided by an embodiment of the present invention. The method for intelligent analysis of endothelial cell heterogeneity in atherosclerosis includes steps S10 to S13:
[0060] S10, obtaining sample data of endothelial cells in the atherosclerosis of the patient; the sample data includes original sequencing data of gene expression information of single endothelial cells, and the original sequencing data is obtained by pre-analyzing the transcriptome level of endothelial cells using single-cell sequencing technology;
[0061] S11, extracting features from the sample data of the endothelial cells to obtain gene expression features of the endothelial cells;
[0062] S12. According to the gene expression characteristics of the endothelial cells, the endothelial cells are divided into different subpopulations according to gene expression similarity to determine the marker genes for distinguishing different subpopulations, and the characteristics of each subpopulation are defined and characterized by analyzing the expression patterns of the marker genes;
[0063] S13. Compare the gene expression differences between different endothelial cell subpopulations, and identify the differentially expressed genes that are significantly up-regulated or down-regulated in the relevant endothelial cell subpopulations during the process of atherosclerosis. The differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells.
[0064] Compared with the prior art, the embodiments of the present invention have the following beneficial effects:
[0065] Based on single-cell sequencing technology, the gene expression information of endothelial cells of patients with atherosclerosis is obtained. Through feature extraction and analysis of these original sequencing data, the endothelial cells are divided into different subpopulations according to gene expression similarity, and the marker genes and their expression patterns of each subpopulation are determined. Subsequently, by comparing the gene expression differences between different subpopulations, the differentially expressed genes that are significantly up-regulated or down-regulated during the process of atherosclerosis are identified as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells. Compared with the prior art that cannot accurately reveal the complex heterogeneity and dynamic changes inside endothelial cells during the process of atherosclerosis, the embodiments of the present invention can accurately analyze endothelial cell heterogeneity through in-depth analysis at the single-cell resolution, so as to reveal the functional changes and potential molecular markers at the single-cell level, and finally achieve a comprehensive understanding of the molecular mechanism of the disease and the discovery of new therapeutic targets.
[0066] As an example, in step S10, the sample acquisition and processing process is as follows:
[0067] Clinical sample collection: Select patients clinically diagnosed with atherosclerosis. During vascular surgery (such as carotid endarterectomy or coronary artery bypass grafting), experienced surgeons use sterile and pyrogen-free delicate surgical instruments to precisely excise small samples of diseased arterial tissue and adjacent relatively normal tissue, and the volume of each sample is about 3-8 cubic millimeters. After the sample is collected, it is immediately placed in pre-cooled physiological saline containing antibiotics (penicillin 100 U / mL and streptomycin 100 μg / mL) and antifungal agent (amphotericin B 0.25 μg / mL), and quickly transported to the cell processing laboratory on ice, and the whole process is controlled within 20 minutes.
[0068] Enzymatic digestion and cell separation: A special composite enzyme solution composed of 0.12% collagenase I, 0.08% trypsin, and 0.03% hyaluronidase was used. The arterial tissue sample was placed in this enzyme solution and gently shaken at a speed of 80 - 100 rpm for enzymatic digestion for 25 - 35 minutes in a constant temperature water bath shaker at 37°C with saturated humidity and 5% CO2.
[0069] After the enzymatic digestion was completed, an equal volume of DMEM / F12 medium containing 15% fetal bovine serum was added to terminate the enzyme reaction. Then, it was filtered through a cell filter with a pore size of 30 μm, and the filtrate was collected. The filtrate was centrifuged at a centrifugal force of 400×g for 8 minutes, the supernatant was discarded, and the precipitate was washed twice with pre-cooled PBS buffer supplemented with 2 mM EDTA. Finally, the cells were resuspended in an appropriate amount of single-cell separation buffer containing 0.04% BSA, and the cell concentration was adjusted to approximately 8×10 5 -1.2×10 6 cells / ml to obtain a high-quality single-cell suspension.
[0070] A single-cell capture system based on laser-induced fluorescence detection technology was used. The single-cell suspension was mixed with reverse transcription primers with unique molecular identifiers (UMIs), fluorescent labels, and reaction reagents using single-cell analysis technology (scRNA-seq). In the microchannels of the microfluidic chip, single cells were identified by laser-induced fluorescence detection, bound to the reverse transcription primers, and then encapsulated by oil droplets. Inside the oil droplets, the cells lysed, and the released mRNA hybridized with the reverse transcription primers and underwent reverse transcription reactions to generate barcoded and fluorescently labeled cDNA.
[0071] After collecting the oil droplets and demulsifying, fluorescence-activated cell sorting technology (FACS) was used to screen and enrich the microbeads with cDNA according to the fluorescent labels, removing impurities and unreacted microbeads. Then, the cDNA was purified and amplified. Using a specially designed ultra-high-precision library preparation kit, the amplified cDNA was first fragmented, with the fragment length controlled at 220 - 480 bp, and then end repair, A-tail addition, adapter ligation, and PCR amplification were carried out in sequence. Finally, a library suitable for ultra-deep high-throughput sequencing was constructed. Sequencing was performed using the Illumina HiSeq X Ten sequencing platform, and a paired-end sequencing mode with a sequencing read length of 180 - 220 bp was set to ensure sufficient gene expression information was obtained.
[0072] Data cleaning and filtering: Using self-written data processing scripts, first remove the reads in which the proportion of bases with a quality value lower than 25 exceeds 15% in the sequencing data, and the reads with a proportion of ambiguous bases (N) greater than 5%. For the problem of adapter sequence contamination, by precisely aligning with a known adapter sequence database, remove the matching adapter sequence part in the reads. Genes with extremely low gene expression levels (total UMI count less than 8 in all cells) are excluded to reduce noise interference in the data.
[0073] As an improvement of the above embodiment, dividing the endothelial cells into different subpopulations according to the gene expression characteristics of the endothelial cells and according to gene expression similarity to determine the marker genes for distinguishing different subpopulations, and defining and characterizing the characteristics of each subpopulation by analyzing the expression patterns of the marker genes, including:
[0074] Similarity measurement and probability calculation steps: Determine the vectors of cells in the high-dimensional gene expression feature space based on the obtained gene expression characteristics, calculate the Gaussian kernel similarity between cells, and calculate the conditional probability and joint probability based on this similarity;
[0075] Low-dimensional vector layout optimization steps: In a low-dimensional space with a preset dimension, set low-dimensional vectors corresponding to the high-dimensional cell vectors, and determine the low-dimensional vector layout by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix and joint probability based on the gene expression feature data, then randomly initialize the low-dimensional vectors so that they are evenly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. Calculate the gradient of the objective function with respect to the low-dimensional vectors at each iteration and update the low-dimensional vectors accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule according to the change of the objective function during the iteration, and the number of iterations is set within a predetermined range and it is judged whether to converge according to the comparison of the change amount of the objective function with the preset convergence threshold;
[0076] Density peak-based clustering steps: After obtaining the optimized low-dimensional vectors, use a density peak-based clustering algorithm in the low-dimensional space to calculate the local density of cells and the minimum distance to cells with higher density. The local density calculation is based on the Euclidean distance between cells and combines the statistical parameters of the distance matrix to determine the calculation scale. For the case where there are no cells with higher density, set the minimum distance to a predetermined large value. Construct a decision graph based on the local density and the minimum distance and select cells with predetermined characteristics as the clustering centers, and assign the remaining cells to the clusters belonging to the nearest clustering centers, thereby realizing the division of endothelial cell subpopulations;
[0077] Marker gene determination step: For each subgroup after division, calculate the change multiple of the average expression level ratio of the genes within the subgroup to the genes of other subgroups. Set a threshold based on this ratio and the average expression level of the genes within the subgroup to screen for marker genes, and use the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
[0078] In this embodiment, focusing on the intelligent analysis of endothelial cell heterogeneity in atherosclerosis, first, through the similarity measurement and probability calculation steps, cell vectors are constructed in a high-dimensional space based on the gene expression characteristics of endothelial cells, and the Gaussian kernel similarity and corresponding probabilities are calculated, laying a foundation for subsequent analysis and aiming to accurately capture the associations between cells. Then, in the low-dimensional vector layout optimization step, corresponding vectors are set in a preset low-dimensional space, and a reasonable layout is determined by minimizing the objective function and through iterative optimization to ensure that the similarity structure of the high-dimensional space can be effectively restored. Subsequently, using the density peak-based clustering algorithm, the local density and relative distance of cells are calculated based on the optimized low-dimensional vectors, and the clustering centers are selected accordingly to achieve the division of endothelial cell subgroups. Finally, in the marker gene determination step, the average expression levels of genes in each subgroup are compared, thresholds are set to screen and verify the marker genes with the help of gene set enrichment analysis, so as to accurately characterize the characteristics of each subgroup. This embodiment can achieve accurate subgroup division of endothelial cells in atherosclerosis according to gene expression similarity, effectively screen out reliable marker genes that can distinguish subgroups, overcome the deficiencies of traditional analysis methods in analyzing endothelial cell heterogeneity, and provide strong data support for in-depth exploration of the disease mechanism of atherosclerosis and subsequent related research.
[0079] For ease of understanding, the working process of this embodiment is specifically described as follows:
[0080] 1. Similarity measurement and probability calculation steps:
[0081] Let the gene expression feature matrix of endothelial cells be F, with dimensions of n×m, where n is the total number of cells, that is, the total number of endothelial cells under study, and m is the gene expression feature dimension, representing the number of different aspects describing the gene expression characteristics of cells. For any two cells i and j, their vectors in the high-dimensional gene expression feature space are respectively denoted as x i and x j ( that is, each vector has m dimensions, and the values of these dimensions correspond to the specific values of the cells in each gene expression feature aspect).
[0082] Calculate the Gaussian kernel similarity S ij between cells, and the formula is Here represents the Euclidean distance between cells i and j in the high-dimensional space, measuring the degree of difference between their gene expression feature vectors, x ik and xjk They are the values of the k-th dimension in the gene expression feature vectors of cells i and j respectively; σ is the scale parameter, ( is the mean of the cell vectors, ), and σ determines the similarity measurement range according to the gene expression distribution of the cell population.
[0083] Based on S ij Calculate the conditional probability which reflects the probability of cell j appearing when cell i is known, and is obtained by normalizing the Gaussian kernel similarity between cells and other cells. S ik represents the Gaussian kernel similarity between cell i and cell k. Then calculate the joint probability p ij : (n is the total number of cells, p i|j is similar to p j|i and is also a conditional probability calculated based on the same Gaussian kernel similarity. The calculation formula is: It reflects the probability of cell i appearing when cell j is known), comprehensively describing the relationship between cells.
[0084] 2. Steps for optimizing the low-dimensional vector layout:
[0085] Set the dimension d of the low-dimensional space, and for the high-dimensional cell vector x i Set the corresponding low-dimensional vector y in the low-dimensional space i Randomly initialize y i to be uniformly distributed in a specific range (such as [-0.03, 0.03]).
[0086] Construct the objective function C: where q ij is the expected cell distribution probability, and the calculation formula is used to compare and measure the difference between the actual and expected cell distributions with p ij ; is related to the vectors of cells i and k in the low-dimensional space, reflecting the influence of their distance relationship pairs in the low-dimensional space; The first term of C measures the difference between the actual and expected cell distributions, the second term is the regularization term to control the smoothness of the low-dimensional space, λ (such as 0.5) is the regularization parameter to adjust its influence degree, and α (such as 0.1), β (such as 0.2) affect q ij and thus act on the objective function; y jl represents the value of the l-th dimension of the low-dimensional vector y j ; y il represents the value of the l-th dimension of the low-dimensional vector y i ; y kl represents the value of the l-th dimension of the low-dimensional vector y i ; yk is the vector in the low-dimensional space corresponding to the high-dimensional cell vector x k where k is used to identify the cell and l is used to identify the dimension number.
[0087] Iteratively optimize C using the gradient descent method. At each iteration t, calculate the gradient and then update y i : is the value of the low-dimensional vector y i at the t-th iteration; is the value of the low-dimensional vector y i at the (t + 1)-th iteration; η (initially set to 0.02) is the learning rate, which decays according to a certain rule (e.g., every 5 iterations, η = η × 0.9), the number of iterations is set to 80 - 100 times, and the convergence threshold is set to 2 × 10 -5 .
[0088] 3. Steps for density peak clustering:
[0089] Use the optimized y i to calculate the local density ρx of the cell: where is the Euclidean distance between cells in the low-dimensional space, construct the distance matrix D (D ij = d ij ), r = μD + 1.3σ D (μ D is the mean of D, σ D is the standard deviation, and d is the dimension of the low-dimensional space), and ρ i reflects the degree of aggregation around cell i.
[0090] Calculate the minimum distance δ between cell i and cells with higher density i : If there are no cells with higher density, δ i is set to the maximum value of D. Based on ρ i and δ i construct a ρ-δ decision graph (constructed using existing decision graph techniques), select cells with large ρ i and δ i as the clustering centers (with ρ i as the abscissa and δ i as the ordinate, select cells with both large ρ i and δ i as the clustering centers), and assign the remaining cells to the corresponding clusters according to the distance from the clustering centers to achieve subpopulation division.
[0091] 4. Steps for determining marker genes:
[0092] For the divided subpopulations, assume the average expression level of gene g in subpopulation A (yg,c is the normalized expression level of gene g in cell c, n A is the number of cells in subgroup A), and the average expression levels of other subgroups are obtained in the same way. Calculate the fold change of the average expression level ratio (μ g,B is the average expression level of gene g in another subgroup B, FC g is used to measure the degree of expression difference of gene g between different subgroups (such as A and B)), and according to FC g and the average gene expression level within the subgroup, set a threshold (such as FC g > 3 and the average expression level is greater than 15) to screen marker genes, and then use gene set enrichment analysis to verify their enrichment in biological function pathways. Specific gene set enrichment analysis can refer to the existing gene set enrichment analysis method (GSEA).
[0093] As an improvement of the above embodiment, comparing the gene expression differences between different endothelial cell subgroups, identifying the differentially expressed genes that are significantly up-regulated or down-regulated in the relevant endothelial cell subgroups during atherosclerosis, and using the differentially expressed genes as key molecular markers to understand the functional changes and heterogeneity of endothelial cells, including:
[0094] Subgroup gene expression level calculation sub-step: Based on the characteristics of the marker genes of different subgroups, for the selected two endothelial cell subgroups, count the number of cells in each subgroup respectively, and then calculate the expression level of the gene in the corresponding subgroup, which is used as the basis data preparation for subsequent model construction and analysis;
[0095] Negative binomial distribution-based model construction and parameter estimation sub-step: Construct a negative binomial distribution model suitable for describing the gene expression level distribution. The model expression contains terms associated with cell characteristic variables, as well as model parameters to be determined and discrete parameters unique to the negative binomial distribution; Use the maximum likelihood estimation method to carry out model parameter estimation work, providing a basis for subsequent gene expression level estimation;
[0096] Support vector machine regression model construction sub-step: Set the training data to consist of cell characteristic vectors and gene expression levels, construct a support vector machine regression model. The objective function of the support vector machine regression model comprehensively considers the relevant norms of the weight vector and the terms related to the slack variables, and at the same time sets a constraint condition system constructed based on cell characteristic vectors, weight vectors, bias terms, slack variables, and insensitive loss parameters; By solving the corresponding quadratic programming problem, determine the parameters of the support vector machine regression model;
[0097] Differential expression analysis and screening sub-step: The gene expression data estimated by the negative binomial distribution model is used as the input information of the support vector machine regression model, and the support vector machine regression model is used to deeply learn and analyze the gene expression difference relationship between cell subsets; based on the adjusted expression data output by the support vector machine regression model, calculate the differential expression fold of the gene between different subsets.
[0098] Testing sub-step: Use a testing method based on the approximation of the normal distribution to perform a significance test on the calculated differential expression fold, and screen out genes with significant differential expression according to the preset significance level threshold.
[0099] In this embodiment, first, through the sub-step of calculating the gene expression of subsets, the cell numbers are counted according to the characteristics of the marker genes, and the average expression of the genes in each subset is calculated, providing basic data for subsequent model construction. Then, in the sub-step of model construction and parameter estimation based on the negative binomial distribution, a model is constructed and parameters are estimated according to the distribution characteristics of gene expression, laying a foundation for the accurate estimation of gene expression. Subsequently, in the sub-step of constructing the support vector machine regression model, a model is constructed based on specific training data, and the model parameters are determined by solving the quadratic programming problem, enhancing the analysis ability of gene expression data. In the subsequent differential expression analysis and screening sub-step, the gene expression data estimated by the negative binomial distribution model is input into the support vector machine regression model to analyze the gene expression difference relationship between cell subsets and calculate the differential expression fold. Finally, in the testing sub-step, a testing method based on the approximation of the normal distribution is used to perform a significance test on the differential expression fold, and genes with significant differential expression are screened out as the key molecular markers for understanding the functional changes and heterogeneity of endothelial cells, providing key information for the research, diagnosis and treatment of atherosclerotic diseases. This embodiment overcomes the deficiencies of traditional methods in identifying differentially expressed genes in endothelial cell subsets, can accurately find the differentially expressed genes with significant changes between endothelial cell subsets in atherosclerosis, and provides data support for deeply exploring the internal molecular mechanism of atherosclerotic diseases, developing new diagnostic methods and optimizing treatment strategies.
[0100] For ease of understanding, the working process of this embodiment is specifically described as follows:
[0101] 1. Sub-step of calculating the gene expression of subsets:
[0102] Let two endothelial cell subsets be M and N. For gene k, use n M to represent the number of cells of gene k in subset M, that is, the total number of cells containing gene k in subset M; n N to represent the number of cells of gene k in subset N. The expression of gene k in cell d of subset M is denoted as x k,d(This expression level, after pre - processing according to the standard, can truthfully reflect the actual expression situation). Then the average expression level of gene k in subgroup M The calculation formula is: Similarly, the average expression level of gene k in subgroup N Through such calculations, the average expression levels of genes in different subgroups can be clearly known, preparing the basic data for subsequent model applications.
[0103] 2. Sub - steps for building a model based on the negative binomial distribution and parameter estimation operator:
[0104] Build a negative binomial distribution model Here, y k,d represents the expression level of gene k in cell d. μ k,d represents the expected expression level of gene k in cell d, and its calculation formula is μ k,d = exp(α 0k + α 1k v d1 + α 2k v d2 +…+ α qk v dq ), where v d1 , v d2 ,…, v dq are the characteristic variables of cell d, such as the specific subgroup where the cell is located, specific metabolic indicators of the cell, etc. They describe the cell characteristics from multiple aspects and affect the calculation of the expected expression level. α 0k , α 1k ,…, α qk are the model parameters. Different combinations of values affect the calculation result of μ k,d and determine the fitting degree of the model to the gene expression level distribution. is a discrete parameter unique to the negative binomial distribution, used to reflect the discrete characteristics of the expression level distribution of gene k. Different genes have different values of this parameter due to their own expression characteristics.
[0105] Build a likelihood function L to carry out parameter estimation. Its expression is: In the formula, m represents the total number of cells participating in the analysis, covering all the endothelial cells of concern. Γ is the gamma function, which is used to handle complex operations involving factorials, integrals, etc., and helps to build a likelihood function that conforms to the gene expression level distribution characteristics. By taking the derivative of L and setting the derivative to 0, the values of α 0k , α 1k ,…, α qk and are determined, and the specific form of the negative binomial distribution model is determined, providing a basis for subsequent gene expression level estimation.
[0106] 3. Sub - steps for building a support vector machine regression model:
[0107] Set the training data as (u i , z i ), where u i is the feature vector of cell i, integrating various information such as the subpopulation to which the cell belongs and indicators related to cell function. Its dimension depends on the number of cell features selected, and the values of each dimension reflect the corresponding feature situation, providing input content for model analysis. z i represents the gene expression level (similar to the y k,d situation mentioned earlier), which is the target variable to be predicted and analyzed by the model.
[0108] Construct a support vector machine regression model, and set its objective function J as: In J, w is the weight vector, whose dimension is associated with the dimension of the cell feature vector u i . Each element corresponds to the weight of different cell features in the model, determining the influence degree of each feature on the result; p represents the number of samples of the training data, that is, the number of cell-gene expression data pairs (u i , z i ) used to train the support vector machine regression model. ||w|| 2 is the square of the norm of the weight vector, which controls the complexity of the model during model training to prevent overfitting. n w is the dimension of the weight vector w, w j is the j-th element of the weight vector w, corresponding to the weight of the j-th feature in the cell feature vector u i in the model, which determines the influence degree of this feature on the model prediction result. C is the penalty parameter, used to adjust the degree of emphasis on the training error. The larger the value, the heavier the penalty for the error, prompting the model to fit the training data better. η i and are slack variables, used to handle the situation where some data points in the actual data are difficult to meet the constraint conditions, and appropriately relax the restrictions to make the model adapt to the real data distribution.
[0109] At the same time, set the constraint condition: z i -w T u i -b ≤ ε + η i , w T is the transpose of the weight vector w.
[0110] Solve this quadratic programming problem constructed based on the objective function and constraint conditions through an existing appropriate optimization algorithm to determine the parameters w, b, etc. of the support vector machine regression model, so that the model can accurately analyze the gene expression level based on the input cell feature vector and achieve an effective mapping from input to output.
[0111] 4. Sub-steps of differential expression analysis and screening:
[0112] Use the gene expression data estimated by the negative binomial distribution model (i.e., the estimated values of μ calculated according to the previously determined model and parameters) as the input features of the support vector machine regression model, and input them into the support vector machine regression model with determined parameters. Based on the input feature information, the model deeply explores and analyzes the gene expression difference relationship between cell subsets, and outputs the adjusted expression data r k,d (r k,d (r k,d is the adjusted estimated value of the expression level of gene k in cell d, which integrates the comprehensive analysis results of the model on cell characteristics and gene expression).
[0113] For the selected two subsets M and N, calculate the differential expression fold change DF of gene k between different subsets (M and N) k , and the calculation formula is n M represents the number of cells of gene k within subset M, that is, the total number of cells containing gene k in subset M; n N represents the number of cells of gene k within subset N. Through this calculation, the difference in the expression levels of genes in different subsets can be intuitively and quantitatively seen, providing a key quantitative index for screening significantly differentially expressed genes, and helping to accurately locate those genes that have important indicative effects on endothelial cell function changes and heterogeneity.
[0114] 5. Inspection sub-step:
[0115] Use the t-test method based on the normal distribution approximation to conduct a significance test on the differential expression fold change DF k . First, calculate the test statistic T according to DF k and the corresponding sample data (data related to the gene expression levels of each subset) (calculated by the conventional mathematical formula of the t-test, integrating relevant data such as the sample mean and sample standard deviation). Then, compare the T value with the preset significance level threshold (for example, set to 0.05, corresponding to the critical value of the corresponding t-distribution), and screen out the genes with a p-value less than this threshold (such as p < 0.05). These genes are the genes with significant differential expression, which are of great significance for deeply understanding the functional changes and heterogeneity characteristics of endothelial cells during the process of atherosclerosis, and can be used as key molecular markers to provide important reference basis for subsequent disease mechanism research, diagnosis and treatment, etc.
[0116] As an improvement of the above embodiment, after comparing the gene expression differences between different endothelial cell subsets, identifying the differentially expressed genes that are significantly up-regulated or down-regulated in the relevant endothelial cell subsets during the process of atherosclerosis, and using the differentially expressed genes as key molecular markers for understanding endothelial cell function changes and heterogeneity, the method further includes:
[0117] Select the protein products corresponding to the significantly differentially expressed genes as potential drug targets, and perform in silico functional validation on the selected potential drug targets.
[0118] Specifically, the selection of the protein products corresponding to the significantly differentially expressed genes as potential drug targets and the in silico functional validation of the selected potential drug targets include:
[0119] Differentially expressed gene screening and ranking sub-step: For the obtained set of differentially expressed genes, calculate the fold change of differential expression and statistical test values of each gene through statistical test methods, and construct a comprehensive scoring function in combination with the gene function importance index. This function importance index is obtained based on a multi-factor comprehensive evaluation of the key degree of the gene participating in biological pathways and the degree of association with known disease-related genes. Rank the differentially expressed genes in descending order according to the comprehensive scoring function;
[0120] Protein product analysis and initial target screening sub-step: Use bioinformatics tools to obtain the information of the protein product set corresponding to the ranked differentially expressed genes. Based on the functional characteristics of the protein, such as its participation in cell signal transduction, interaction with lipid metabolism and inflammatory factors, and the localization characteristics of membrane proteins in subcellular compartments, construct a screening function that includes the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor. Screen out the candidate set of potential drug targets through this screening function in combination with a set threshold. The weight of each screening factor is set according to its relative importance in drug target screening, and the scoring function determines the value according to the matching degree of the protein characteristics and the screening criteria, so as to achieve the initial screening of protein products to determine potential drug target candidates;
[0121] Binding affinity calculation sub-step based on molecular docking: For each potential drug target, select candidate ligand molecules from a specific molecular library. First, perform structural preprocessing on the target protein and ligand molecules, and then perform conformational search and docking in the active pocket region of the target protein. Use a molecular docking algorithm formula that includes the number of atoms of the ligand and target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the parameter controlling the angle range to calculate the binding energy. Determine the partial charges and van der Waals interaction parameters according to the atomic chemical element type, set the electrostatic interaction adjustment parameter according to experience, and set the angle range parameter according to the characteristics of the protein active pocket. Screen out the ligand-target combinations with strong binding affinity by comparing the binding energy with a preset threshold to preliminarily verify the druggability of the target;
[0122] Sub-steps for target stability analysis based on kinetic simulation: For ligand-target combinations with good binding energy, construct a simulation system containing the ligand, target protein, and solvent molecules. Determine the atomic mass parameters according to the atomic chemical composition, set the calculation of the atomic potential energy function and the intermolecular forces according to the molecular mechanics force field parameters, determine the external forces according to the simulation environment temperature and pressure conditions, perform the simulation using the molecular dynamics simulation algorithm formula, set the simulation time and step size, monitor the root mean square deviation (RMSD), an indicator of the target protein structure stability, and judge the structural stability of the target protein after binding to the ligand based on the comparison between the RMSD and the preset range to further verify the effectiveness of the target.
[0123] Sub-steps for verifying the target regulatory relationship based on network analysis: Construct a gene regulatory network containing differentially expressed genes, potential drug targets, and other related genes. The network nodes are genes and targets, and the edges represent regulatory relationships with weight assignments. For each potential drug target, determine the gene set with regulatory relationships with the target and the shortest path length between the gene and the target according to the network topology structure, calculate the importance index of the target in the gene regulatory network using the network analysis algorithm formula, and judge the key degree of the target in the gene regulatory network based on the comparison between this index and the preset threshold, so as to verify the rationality of the target from the gene regulation level.
[0124] In this embodiment, after identifying the differentially expressed genes of endothelial cell subsets during the atherosclerotic process, potential targets for drug development are further explored. First, through the sub-step of differentially expressed gene screening and ranking, considering the degree of differential expression and functional importance of genes comprehensively, the differentially expressed genes are ranked by priority to focus on key genes. Then, in the sub-step of protein product analysis and initial target screening, based on various functional characteristics and localization features of proteins, potential drug target candidates are initially screened using a carefully constructed screening function. Next, through the sub-step of binding affinity calculation based on molecular docking, using complex molecular docking algorithm formulas, the binding ability of ligands to targets is evaluated at the molecular level, and combinations with strong binding affinity are screened to preliminarily verify the druggability of the targets. After that, in the sub-step of target stability analysis based on kinetic simulation, a simulation system is constructed and according to the kinetic simulation algorithm formula, by monitoring the target protein structure stability index, the stability of the target after binding to the ligand is further investigated to enhance the reliability of target verification. Finally, in the sub-step of target regulatory relationship verification based on network analysis, a gene regulatory network is constructed and using network analysis algorithm formulas, the key degree of the target is explored at the gene regulatory level to comprehensively verify the rationality of the target. This embodiment can start from the gene differential expression data of atherosclerotic endothelial cells, efficiently and accurately screen out targets with potential drug development value, and through multi-dimensional and multi-level simulation verification methods, improve the accuracy and reliability of target screening, provide high-quality target information for subsequent drug development against atherosclerosis, and strongly promote the transformation process from basic research to clinical drug development.
[0125] For the sake of easy understanding, the working process of this embodiment is specifically described as follows:
[0126] 1. Sub-step of differentially expressed gene screening and ranking:
[0127] Let the set of obtained differentially expressed genes be G = {g1, g2,..., g n}, for each gene g i , calculate its fold change of differential expression FC i and statistical test value p i through a specific statistical test method (such as t-test or analysis of variance, etc.). At the same time, based on the key degree C(g i ) of the gene participating in biological pathways (the key degree of gene g i participating in biological pathways, which is a factor to evaluate the functional importance of genes, and the larger the value, the more critical the gene is in biological pathways), the association degree A(g i ) with known disease-related genes (gene g iThe degree of association with genes related to known diseases, which is used to measure the correlation between genes and diseases. The higher the value, the closer the connection with genes related to diseases), and other factors are used to comprehensively evaluate the gene function importance index I(g i )(The gene function importance index of gene g i ), for example, I(g i ) = αC(g i ) + βA(g i ), where α and β are weight coefficients preset according to the importance of each factor (assuming α = 0.6 and β = 0.4). Construct the comprehensive scoring function S i of gene g i = log(FC i ) × (-log(p i )) × I(g i ). According to S i , the differentially expressed genes are sorted in descending order to determine those genes that are significantly changed in expression level and may play a key role in function, providing a priority order for subsequent target screening.
[0128] 2. Protein product analysis and initial target screening sub - steps:
[0129] Using bioinformatics tools such as UniProt, obtain the set of protein products P = {p1, p2,..., p m} corresponding to the sorted differentially expressed genes. According to the participation of the protein in cell signal transduction S(p j )(The participation of protein p j in cell signal transduction is a factor for screening potential drug targets; p j is the j - th protein product in the set of protein products P), the interaction with lipid metabolism L(p j )(The interaction of protein p j with lipid metabolism is used to evaluate the association between the protein and lipid metabolism), the interaction with inflammatory factors F(p j )(The interaction of protein p j with inflammatory factors reflects the relationship between the protein and the inflammatory process), and the membrane protein localization feature M(p j )(The membrane protein localization feature of protein p j , such as whether it is a transmembrane protein, etc., is an important basis for screening) and other screening factors in many aspects, construct the screening function where K is the number of screening factors (here K = 4), w k is the weight of the k - th screening factor (for example, w1 = 0.3 corresponds to signal transduction, w2 = 0.25 corresponds to lipid metabolism, w3 = 0.25 corresponds to inflammatory factors, w4 = 0.2 corresponds to membrane protein localization), f k (pj ) is protein p j The scoring function for the k-th screening factor (for example, for membrane protein localization, if it is a transmembrane protein, then f4(p j ) = 1, otherwise f4(p j ) = 0). Set a threshold T (such as T = 0.6), and select the protein products with F(p j ) ≥ T as the candidate set T = {t1, t2, …, t l} of potential drug targets. In this way, proteins with the characteristics of potential drug targets are initially screened out.
[0130] 3. Sub-steps for calculating binding affinity based on molecular docking:
[0131] For each potential drug target t s ∈T, select a series of candidate ligand molecules L = {l1, l2, …, l u} from the known drug molecule library or the virtual designed small molecule compound library. For the target protein t s and the ligand molecule l v , first perform structure preprocessing, including operations such as adding hydrogen atoms and optimizing charge distribution. Then, perform conformational search and docking in the active pocket region of the target protein, and calculate the binding energy E bind using the following molecular docking algorithm formula: Among them, E bind is the binding energy between the ligand and the target protein, N is the number of ligand atoms, M is the number of target protein atoms, q i and q j are the partial charges of ligand atom i and protein atom j respectively, and their values are determined by quantum chemical calculations or empirical databases according to the chemical element types of the atoms (for example, for carbon atoms, the q value may fluctuate within a certain range according to its chemical environment), r ij is the distance between atom i and atom j, A ij and B ij are the van der Waals interaction parameters (depending on the atom type, and different atom pairs have specific empirical values), C is the electrostatic interaction adjustment parameter (set to 0.5 according to experience to balance the contributions of van der Waals forces and electrostatic interactions), θ ij is the angle between the chemical bond vector connecting atom i and atom j and the vector connecting their geometric centers, and σ is the parameter controlling the action range of the angle-related term (set to 0.3 nanometers according to the size and shape characteristics of the protein active pocket). By calculating the binding energies of different ligands and the target, select the ligand-target combinations with lower binding energies (such as E bind < -5.0 kcal / mol), indicating that they have strong binding affinities and initially verifying the druggability of the target.
[0132] 4. Sub - steps for target stability analysis based on kinetic simulation:
[0133] For the ligand - target combinations with good binding energy above, construct a simulation system including the ligand, target protein, and surrounding solvent molecules (such as water molecules). Determine the atomic mass parameter m according to the chemical composition of the atoms i (Based on the relative atomic mass of the atoms in the periodic table and appropriately corrected considering the isotope distribution), and calculate the potential energy function U(r i ) of atom i according to the parameter settings of the molecular mechanics force field (such as AMBER, CHARMM, etc) (including bond energy, angle energy, dihedral angle energy, van der Waals energy, and electrostatic energy, etc. Each energy term is accurately calculated according to the force field parameters and the geometric relationship between atoms) and the intermolecular force F ij (Determined by the calculation formula of the intermolecular interaction in the force field, considering factors such as the distance between atoms, charge, and chemical bond type). Determine the external force F according to the simulated temperature (such as 300K) and pressure (such as 1 atmosphere) conditions ext (Calculated through the ideal gas state equation and thermodynamic principles, which is determined according to the temperature and pressure conditions of the simulation environment. For example, in a simulated physiological environment, temperature and pressure will have a certain effect on molecules, thus affecting the structure and kinetic behavior of molecules). Use the molecular dynamics simulation algorithm formula for simulation: Set the simulation time to 20 - 50 nanoseconds and the time step to 1 - 2 femtoseconds. During the simulation, monitor the structural stability indicators of the target protein, such as the root - mean - square deviation (RMSD): (where r i (t) is the position of atom i at time t, and r i (0) is the initial position). If during the simulation, the RMSD of the target protein remains within a small range (such as RMSD < 0.2 nm), it indicates that the structure of the target protein is stable after binding to the ligand, further verifying the effectiveness of this target.
[0134] 5. Sub - steps for verifying the target regulatory relationship based on network analysis:
[0135] Construct a gene regulatory network including differentially expressed genes, potential drug targets, and other related genes. The nodes in the network are genes and targets, and the edges represent the regulatory relationships between genes or between genes and targets (such as the binding relationship between transcription factors and target genes, the interaction between upstream and downstream molecules in signal pathways, etc.). The weight of the edge can be assigned according to the credibility of experimental evidence or bioinformatics prediction. For each potential drug target t s , determine the gene set with direct or indirect regulatory relationships with the target according to the network topology and the shortest path length d(g, t s)(Calculated by the shortest path algorithm in graph theory, such as Dijkstra's algorithm). Use the network analysis algorithm formula to calculate the importance index I(t) of the target in the gene regulatory network s ) where FC g , p g and I(g) are defined as in the previous differential expression gene screening. According to the comparison between I(t s ) and a preset threshold (such as 0.8), judge the key degree of the target in the gene regulatory network. If I(t s ) is relatively high, it indicates that the target t s is in a key position in the gene regulatory network and has a close regulatory relationship with multiple important differentially expressed genes, verifying the rationality and importance of the target from the gene regulatory level. Based on the above simulation function verification results, determine the most potential drug targets for subsequent drug R & D work.
[0136] See Figure 2 , which is a schematic structural diagram of an intelligent analysis device for endothelial cell heterogeneity in atherosclerosis provided by an embodiment of the present invention. The intelligent analysis system for endothelial cell heterogeneity in atherosclerosis includes:
[0137] An acquisition module 10 for acquiring sample data of endothelial cells in atherosclerosis of a patient; the sample data includes the original sequencing data of the gene expression information of single endothelial cells, and the original sequencing data is obtained by pre-analyzing the transcriptome level of endothelial cells using single-cell sequencing technology;
[0138] A feature extraction module 11 for extracting features from the sample data of the endothelial cells to obtain the gene expression features of the endothelial cells;
[0139] A subpopulation division module 12 for dividing the endothelial cells into different subpopulations according to the gene expression features of the endothelial cells and according to gene expression similarity, to determine the marker genes for distinguishing different subpopulations, and to define and characterize the features of each subpopulation by analyzing the expression patterns of the marker genes;
[0140] A differential analysis module 13 for comparing the gene expression differences between different endothelial cell subpopulations, and identifying the differentially expressed genes that are significantly up-regulated or down-regulated in the relevant endothelial cell subpopulations during the process of atherosclerosis, and the differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells.
[0141] Compared with the prior art, the embodiments of the present invention have the following beneficial effects:
[0142] Based on single-cell sequencing technology, obtain the gene expression information of endothelial cells in patients with atherosclerosis. Through feature extraction and analysis of these original sequencing data, divide endothelial cells into different subpopulations according to gene expression similarity, and determine the marker genes and their expression patterns of each subpopulation. Subsequently, by comparing the gene expression differences between different subpopulations, identify the differentially expressed genes that are significantly upregulated or downregulated during the process of atherosclerosis, and use them as key molecular markers to understand the functional changes and heterogeneity of endothelial cells. Compared with the prior art that cannot accurately reveal the complex heterogeneity and dynamic changes within endothelial cells during the process of atherosclerosis, the embodiments of the present invention can accurately analyze endothelial cell heterogeneity through in-depth analysis at the single-cell resolution, thereby being able to reveal the functional changes and potential molecular markers at the single-cell level, and ultimately achieving a comprehensive understanding of the molecular mechanism of the disease and discovering new therapeutic targets.
[0143] As an improvement of the above solution, the subpopulation division module is specifically used for:
[0144] Determine the vector of cells in the high-dimensional gene expression feature space based on the obtained gene expression features, calculate the Gaussian kernel similarity between cells, and calculate the conditional probability and joint probability based on this similarity;
[0145] In the low-dimensional space of a preset dimension, set low-dimensional vectors corresponding to the high-dimensional cell vectors, and determine the layout of the low-dimensional vectors by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix, and joint probability based on the gene expression feature data, then randomly initialize the low-dimensional vectors so that they are evenly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. Calculate the gradient of the objective function with respect to the low-dimensional vectors at each iteration and update the low-dimensional vectors accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule based on the change of the objective function during the iteration process. The number of iterations is set within a predetermined range and it is judged whether to converge based on the comparison between the change amount of the objective function and the preset convergence threshold;
[0146] After obtaining the optimized low-dimensional vectors, use the density peak-based clustering algorithm in the low-dimensional space to calculate the local density of cells and the minimum distance to cells with higher density. The local density calculation is based on the Euclidean distance between cells and combines the statistical parameters of the distance matrix to determine the calculation scale. For the case where there are no cells with higher density, set the minimum distance to a predetermined large value. Construct a decision graph based on the local density and the minimum distance and select cells with predetermined features as cluster centers, and assign the remaining cells to the clusters belonging to the nearest cluster centers, thereby realizing the division of endothelial cell subpopulations;
[0147] For each subgroup after division, calculate the fold change of the average expression level ratio of the genes within the subgroup to the genes in other subgroups. Set a threshold to screen for marker genes based on this ratio and the average expression level of the genes within the subgroup, and use the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
[0148] As an improvement to the above solution, the differential analysis module is specifically used for:
[0149] Based on the characteristics of the marker genes of different subgroups, for the two selected endothelial cell subgroups, respectively count the number of cells within each subgroup, and then calculate the expression level of the genes in the corresponding subgroup, which serves as the basis data preparation for subsequent model construction and analysis;
[0150] Construct a negative binomial distribution model suitable for describing the gene expression level distribution. The model expression contains terms associated with cell characteristic variables, as well as model parameters to be determined and discrete parameters unique to the negative binomial distribution; Use the maximum likelihood estimation method to carry out model parameter estimation work to provide a basis for subsequent gene expression level estimation;
[0151] Set the training data to consist of cell characteristic vectors and gene expression levels, and construct a support vector machine regression model. The objective function of the support vector machine regression model comprehensively considers the relevant norms of the weight vector and the terms related to the slack variables, and at the same time sets a constraint condition system constructed based on the cell characteristic vector, weight vector, bias term, slack variables, and insensitive loss parameters; By solving the corresponding quadratic programming problem, determine the parameters of the support vector machine regression model;
[0152] Use the gene expression level data estimated by the negative binomial distribution model as the input information of the support vector machine regression model, and use the support vector machine regression model to deeply learn and analyze the gene expression difference relationship between cell subgroups; Calculate the differential expression fold change of genes between different subgroups based on the adjusted expression level data output by the support vector machine regression model;
[0153] Adopt a test method based on the approximation of the normal distribution to perform a significance test on the calculated differential expression fold change, and screen out genes with significantly differential expression according to the preset significance level threshold.
[0154] As an improvement to the above solution, the system further includes:
[0155] A drug target verification module, which is used to select the protein products corresponding to the significantly differentially expressed genes as potential drug targets, and perform simulated functional verification on the selected potential drug targets.
[0156] As an improvement to the above solution, for the obtained differentially expressed gene set, the fold change of differential expression and the statistical test value of each gene are calculated through statistical test methods, and a comprehensive scoring function is constructed by combining the gene function importance index. This function importance index is obtained based on a multi-factor comprehensive evaluation of the key degree of genes participating in biological pathways and the degree of association with known disease-related genes. The differentially expressed genes are sorted in descending order according to the comprehensive scoring function;
[0157] Using bioinformatics tools to obtain the information of the protein product set corresponding to the differentially expressed genes after sorting, based on the functional characteristics of proteins such as their participation in cell signal transduction, interaction with lipid metabolism and inflammatory factors, and subcellular localization characteristics such as membrane protein localization, multiple screening factors are considered. A screening function is constructed that includes the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor. Through this screening function and a set threshold, a candidate set of potential drug targets is screened out. The weight of each screening factor is set according to its relative importance in drug target screening, and the scoring function determines the value according to the matching degree of protein characteristics and screening criteria, so as to achieve the preliminary screening of protein products to determine potential drug target candidates;
[0158] For each potential drug target, candidate ligand molecules are selected from a specific molecular library. First, the structures of the target protein and ligand molecules are preprocessed, and then conformational search and docking are performed in the active pocket region of the target protein. The binding energy is calculated using a molecular docking algorithm formula that includes parameters such as the number of atoms of the ligand and target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the control angle range parameter. The partial charge and van der Waals interaction parameters are determined according to the atomic chemical element type, the electrostatic interaction adjustment parameter is set according to experience, and the angle range parameter is set according to the characteristics of the protein active pocket. Ligand-target combinations with strong binding affinity are screened out by comparing the binding energy with a preset threshold to preliminarily verify the druggability of the target;
[0159] For ligand-target combinations with good binding energy, a simulation system including the ligand, target protein, and solvent molecules is constructed. The atomic mass parameters are determined according to the atomic chemical composition, the atomic potential energy function and interatomic interaction forces are calculated according to the molecular mechanics force field parameters, the external forces are determined according to the temperature and pressure conditions of the simulation environment, and the simulation is carried out using the molecular dynamics simulation algorithm formula. The simulation time and step size are set, and the root mean square deviation, an indicator of the structural stability of the target protein, is monitored. The structural stability of the target protein after binding to the ligand is judged by comparing the root mean square deviation with a preset range to further verify the effectiveness of the target;
[0160] Construct a gene regulatory network that includes differentially expressed genes, potential drug targets, and other related genes. The network nodes are genes and targets, and the edges represent regulatory relationships with weight assignments. For each potential drug target, determine the gene set that has a regulatory relationship with the target and the shortest path length between the gene and the target based on the network topology. Use the network analysis algorithm formula to calculate the importance index of the target in the gene regulatory network. Determine the key degree of the target in the gene regulatory network based on the comparison between this index and a preset threshold, so as to verify the rationality of the target from the gene regulation level.
[0161] It should be noted that the system embodiments described above are merely illustrative. The units described as separate components may or may not be physically separated. The components shown as units may or may not be physical units, that is, they may be located in one place or distributed to multiple network units. Some or all of the modules can be selected according to actual needs to achieve the purpose of the solution of this embodiment. In addition, in the accompanying drawings of the system embodiments provided by the present invention, the connection relationships between the modules indicate that they have communication connections, which can be specifically implemented as one or more communication buses or signal lines. Those of ordinary skill in the art can understand and implement without creative efforts.
[0162] The above is the preferred embodiment of the present invention. It should be pointed out that for those of ordinary skill in the art, without departing from the principle of the present invention, several improvements and refinements can be made, and these improvements and refinements are also regarded as the protection scope of the present invention.
Claims
1. An intelligent analysis method for endothelial cell heterogeneity in atherosclerosis, characterized in that, Including: Obtaining sample data of endothelial cells in atherosclerosis of a patient; the sample data includes raw sequencing data of gene expression information of single endothelial cells, and the raw sequencing data is obtained by pre - analyzing the transcriptome level of endothelial cells using single - cell sequencing technology; Performing feature extraction on the sample data of the endothelial cells to obtain gene expression features of the endothelial cells; According to the gene expression features of the endothelial cells and based on gene expression similarity, dividing the endothelial cells into different sub - populations to determine marker genes for distinguishing different sub - populations, and defining and characterizing the features of each sub - population by analyzing the expression patterns of the marker genes; Comparing the gene expression differences between different endothelial cell sub - populations, and identifying differentially expressed genes that are significantly up - regulated or down - regulated in the process of atherosclerosis in related endothelial cell sub - populations, and the differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells; Among them, the comparing the gene expression differences between different endothelial cell sub - populations, and identifying differentially expressed genes that are significantly up - regulated or down - regulated in the process of atherosclerosis in related endothelial cell sub - populations, and the differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells, includes: Sub - population gene expression quantity calculation sub - step: Based on the characteristics of marker genes of different sub - populations, for two selected endothelial cell sub - populations, respectively count the number of cells in each sub - population, and then calculate the expression quantity of genes in the corresponding sub - populations, which is used as the basic data preparation for subsequent model construction and analysis; Model construction and parameter estimation sub - step based on negative binomial distribution: Construct a negative binomial distribution model suitable for describing the distribution of gene expression quantities. The model expression contains terms associated with cell characteristic variables, model parameters to be determined, and discrete parameters unique to the negative binomial distribution; use the maximum likelihood estimation method to carry out model parameter estimation work to provide a basis for subsequent estimation of gene expression quantities; Support vector machine regression model construction sub - step: Set the training data to consist of cell characteristic vectors and gene expression quantities, construct a support vector machine regression model. The objective function of the support vector machine regression model comprehensively considers the relevant norm of the weight vector and the terms related to slack variables, and at the same time sets a constraint condition system constructed based on cell characteristic vectors, weight vectors, bias terms, slack variables, and insensitive loss parameters; by solving the corresponding quadratic programming problem, determine the parameters of the support vector machine regression model; Differentially expressed analysis and screening sub - step: Use the gene expression quantity data estimated by the negative binomial distribution model as the input information of the support vector machine regression model, and use the support vector machine regression model to deeply learn and analyze the gene expression difference relationship between cell sub - populations; calculate the differential expression fold change of genes between different sub - populations based on the adjusted expression quantity data output by the support vector machine regression model; Testing sub - step: Use a testing method based on normal distribution approximation to perform a significance test on the calculated differential expression fold change, and screen out genes with significant differential expression according to a preset significance level threshold.
2. The intelligent analysis method for endothelial cell heterogeneity in atherosclerosis according to claim 1, wherein, The endothelial cells are divided into different subpopulations according to the gene expression characteristics of the endothelial cells and gene expression similarity, so as to determine marker genes for distinguishing different subpopulations, and the characteristics of each subpopulation are defined and characterized by analyzing the expression patterns of the marker genes, including: Similarity measurement and probability calculation steps: determining the vectors of cells in the high-dimensional gene expression feature space based on the obtained gene expression characteristics, calculating the Gaussian kernel similarity between cells, and calculating the conditional probability and joint probability based on the similarity; Low-dimensional vector layout optimization steps: in a low-dimensional space of a preset dimension, setting low-dimensional vectors corresponding to the high-dimensional cell vectors, and determining the low-dimensional vector layout by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix and joint probability according to the gene expression feature data, then randomly initialize the low-dimensional vectors so that they are evenly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. Calculate the gradient of the objective function with respect to the low-dimensional vectors at each iteration and update the low-dimensional vectors accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule according to the change of the objective function during the iteration, and the number of iterations is set within a predetermined range and convergence is judged according to the comparison between the change amount of the objective function and the preset convergence threshold; Density peak-based clustering steps: after obtaining the optimized low-dimensional vectors, using a density peak-based clustering algorithm in the low-dimensional space, calculating the local density of cells and the minimum distance to higher density cells, calculating the local density based on the Euclidean distance between cells and combining the statistical parameters of the distance matrix to determine the calculation scale, setting the minimum distance to a predetermined large value for the case where there are no higher density cells, constructing a decision graph based on the local density and the minimum distance, and selecting cells with predetermined characteristics as cluster centers, and assigning the remaining cells to the clusters belonging to the nearest cluster centers, thereby realizing the division of endothelial cell subpopulations; Marker gene determination steps: for each subpopulation after division, calculating the change multiple of the average expression level ratio of genes within the subpopulation to genes in other subpopulations, setting a threshold according to the ratio and the average expression level of genes within the subpopulation to screen marker genes, and using the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
3. The intelligent analysis method for endothelial cell heterogeneity in atherosclerosis according to claim 2, characterized in that After comparing the gene expression differences between different endothelial cell subpopulations and identifying the differentially expressed genes that are significantly up-regulated or down-regulated in the process of atherosclerosis in the relevant endothelial cell subpopulations, and using the differentially expressed genes as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells, the method further includes: Selecting the protein products corresponding to the significantly changed differentially expressed genes as potential drug targets, and performing simulated functional verification on the selected potential drug targets.
4. The intelligent analysis method for endothelial cell heterogeneity in atherosclerosis according to claim 3, wherein, The selection of the protein products corresponding to the significantly changed differentially expressed genes as potential drug targets and the performance of simulated functional verification on the selected potential drug targets include: Sub-step for screening and ranking differentially expressed genes: For the obtained set of differentially expressed genes, calculate the fold change of differential expression and statistical test values for each gene through statistical test methods, and construct a comprehensive scoring function in combination with the gene function importance index. This function importance index is obtained based on a multi-factor comprehensive evaluation of the degree of gene participation in key biological pathways and the degree of association with known disease-related genes. Rank the differentially expressed genes in descending order according to the comprehensive scoring function; Sub-step for protein product analysis and preliminary target screening: Use bioinformatics tools to obtain the information of the protein product set corresponding to the ranked differentially expressed genes. Based on multiple screening factors including the participation in cell signal transduction, the interaction with lipid metabolism and inflammatory factors, and the localization characteristics of membrane proteins in subcellular compartments of the protein's functional characteristics, construct a screening function that includes the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor. Screen out the candidate set of potential drug targets through this screening function in combination with a set threshold. The weight of each screening factor is set according to its relative importance in drug target screening, and the scoring function determines the value according to the matching degree between the protein characteristics and the screening criteria, so as to achieve the preliminary screening of protein products to determine potential drug target candidates; Sub-step for calculating binding affinity based on molecular docking: For each potential drug target, select candidate ligand molecules from a specific molecular library. First, perform structural preprocessing on the target protein and ligand molecules, then conduct conformational search and docking in the active pocket region of the target protein. Use the molecular docking algorithm formula that includes the number of atoms of the ligand and target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the control angle range parameters to calculate the binding energy. Determine the partial charges and van der Waals interaction parameters according to the atomic chemical element type, set the electrostatic interaction adjustment parameter according to experience, and set the angle range parameter according to the characteristics of the protein active pocket. Screen out the ligand-target combinations with strong binding affinity by comparing the binding energy with a preset threshold to preliminarily verify the druggability of the target; Sub-step for target stability analysis based on molecular dynamics simulation: For the ligand-target combinations with good binding energy, construct a simulation system that includes the ligand, target protein, and solvent molecules. Determine the atomic mass parameters according to the atomic chemical composition, set the calculation of the atomic potential energy function and the interatomic interaction force according to the molecular mechanics force field parameters, determine the external force according to the simulation environment temperature and pressure conditions, use the molecular dynamics simulation algorithm formula to conduct the simulation, set the simulation time and step size, monitor the root mean square deviation of the target protein structure stability index, and judge the structural stability of the target protein after binding to the ligand according to the comparison between the root mean square deviation and the preset range to further verify the effectiveness of the target; Sub-steps for verifying the target regulation relationship based on network analysis: Construct a gene regulation network containing differentially expressed genes, potential drug targets, and other related genes. The network nodes are genes and targets, and the edges represent regulation relationships with weight assignments. For each potential drug target, determine the gene set with regulation relationships with the target and the shortest path length between the gene and the target according to the network topology structure. Use the network analysis algorithm formula to calculate the importance index of the target in the gene regulation network. Judge the key degree of the target in the gene regulation network based on the comparison between this index and the preset threshold, so as to verify the rationality of the target from the gene regulation level.
5. An intelligent analysis system for endothelial cell heterogeneity in atherosclerosis, characterized in that, Including: An acquisition module for acquiring sample data of endothelial cells in atherosclerosis of a patient; the sample data includes the original sequencing data of the gene expression information of a single endothelial cell, and the original sequencing data is obtained by pre-analyzing the transcriptome level of endothelial cells using single-cell sequencing technology; A feature extraction module for extracting features from the sample data of the endothelial cells to obtain the gene expression features of the endothelial cells; A subpopulation division module for dividing the endothelial cells into different subpopulations according to the gene expression features of the endothelial cells and according to gene expression similarity, to determine the marker genes for distinguishing different subpopulations, and to define and characterize the features of each subpopulation by analyzing the expression patterns of the marker genes; A differential analysis module for comparing the gene expression differences between different endothelial cell subpopulations, and identifying the differentially expressed genes that are significantly up-regulated or down-regulated in the process of atherosclerosis in the relevant endothelial cell subpopulations. The differentially expressed genes are used as key molecular markers for understanding the functional changes and heterogeneity of endothelial cells; Among them, the differential analysis module is specifically used for: Based on the characteristics of the marker genes of different subpopulations, for two selected endothelial cell subpopulations, respectively count the number of cells in each subpopulation, and then calculate the expression level of the gene in the corresponding subpopulation, which is used as the basic data preparation for subsequent model construction and analysis; Construct a negative binomial distribution model suitable for describing the gene expression level distribution. The model expression contains terms associated with cell characteristic variables, as well as model parameters to be determined and discrete parameters unique to the negative binomial distribution; use the maximum likelihood estimation method to carry out model parameter estimation work to provide a basis for subsequent estimation of gene expression levels; Set the training data to consist of cell characteristic vectors and gene expression levels, and construct a support vector machine regression model. The objective function of the support vector machine regression model comprehensively considers the relevant norms of the weight vector and the terms related to the slack variables, and at the same time sets a constraint condition system constructed based on cell characteristic vectors, weight vectors, bias terms, slack variables, and insensitive loss parameters; by solving the corresponding quadratic programming problem, determine the parameters of the support vector machine regression model; Use the gene expression level data estimated by the negative binomial distribution model as the input information of the support vector machine regression model, and use the support vector machine regression model to deeply learn and analyze the gene expression difference relationship between cell subpopulations; calculate the differential expression fold of the gene between different subpopulations based on the adjusted expression level data output by the support vector machine regression model; A test method based on the approximation of the normal distribution is used to perform a significance test on the calculated fold change of differential expression, and genes with significant differential expression are screened according to a preset significance level threshold.
6. The intelligent analysis system for endothelial cell heterogeneity in atherosclerosis according to claim 5, wherein The subpopulation division module is specifically used for: Determining the vectors of cells in the high-dimensional gene expression feature space based on the obtained gene expression features, calculating the Gaussian kernel similarity between cells, and calculating the conditional probability and joint probability based on this similarity; In the low-dimensional space of a preset dimension, low-dimensional vectors are set corresponding to the high-dimensional cell vectors, and the layout of the low-dimensional vectors is determined by minimizing a specific objective function; during the optimization process, first calculate the Gaussian kernel similarity matrix, conditional probability matrix and joint probability according to the gene expression feature data, then randomly initialize the low-dimensional vectors so that they are evenly distributed in the low-dimensional space and within a specific numerical range, and then iteratively optimize the objective function by the gradient descent method. Calculate the gradient of the objective function with respect to the low-dimensional vectors at each iteration and update the low-dimensional vectors accordingly, where the learning rate is initially set to a specific value and decays according to a predetermined rule according to the change of the objective function during the iteration process, and the number of iterations is set within a predetermined range and whether to converge is judged by comparing the change amount of the objective function with the preset convergence threshold; After obtaining the optimized low-dimensional vectors, use the density peak-based clustering algorithm in the low-dimensional space to calculate the local density of cells and the minimum distance to cells with higher density. The local density calculation is based on the Euclidean distance between cells and combines the statistical parameters of the distance matrix to determine the calculation scale. For the case where there are no cells with higher density, the minimum distance is set to a predetermined large value. Construct a decision graph based on the local density and the minimum distance and select cells with predetermined features as cluster centers, and assign the remaining cells to the clusters belonging to the nearest cluster centers, so as to realize the division of endothelial cell subpopulations; For each subpopulation after division, calculate the change multiple of the average expression quantity ratio of the genes within the subpopulation to the genes in other subpopulations, set a threshold to screen marker genes according to this ratio and the average expression quantity of the genes within the subpopulation, and use the gene set enrichment analysis method to perform biological function pathway enrichment analysis on the marker genes for verification.
7. The intelligent analysis system for endothelial cell heterogeneity in atherosclerosis according to claim 6, wherein It further includes: A drug target verification module, which is used to select the protein products corresponding to the significantly differentially expressed genes as potential drug targets and perform simulation function verification on the selected potential drug targets.
8. The intelligent analysis system for endothelial cell heterogeneity in atherosclerosis according to claim 7, wherein The drug target verification module is specifically used for: For the obtained set of differentially expressed genes, calculate the change multiple of the differential expression of each gene and the statistical test value through a statistical test method, and construct a comprehensive scoring function in combination with the gene function importance index. This function importance index is based on a multi-factor comprehensive evaluation of the degree of participation of genes in biological pathways and the degree of association with known disease-related genes. Sort the differentially expressed genes in descending order according to the comprehensive scoring function; Use bioinformatics tools to obtain the information of the set of protein products corresponding to the differentially expressed genes after sorting. Based on the functional characteristics of the proteins, including their participation in cell signaling transduction, interaction with lipid metabolism and inflammatory factors, and subcellular localization characteristics such as membrane protein localization, a screening function is constructed that includes the number of screening factors, the weight of each factor, and the scoring function of the protein for each factor. Through this screening function and in combination with a set threshold, a candidate set of potential drug targets is screened out. Among them, the weights of each screening factor are set according to their relative importance in drug target screening, and the scoring function determines the value according to the matching degree between the protein characteristics and the screening criteria, so as to achieve the preliminary screening of protein products to determine the candidates of potential drug targets; For each potential drug target, select candidate ligand molecules from a specific molecular library. First, perform structural preprocessing on the target protein and the ligand molecules, and then conduct conformational search and docking in the active pocket region of the target protein. Use the molecular docking algorithm formula that includes the number of atoms of the ligand and the target protein, atomic partial charges, interatomic distances, van der Waals interaction parameters, electrostatic interaction adjustment parameters, the angle between the interatomic chemical bond vector and the geometric center vector, and the parameter for controlling the angle range to calculate the binding energy. Determine the partial charges and van der Waals interaction parameters according to the atomic chemical element type, set the electrostatic interaction adjustment parameter according to experience, and set the angle range parameter according to the characteristics of the protein active pocket. By comparing the binding energy with a preset threshold, select the ligand-target combinations with strong binding affinity to preliminarily verify the druggability of the target; For the ligand-target combinations with good binding energy, construct a simulation system that includes the ligand, the target protein, and solvent molecules. Determine the atomic mass parameters according to the atomic chemical composition, set the calculation of the atomic potential energy function and the interatomic interaction force according to the molecular mechanics force field parameters, determine the external force according to the simulation environment temperature and pressure conditions, and use the molecular dynamics simulation algorithm formula to perform the simulation. Set the simulation time and step size, and monitor the root mean square deviation of the target protein structure stability index. Judge the structural stability of the target protein after binding to the ligand according to the comparison between the root mean square deviation and the preset range to further verify the effectiveness of the target; Construct a gene regulatory network that includes differentially expressed genes, potential drug targets, and other related genes. The network nodes are genes and targets, and the edges represent regulatory relationships and have weight assignments. For each potential drug target, determine the gene set with regulatory relationships with the target and the shortest path length between the gene and the target according to the network topology structure. Use the network analysis algorithm formula to calculate the importance index of the target in the gene regulatory network. Judge the key degree of the target in the gene regulatory network according to the comparison between this index and the preset threshold, so as to verify the rationality of the target from the gene regulatory level.
Citation Information
Patent Citations
Analysis method based on 10X unicell transcriptome sequencing data
CN109979538A