Full-path signal flow analysis method, device, equipment and storage medium

Abstract

Embodiments of the present application provide a full-path signal flow analysis method, device, equipment and storage medium, and relate to the technical field of bioinformatics. A transcription factor activity matrix is determined according to a single-cell transcriptome expression matrix and prior regulatory network data; at least one biological pathway to be verified is determined from the single-cell transcriptome expression matrix based on a prior pathway database, the biological pathway to be verified is traversed, and a pathway activation parameter of the biological pathway to be verified in a single cell is determined based on the single-cell transcriptome expression matrix, the prior regulatory network data and the transcription factor activity matrix; and the biological pathway to be verified is subjected to significance verification based on the pathway activation parameter, so as to determine a target signal flow pathway. The biological process of signal from extracellular recognition to intracellular transcription regulation is restored through a complete signal transduction link, and the false positive probability is reduced. The cell sample analysis in a fluid environment such as blood can be adapted, and the integrity, accuracy and scene adaptability of the cell communication analysis result are improved.

Description

Technical Field

[0001] This application relates to the field of bioinformatics technology, and in particular to methods, apparatus, devices and storage media for full-pathway signal flow analysis. Background Technology

[0002] In the field of life science research, cell communication analysis at the single-cell level has become an important approach for discovering key biological mechanisms. For example, in tumor immunology, communication analysis is widely used to discover how cancer cells reprogram infiltrating immune cells by secreting specific immune checkpoint ligands, thereby achieving immune evasion.

[0003] Related technologies typically utilize cell communication analysis models to calculate communication probabilities based on the expression levels of ligands in sending cell clusters and receptors in receiving cell clusters using statistical methods. Alternatively, inference models based on ligand-target gene associations can be used to establish associations from extracellular ligands to intracellular target genes, predicting genes that exhibit differential expression driven by signals. However, cell communication analysis models focus solely on the ligand-receptor relationship between cells, failing to verify the regulation of signals within recipient cells. Furthermore, these models often rely on static tissue assumptions, neglecting the kinetics and nonlinear responses of ligand diffusion in fluid environments such as blood. Inference models, which predict cell communication directly from the originating ligand to the recipient cell's target gene, lack quantitative characterization of the activity of key intermediate transcription factors, resulting in an incomplete intracellular signaling pathway. Summary of the Invention

[0004] The main objective of this application is to propose a method, apparatus, device, and storage medium for full-pathway signal flow analysis, thereby improving the completeness of cellular pathway signal flow analysis.

[0005] To achieve the above objectives, a first aspect of this application proposes a method for analyzing signal flow across an entire path, comprising: The single-cell transcriptome expression matrix and prior regulatory network data of the sample are obtained, and the transcription factor activity matrix is ​​determined based on the single-cell transcriptome expression matrix and the prior regulatory network data. Based on a priori pathway database, at least one experimental biological pathway is identified from the single-cell transcriptome expression matrix, and the experimental biological pathway includes at least a ligand, a receptor, a transcription factor, and a target gene. The biological pathways to be tested are traversed, and based on the single-cell transcriptome expression matrix, the prior regulatory network data, and the transcription factor activity matrix, the pathway activation parameters of the biological pathways to be tested in single cells are determined; the pathway activation parameters include the signal flow probability intensity from the receptor to the transcription factor and the signal regulation intensity from the transcription factor to the target gene. Based on the signal modulation intensity and the signal flow probability intensity, the biological pathway to be tested is validated for significance in order to determine the target signal flow pathway.

[0006] In some embodiments, the single-cell transcriptome expression matrix includes the expression levels of different genes in each cell; the prior regulatory network data includes the regulatory weights between different genes and each transcription factor; determining the transcription factor activity matrix based on the single-cell transcriptome expression matrix and the prior regulatory network data includes: Obtain the intersection genes between the single-cell transcriptome expression matrix and the prior regulatory network data, construct an intersection expression matrix based on the expression levels of the intersection genes, and construct an intersection weight matrix based on the regulatory weights of the intersection genes. The transcription factor activity matrix is ​​calculated based on the transpose of the intersection expression matrix and the intersection weight matrix; the transcription factor activity matrix contains the activity value of each transcription factor in different cells.

[0007] In some embodiments, the determination of at least one experimental biological pathway from the single-cell transcriptome expression matrix based on a priori pathway database includes: Based on a priori pathway database, at least one initial biological pathway is obtained, wherein the initial biological pathway includes at least a ligand, a receptor, a transcription factor, and a target gene; From the single-cell transcriptome expression matrix, obtain the expression levels of the ligands, receptors, transcription factors, and target genes in the initial biological pathway in different cells; Determine whether the proportion of cells with non-zero expression levels in the ligand, receptor, transcription factor, and target gene is greater than or equal to a preset threshold. If so, the initial biological pathway is taken as the biological pathway to be tested.

[0008] In some embodiments, determining the pathway activation parameters of the biological pathway to be tested in a single cell based on the single-cell transcriptome expression matrix, the prior regulatory network data, and the transcription factor activity matrix includes: Based on the single-cell transcriptome expression matrix, the mean expression levels of the receptor, transcription factor, and target gene in each cell were determined, and the receptor half-saturation coefficient, transcription factor half-saturation coefficient, and target gene half-saturation coefficient were obtained. For the receptor and the transcription factor, the probability intensity of the signal flow is determined based on the single-cell transcriptome expression matrix, the prior regulatory network data, the receptor half-saturation coefficient, the transcription factor half-saturation coefficient, and the transcription factor activity matrix; The signal regulation intensity is determined for the transcription factor and the target gene based on the single-cell transcriptome expression matrix, the prior regulatory network data, the half-saturation coefficient of the transcription factor, the half-saturation coefficient of the target gene, and the transcription factor activity matrix.

[0009] In some embodiments, determining the signal flow probability intensity based on the single-cell transcriptome expression matrix, the prior regulatory network data, the receptor half-saturation coefficient, the transcription factor half-saturation coefficient, and the transcription factor activity matrix includes: Obtain the receptor expression level corresponding to the receptor in the single-cell transcriptome expression matrix, and the transcription factor activity value corresponding to the transcription factor in the transcription factor activity matrix; The receptor activation percentage is obtained by calculating the proportion of the receptor expression level in the sum of the receptor expression level and the receptor half-saturation coefficient. The transcription factor activation percentage is obtained by calculating the proportion of the transcription factor activity value in the sum of the transcription factor activity value and the transcription factor half-saturation coefficient. The physical activation probability of the receptor activating the transcription factor is obtained based on the receptor activation percentage and the transcription factor activation percentage. The signal flow probability intensity is obtained based on the physical activation probability and the regulatory weights of the transcription factors in the prior regulatory network data.

[0010] In some embodiments, determining the signal regulation intensity based on the single-cell transcriptome expression matrix, the prior regulatory network data, the transcription factor half-saturation coefficient, the target gene half-saturation coefficient, and the transcription factor activity matrix includes: Obtain the expression levels of the target genes corresponding to the target genes in the single-cell transcriptome expression matrix and the activity values ​​of the transcription factors corresponding to the transcription factors in the transcription factor activity matrix; The second-order activity value is obtained from the transcription factor activity value, and the second-order half-saturation coefficient is obtained from the transcription factor half-saturation coefficient. The target gene intensity ratio is obtained by calculating the proportion of the target gene expression level in the sum of the target gene expression level and the target gene half-saturation coefficient; and the transcription factor intensity ratio is obtained by calculating the proportion of the second-order activity value in the sum of the second-order activity value and the second-order half-saturation coefficient. The physical response intensity of the transcription factor activating the target gene is obtained based on the target gene intensity ratio and the transcription factor intensity ratio. The signal regulation intensity is obtained based on the physical response intensity and the regulatory weights corresponding to the transcription factors in the prior regulatory network data.

[0011] In some embodiments, the significance verification of the biological pathway under test based on the signal modulation intensity and the signal flow probability intensity includes: The pathway activation parameters of a single cell are obtained by multiplying the signal modulation intensity and the signal flow probability intensity. The true observed signal intensity of the biological pathway under test is obtained by calculating the average value of the pathway activation parameters in all cells within a specific cell population. Randomly shuffle the pairings of the receptor and the transcription factor or the transcription factor and the target gene in the biological pathway to be tested to obtain the corresponding randomized biological pathway; For each of the random biological pathways, a corresponding virtual signal intensity is calculated. The significance P-value of the biological pathway to be tested is calculated based on the actual observed signal intensity and the virtual signal intensity. The significance verification result is obtained based on the significance P-value.

[0012] To achieve the above objectives, a second aspect of this application provides an omnichannel signal flow analysis device, comprising: Activity matrix calculation module: used to obtain the single-cell transcriptome expression matrix and prior regulatory network data of the sample, and to determine the transcription factor activity matrix based on the single-cell transcriptome expression matrix and the prior regulatory network data; Initial Pathway Filtering Module: Used to identify at least one experimental biological pathway from the single-cell transcriptome expression matrix based on a priori pathway database, wherein the experimental biological pathway includes at least a ligand, receptor, transcription factor, and target gene; Quantitative calculation module: used to traverse the biological pathway to be tested, and determine the pathway activation parameters of the biological pathway to be tested in a single cell based on the single-cell transcriptome expression matrix, the prior regulatory network data, and the transcription factor activity matrix; the pathway activation parameters include the signal flow probability intensity from the receptor to the transcription factor and the signal regulation intensity from the transcription factor to the target gene; Verification output module: used to perform saliency verification on the biological pathway to be tested based on the signal modulation intensity and the signal flow probability intensity, so as to determine the target signal flow pathway. To achieve the above objective, a third aspect of this application provides an electronic device, the electronic device including a memory and a processor, the memory storing a computer program, and the processor executing the computer program to implement the method described in the first aspect.

[0013] To achieve the above objectives, a fourth aspect of the present application provides a storage medium storing a computer program that, when executed by a processor, implements the method described in the first aspect.

[0014] The method, apparatus, device, and storage medium for full-pathway signal flow analysis proposed in this application acquire single-cell transcriptome expression matrices and prior regulatory network data of samples, and determine transcription factor activity matrices based on these data. Based on a prior pathway database, at least one experimental biological pathway is identified from the single-cell transcriptome expression matrix. This experimental biological pathway includes at least ligands, receptors, transcription factors, and target genes. The experimental biological pathway is traversed, and based on the single-cell transcriptome expression matrix, prior regulatory network data, and transcription factor activity matrix, pathway activation parameters of the experimental biological pathway in single cells are determined. These activation parameters include the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene. Based on the signal regulation intensity and signal flow probability intensity, the experimental biological pathway is significance-verified to determine the target signal flow pathway. This application constructs a complete signal transduction pathway including receptors, transcription factors, and target genes by quantitatively calculating transcription factor activity, supplementing the intermediate links in intracellular signal transduction, and restoring the entire process from extracellular recognition to intracellular transcriptional regulation. The resulting analysis results are more consistent with real biological processes, improving the completeness of cell communication analysis pathways. Simultaneously, by quantitatively calculating the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene, layer-by-layer verification of signal transduction is achieved. This effectively eliminates false-positive communication events where only receptor-ligand co-expression exists but downstream transcriptional regulation is not actually triggered, improving the accuracy of communication inference. Furthermore, quantitative calculations of signal intensity based on transcriptome expression characteristics do not rely on fixed cell spatial location information, breaking free from the constraints of static tissue spatial assumptions. This allows for adaptation to the physiological characteristics of dynamic cell flow and dispersed ligand distribution in fluid environments such as blood and lymph, enhancing the rationality and applicability of cell sample communication analysis results in fluid environments. Finally, after statistical significance verification, the output signal flow pathways possess statistical significance support, improving the completeness, accuracy, and scenario adaptability of cell communication analysis results. Attached Figure Description

[0015] Figure 1 This is a flowchart of the full-path signal flow analysis method provided in the embodiments of this application.

[0016] Figure 2 This is a flowchart provided in the embodiments of this application for determining the transcription factor activity matrix based on the single-cell transcriptome expression matrix and prior regulatory network data.

[0017] Figure 3 This is a flowchart provided in this application embodiment for determining at least one experimental biological pathway from a single-cell transcriptome expression matrix based on a priori pathway database.

[0018] Figure 4This is a flowchart provided in this application embodiment for determining the pathway activation parameters of a biological pathway to be tested in a single cell based on a single-cell transcriptome expression matrix, prior regulatory network data, and transcription factor activity matrix.

[0019] Figure 5 This is a flowchart provided in the embodiments of this application for determining the probability intensity of signal flow based on single-cell transcriptome expression matrix, prior regulatory network data, receptor half-saturation coefficient, transcription factor half-saturation coefficient, and transcription factor activity matrix.

[0020] Figure 6 This is a flowchart provided in the embodiments of this application for determining the intensity of signal regulation based on single-cell transcriptome expression matrix, prior regulatory network data, transcription factor half-saturation coefficient, target gene half-saturation coefficient, and transcription factor activity matrix.

[0021] Figure 7 This is a schematic diagram comparing the nonlinear response characteristics of the first-order Hill function and the second-order Hill function provided in the embodiments of this application.

[0022] Figure 8 This is a flowchart provided in this application embodiment for significant verification of the biological pathway under test based on signal modulation intensity and signal flow probability intensity.

[0023] Figure 9 This is a schematic diagram showing the comparison of verification results provided in the embodiments of this application.

[0024] Figure 10 This is a diagram illustrating the visualization results.

[0025] Figure 11 This is a schematic diagram of the cell signaling flow results.

[0026] Figure 12 This is a structural block diagram of an all-path signal flow analysis device provided in another embodiment of this application.

[0027] Figure 13 This is a schematic diagram of the hardware structure of the electronic device provided in the embodiments of this application. Detailed Implementation

[0028] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the scope of this application.

[0029] It should be noted that although functional modules are divided in the device schematic diagram and the logical order is shown in the flowchart, in some cases, the steps shown or described may be performed in a different order than the module division in the device or the order in the flowchart.

[0030] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application belongs. The terminology used herein is for the purpose of describing embodiments of this application only and is not intended to limit this application.

[0031] With the rapid development of single-cell transcriptome sequencing technology, researchers' understanding of cellular heterogeneity has shifted from static cell typing to dynamic functional synergy, and the analysis of intercellular communication has become crucial for revealing the homeostasis and disease evolution of complex biological systems. In the field of life science research, single-cell-level cell communication analysis has become an important approach to discovering key biological mechanisms. For example, in tumor immunology, communication analysis is widely used to discover how cancer cells reprogram infiltrating immune cells by secreting specific immune checkpoint ligands (such as PD-L1), thereby achieving immune escape. In developmental biology, researchers use communication analysis to identify key signal gradients driving stem cell differentiation into specific lineages and to discover the precise spatial and temporal synergy of core microenvironmental signals that maintain tissue homeostasis, such as the Wnt and Notch pathways.

[0032] Cell communication analysis tools in related technologies mainly rely on the co-expression frequency of ligands and receptors in different cell clusters, focusing on answering the initial question of "which cells might be communicating." However, as research progressed, researchers discovered that simple ligand-receptor pair matching does not equate to the activation of biological functions. Due to the sparsity of single-cell sequencing data and strong background noise, many predicted communication pathways proved to be false positives in biological validation. This has led to an urgent technical need in the field: not only to identify the sending and receiving of signals, but also to quantify whether the signals, after reaching the receptor, actually cross the cell membrane, enter the intracellular regulatory network, and ultimately drive the transcriptional expression of the target gene. Based on this need, researchers further developed cell communication software such as NicheNet, but this only links ligands to genes within the cell.

[0033] Furthermore, most communication models in related technologies implicitly assume "relatively fixed organizational structure" or "spatial contact dependence," which shows a significant mismatch when analyzing the highly dynamic blood system. In solid tissues, cell communication is highly dependent on direct physical contact, such as protein anchoring or synaptic interactions. However, in the blood microenvironment, the spatial location of immune cells is constantly in flux, and their communication patterns largely follow paracrine and endocrine mechanisms. In this dynamic physical context, a large number of free signaling molecules, such as interleukins, chemokines, and metabolites, are distributed globally with the blood flow. Immune cells capture these mobile ligands through surface receptors or internalize them via endocytosis, thereby converting extracellular physical signals into intracellular chemical signal cascades. This process further triggers complex signal transduction pathways, driving the nuclear translocation, allosteric changes, and activity alterations of specific transcription factors, and ultimately determining the functional polarization and phenotypic evolution of immune cells such as monocytes through precise regulation of target gene sets.

[0034] Furthermore, research on cell communication in blood is crucial for establishing the correlation between macroscopic blood test indicators and microscopic changes in the PBMC transcriptome in clinical studies. Current clinical biochemical and protein assays mostly focus on describing systemic results, making it difficult to reveal the decision-making mechanisms of immune cell responses to signal stimulation. However, by using cross-scale modeling to couple plasma signal concentration with intracellular signal flux intensity, we can not only uncover clinical molecular biomarkers earlier than traditional indicators, improving the sensitivity of disease diagnosis, but also predict individualized differences in the responses of different patients to systemic signals.

[0035] As mentioned above, CellChat, a cell communication analysis software in related technologies, is mainly based on a statistical model of ligand-receptor co-expression frequency. It integrates a known database of ligand-receptor interactions and calculates the communication probability using statistical methods such as the law of mass action or permutation tests, based on the expression levels of the ligand in the sending cell cluster and the receptor in the receiving cell cluster. In contrast, NicheNet, an inference model based on ligand-target gene associations, integrates signal transduction networks and gene regulatory networks as prior knowledge to infer which ligands may affect observed gene expression changes in target cells, establishing associations from extracellular ligands to intracellular target genes, and predicting genes that undergo differential expression driven by signals.

[0036] However, CellChat's focus on cell communication is limited to ligand-receptor relationships between cells, failing to verify the regulation of signals within recipient cells. Furthermore, its model is largely based on static tissue assumptions, neglecting the kinetics and nonlinear responses of ligand diffusion in fluid environments such as blood, making it difficult to simulate fluid-based signal capture and kinetic response processes. NicheNet's proposed cell communication, which runs directly from the originating ligand to the recipient cell's target gene, lacks quantitative characterization of the activity of key intermediate transcription factors, resulting in an incomplete intracellular signaling pathway.

[0037] Based on this, embodiments of this application provide a method, apparatus, device, and storage medium for full-pathway signal flow analysis. By quantitatively calculating transcription factor activity, a complete signal transduction pathway including receptors, transcription factors, and target genes is constructed, supplementing the intermediate links of intracellular signal transduction and restoring the entire process from extracellular recognition to intracellular transcriptional regulation. The resulting analysis results are more consistent with real biological processes, improving the completeness of the cell communication analysis pathway. Simultaneously, by quantitatively calculating the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene, layer-by-layer verification of signal transduction is achieved, effectively eliminating false-positive communication events where only receptor-ligand co-expression exists but downstream transcriptional regulation is not actually triggered, thus improving the accuracy of communication inference. Furthermore, quantitative calculations related to signal intensity are performed based on transcriptome expression characteristics, eliminating the need for fixed cell spatial location information and overcoming the constraints of static tissue spatial assumptions. This allows for adaptation to the physiological characteristics of dynamic cell flow and dispersed ligand distribution in fluid environments such as blood and lymph, improving the rationality and applicability of cell sample communication analysis results in fluid environments. Finally, after statistical significance verification, the output signal flow pathway has statistical significance support, improving the completeness, accuracy, and scenario adaptability of cell communication analysis results.

[0038] This application provides a method, apparatus, device, and storage medium for full-path signal flow analysis, which are specifically described through the following embodiments. First, the full-path signal flow analysis method in this application embodiment is described.

[0039] The omnichannel signal flow analysis method provided in this application relates to the field of bioinformatics technology. This method can be applied to a terminal, a server, or a computer program running on either a terminal or a server. For example, the computer program can be a native program or software module in an operating system; it can be a native application (APP), i.e., a program that needs to be installed in the operating system to run, such as a client supporting omnichannel signal flow analysis, i.e., a program that only needs to be downloaded to a browser environment to run; or it can be a small program that can be embedded in any APP. In short, the above-mentioned computer program can be any form of application, module, or plugin. The terminal communicates with the server via a network. The omnichannel signal flow analysis method can be executed by the terminal or the server, or by the terminal and the server working together.

[0040] In some embodiments, the terminal can be a smartphone, tablet, laptop, desktop computer, or smartwatch, etc. The server can be a standalone server, or a cloud server providing basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communication, middleware services, domain name services, security services, content delivery networks (CDNs), and big data and artificial intelligence platforms; it can also be a service node in a blockchain system, where the service nodes form a peer-to-peer (P2P) network. The P2P protocol is an application layer protocol running on top of the Transmission Control Protocol (TCP). The terminal and server can connect via Bluetooth, Universal Serial Bus (USB), or a network, etc., and this embodiment does not impose any limitations.

[0041] This application can be used in a wide variety of general-purpose or special-purpose computer system environments or configurations. Examples include: personal computers, server computers, handheld or portable devices, tablet devices, multiprocessor systems, microprocessor-based systems, set-top boxes, programmable consumer electronics, network PCs, minicomputers, mainframe computers, and distributed computing environments including any of the above systems or devices. This application can be described in the general context of computer-executable instructions executed by a computer, such as program modules. Generally, program modules include routines, programs, objects, components, data structures, etc., that perform specific tasks or implement specific abstract data types. This application can also be practiced in distributed computing environments where tasks are performed by remote processing devices connected via a communication network. In distributed computing environments, program modules can reside in local and remote computer storage media, including storage devices.

[0042] The following describes the omnichannel signal flow analysis method in the embodiments of this application.

[0043] Figure 1 This is an optional flowchart of the full-path signal flow analysis method provided in the embodiments of this application. Figure 1 The method may include, but is not limited to, steps 110 to 140. It is also understood that this embodiment... Figure 1 The order of steps 110 to 140 is not specifically limited. The order of steps can be adjusted or some steps can be reduced or added according to actual needs.

[0044] Step 110: Obtain the single-cell transcriptome expression matrix and prior regulatory network data of the sample, and determine the transcription factor activity matrix based on the single-cell transcriptome expression matrix and prior regulatory network data.

[0045] In one embodiment, a single-cell transcriptome expression matrix can be obtained by performing a standard preprocessing procedure using single-cell transcriptome data from the sample. For example, the standard preprocessing procedure could be the Seurat V5 standard procedure. After performing quality control, standardization, and cell annotation based on this procedure, the single-cell transcriptome expression matrix is ​​extracted. Here, the single-cell transcriptome expression matrix is ​​a two-dimensional matrix of genes × cells. For example, rows correspond to different genes, columns correspond to different cells, and the values ​​in the matrix represent the standardized expression levels of each gene in the corresponding cells. The expression levels are used to reflect the transcriptional abundance information of each gene in a single cell.

[0046] Next, prior regulatory network data is obtained based on publicly available prior regulatory databases, such as the NicheNet regulatory network dataset, as a priori reference for signaling pathway inference. The prior regulatory network data employs a two-dimensional matrix structure of genes × transcription factors. Rows correspond to regulated genes, columns correspond to transcription factors that exert regulatory effects, and the values ​​within the matrix represent the corresponding regulatory weights, indicating the strength of the regulatory effect of the corresponding transcription factor on the gene. Higher regulatory weights indicate more sufficient evidence of the transcription factor's regulation of the target gene and a stronger regulatory effect.

[0047] Understandably, with the development of spatial omics, the whole-pathway signal flow analysis method of this embodiment can be extended from suspended cell systems such as PBMCs to spatial transcriptome data. For example, a physical distance attenuation function can be introduced to replace the global distribution assumption of blood flow. In suspended cell systems (such as blood, lymph, etc.), signaling molecules are globally distributed with the flow of body fluids, and the contact probability between ligands and receptors mainly depends on the overall concentration of molecules in the system, with spatial proximity having a relatively weak effect on signal transduction. In this embodiment, the regulatory weight of ligand sensing adopts the global distribution assumption, that is, it is assumed that signaling molecules can act uniformly on cells throughout the system. However, in solid tissues (such as the brain, liver, tumor microenvironment, etc.), signaling molecules usually function in the local microenvironment in a paracrine manner, and their concentration gradient is significantly affected by the physical distance between cells. After being secreted by the sending cell, the ligand diffuses in the interstitial space, and its effective concentration reaching the receptor cell decreases with increasing distance. In this case, if the global distribution assumption is adopted, the probability of signal transduction between distant cells may be overestimated, leading to a large number of false positive pathways. Therefore, when the method of this embodiment is applied to spatial transcriptome data, the regulatory weights of ligand sensing can be adaptively modified. For example, the spatial proximity of spatial transcriptome data can be used as the regulatory weights of ligand sensing to construct prior regulatory network data.

[0048] Specifically, the two-dimensional or three-dimensional spatial coordinates of cells or spot points provided by spatial transcriptome data can be used to calculate the physical distance between any two spatial locations. Based on this physical distance, a distance decay function is introduced to replace the original global distribution assumption. This distance decay function includes, but is not limited to, Gaussian decay function and exponential decay function, etc., and its range is (0,1]. As a spatial weight coefficient, it can be used to recalibrate the regulatory weights in the prior regulatory network data. Based on the above spatial regulatory weight matrix, all subsequent steps in this embodiment, including the calculation of the transcription factor activity matrix, the screening of the biological pathway to be tested, the calculation of the probability intensity of the signal flow from the receptor to the transcription factor, the calculation of the signal regulation intensity from the transcription factor to the target gene, and the significance verification, can follow the original process, thereby achieving a seamless migration of the whole pathway signal flow analysis method from the suspension cell system to spatial transcriptome data. It can restore the concentration gradient distribution of signal molecules in the local microenvironment in solid tissues, solve the problem of quantitative inference of the signal intensity of the local microenvironment in solid tissues, and realize the close correlation between spatial location information and intracellular regulatory logic.

[0049] In one embodiment, after obtaining the single-cell transcriptome expression matrix and prior regulatory network data, to distinguish between them, genes in the single-cell transcriptome expression matrix are called expressing genes, and genes in the prior regulatory network data are called regulatory genes. This allows for the acquisition of the expression levels of expressing genes in different cells and the regulatory weights of different transcription factors corresponding to regulatory genes. (Refer to...) Figure 2 , Figure 2 This is a flowchart provided in this application embodiment for determining the transcription factor activity matrix based on the single-cell transcriptome expression matrix and prior regulatory network data, specifically including the following steps: Step 210: Obtain the intersection genes between the single-cell transcriptome expression matrix and the prior regulatory network data, construct the intersection expression matrix based on the expression levels of the intersection genes, and construct the intersection weight matrix based on the regulatory weights of the intersection genes.

[0050] In one embodiment, the intersection of all expressed genes in the single-cell transcriptome expression matrix and all regulatory genes in the regulatory network is obtained to obtain the intersection genes. These intersection genes represent genes that appear simultaneously in both the transcriptome expression data and the prior regulatory network, allowing for the acquisition of expression levels in different cells and regulatory weights on different transcription factors. Non-intersection genes are directly filtered out to prevent invalid data from participating in subsequent calculations.

[0051] Next, based on the list of intersection genes, the numerical values ​​of all expression levels for corresponding gene rows are extracted from the single-cell transcriptome expression matrix, preserving the original cell column structure to obtain the intersection expression matrix. This matrix still has a gene × cell structure and is a subset of the single-cell transcriptome expression matrix. Simultaneously, the numerical values ​​of all regulatory weights for corresponding gene rows are extracted from the prior regulatory network data, preserving the original transcription factor column structure to obtain the intersection regulatory network matrix. This matrix also has a gene × transcription factor structure and is a subset of the prior regulatory network data.

[0052] Therefore, there is a one-to-one correspondence between the row genes of the intersection expression matrix and the intersection weight matrix. Furthermore, to adapt to large-scale data computation in single cells, the intersection weight matrix can be converted to a sparse matrix format for storage, retaining only the non-zero regulatory weight values, thereby significantly reducing memory usage and improving subsequent computational efficiency.

[0053] Step 220: Calculate the transcription factor activity matrix based on the transpose of the intersection expression matrix and the intersection weight matrix.

[0054] In one embodiment, based on the biological hypothesis, the activity of a transcription factor can be reflected by the weighted sum of the expression levels of its downstream target genes; the higher the expression level of the target genes regulated by the transcription factor, the stronger the functional activity of the corresponding transcription factor. Therefore, in this embodiment, the activity value of each transcription factor in different cells is obtained based on the transpose of the intersection expression matrix and the intersection regulatory network matrix, and the activity values ​​constitute the transcription factor activity matrix.

[0055] The specific calculation formula is as follows:

[0056] in, The intersection expression matrix, The intersection weight matrix, This is the transpose of the intersection weight matrix. After transposing, the dimensions of the transpose matrix are converted to transcription factor × gene. This represents the transcription factor activity matrix, with dimensions of transcription factor × cell. Each element in the transcription factor activity matrix corresponds to the activity value of a transcription factor in a single cell. The transcription factor activity matrix calculation program employs sparse weighted matrix storage and matrix multiplication operations, which reduces memory usage and improves computational efficiency. It also supports multi-process parallel computing and task scripting, enabling efficient processing of large-scale single-cell transcriptome data (tens of thousands of cells) to meet the needs of cell communication analysis for large cohort samples.

[0057] Furthermore, as calculated above, higher regulatory weights and higher target gene expression levels correspond to higher transcription factor activity values. Therefore, the transcription factor activity matrix obtained in this embodiment can quantitatively reflect the functional activity of each transcription factor at single-cell resolution. This embodiment can improve the processing efficiency of large-scale single-cell data while ensuring the rationality of computational biology.

[0058] Step 120: Based on the prior pathway database, identify at least one experimental biological pathway from the single-cell transcriptome expression matrix.

[0059] In one embodiment, the biological pathway to be analyzed is the target of subsequent analysis, which includes ligands, receptors, transcription factors, and target genes. (See reference...) Figure 3 , Figure 3 This is a flowchart provided in this application embodiment for determining at least one experimental biological pathway from a single-cell transcriptome expression matrix based on a priori pathway database, specifically including the following steps: Step 310: Based on the prior pathway database, obtain at least one initial biological pathway.

[0060] In one embodiment, both the Minipath public database and the NicheNet dataset are used as prior pathway databases. Biologically validated signaling pathway node information is extracted from the Minipath public database. Simultaneously, the ligand-receptor interaction network and transcription factor-target gene regulatory network contained in the NicheNet dataset are combined to complete pathway splicing through multi-level upstream and downstream relationship mapping, resulting in multiple initial biological pathways. Each initial biological pathway includes at least a ligand, a receptor, a transcription factor, and a target gene. In the initial biological pathway, the ligand-receptor-transcription factor-target gene pathway follows the order of signal transduction: upstream is the ligand acting as an extracellular signal carrier, followed by the receptor located on the cell membrane that recognizes the ligand, then the transcription factor that receives the intracellular signal, and downstream is the target gene that performs the functional output, together forming a complete signal transduction pathway, represented as: ligand-receptor-transcription factor-target gene pathway.

[0061] Understandably, to improve information processing efficiency and narrow the scope of data processing, subsequent analysis can be conducted by identifying user-specified initial biological pathways from the initial biological pathways, based on the specific research objectives of the scenario. For example, a specific ligand can be designated as the starting point based on the research objective, and all initial biological pathways related to that ligand can be screened out, thereby narrowing the scope of analysis and focusing primarily on the downstream transduction pathways of the target signaling molecule.

[0062] Step 320: Obtain the expression levels of ligands, receptors, transcription factors, and target genes in different cells from the single-cell transcriptome expression matrix.

[0063] In one embodiment, the ligands, receptors, transcription factors and target genes contained in all initial biological pathways are matched one by one with the corresponding gene entries in the single-cell transcriptome expression matrix of the single-cell transcriptome, and the initial expression level of each gene in all cells is extracted.

[0064] Step 330: Determine whether the proportion of cells with non-zero expression levels in ligands, receptors, transcription factors, and target genes is greater than or equal to a preset threshold.

[0065] In one embodiment, for each ligand, receptor, transcription factor, and target gene, the total number of cells with initial non-zero expression levels is counted, and the proportion of non-zero expression cells to the total number of analyzed cells is calculated to reflect the general expression level of the corresponding gene in the sample. Then, this cell proportion is compared with a preset proportion, and genes with an initial expression proportion not lower than the preset proportion are selected as target genes. Ligands, receptors, transcription factors, or target genes with a non-zero expression proportion lower than the preset proportion are judged to lack a basis for general expression and are filtered out. This excludes genes that are expressed only in a very few cells and do not have population biological significance, reducing technical interference from random missed detections in single-cell sequencing and avoiding the introduction of low-expression genes into invalid pathway calculations. For example, the preset proportion can be set to 20%, determining that at least 20% of cells have non-zero expression levels of ligands, receptors, transcription factors, and target genes.

[0066] Step 340: If so, the initial biological pathway is taken as the biological pathway to be tested.

[0067] In one embodiment, an initial biological pathway whose ligands, receptors, transcription factors, and target genes all meet preset threshold conditions is selected as the biological pathway to be tested. That is, only when all ligands, receptors, transcription factors, and target genes contained in a pathway are selected as target genes through a non-zero expression percentage screening, will this initial biological pathway be retained as a biological pathway to be tested. Therefore, if the ligand, receptor, transcription factor, or target gene corresponding to any node in a pathway does not reach the preset expression percentage requirement, the entire pathway will be filtered out. All the biological pathways to be tested ultimately possess a complete signal transduction structure, and the genes at each node have reached the minimum expression percentage requirement in the current sample, ensuring that each retained biological pathway to be tested has the biological expression basis for achieving complete signal transduction.

[0068] In the above embodiments, a uniform expression percentage is used as a threshold to filter all genes in the pathway. This ensures that the retained biological pathways to be tested have genuine expression support in the current sample, while also reducing the data scale for subsequent calculations and improving overall computational efficiency. Furthermore, this screening method does not rely on the spatial location information of cells, making it suitable for suspended cell samples in fluid environments such as blood. It effectively filters out spurious pathways included in the database but not actually expressed in the sample, reducing the probability of false positives in subsequent signal flow inference and ensuring the biological rationality and reliability of cell communication analysis results.

[0069] Step 130: Traverse the biological pathways to be tested, and determine the pathway activation parameters of the biological pathways to be tested in single cells based on the single-cell transcriptome expression matrix, prior regulatory network data, and transcription factor activity matrix.

[0070] In one embodiment, it is necessary to traverse the biological pathways under investigation and analyze each pathway. Pathway activation parameters include the probability intensity of signal flow from receptor to transcription factor and the intensity of signal regulation from transcription factor to target gene. The analysis process of any one of the biological pathways under investigation is described below as an example. (Refer to...) Figure 4 , Figure 4 This is a flowchart provided in this application embodiment for determining the pathway activation parameters of a biological pathway in a single cell based on a single-cell transcriptome expression matrix, prior regulatory network data, and transcription factor activity matrix. The flowchart specifically includes the following steps: Step 410: Based on the single-cell transcriptome expression matrix, determine the mean expression levels of receptors, transcription factors, and target genes in each cell, and obtain the receptor half-saturation coefficient, transcription factor half-saturation coefficient, and target gene half-saturation coefficient.

[0071] In one embodiment, the mean expression level here can be determined based on the single-cell transcriptome expression matrix or the intersection expression matrix obtained in step S220. For the receptor, transcription factor, and target gene in the current experimental biological pathway, the corresponding half-saturation coefficients need to be calculated, namely the receptor half-saturation coefficient, the transcription factor half-saturation coefficient, and the target gene half-saturation coefficient. The half-saturation coefficient refers to the level of biochemical activation of a gene when its expression level or activity equals the half-saturation coefficient, indicating that the gene participates in exactly 50% of the activation.

[0072] The formula for calculating the half-saturation coefficient of genes in this embodiment is:

[0073] in, This indicates the non-zero expression level of the corresponding gene. This represents the total number of non-zero expression levels. In this embodiment, the half-saturation coefficient K is adaptively determined using the gene's own non-zero mean, without relying on a universal fixed value, making the half-saturation coefficient more closely match the tissue environment and physiological state of the research scenario.

[0074] The receptor half-saturation coefficient can be obtained according to the above calculation process. transcription factor half-saturation coefficient and target gene half-saturation coefficient .

[0075] Step 420: For receptors and transcription factors, determine the probability intensity of signal flow based on the single-cell transcriptome expression matrix, prior regulatory network data, receptor half-saturation coefficient, transcription factor half-saturation coefficient, and transcription factor activity matrix.

[0076] In one embodiment, the activation process in biochemistry is simulated using a first-order kinetic equation to calculate the probability intensity of the signal flow. It is understood that the kinetic equation can be the Hill equation, the Michaelis-Menten equation, logistic regression, or an exponential decay function, etc. The Hill equation is used as an example below to describe the process of simulating the saturation effect.

[0077] In one embodiment, reference is made to Figure 5 , Figure 5 This is a flowchart provided in this application embodiment for determining the probability intensity of signal flow based on a single-cell transcriptome expression matrix, prior regulatory network data, receptor half-saturation coefficient, transcription factor half-saturation coefficient, and transcription factor activity matrix. The flowchart specifically includes the following steps: Step 510: Obtain the receptor expression level corresponding to the receptor in the single-cell transcriptome expression matrix, and the transcription factor activity value corresponding to the transcription factor in the transcription factor activity matrix.

[0078] In one embodiment, for each receptor and transcription factor pairing in the biological pathway under investigation, corresponding quantitative data are extracted. The expression level of the receptor-corresponding gene is extracted from the single-cell transcriptome expression matrix as the receptor expression level. The activity values ​​of corresponding transcription factors are extracted from the transcription factor activity matrix as transcription factor activity values. .

[0079] Step 520: Calculate the proportion of receptor expression level in the sum of receptor expression level and receptor half-saturation coefficient to obtain the receptor activation proportion; calculate the proportion of transcription factor activity value in the sum of transcription factor activity value and transcription factor half-saturation coefficient to obtain the transcription factor activation proportion.

[0080] In one embodiment, receptor expression levels are calculated during the first-order Hill equation calculation. and receptor half-saturation coefficient sum Then calculate the receptor expression level. exist The percentage of receptor activation is obtained by taking the proportion of the receptor activation percentage, which is expressed as: And the calculation of transcription factor activity values. and transcription factor half-saturation coefficient sum Then calculate the transcription factor activity value. exist The proportion of transcription factor activation is obtained by taking the percentage of the total number of transcription factors, which is expressed as: .

[0081] In this study, higher receptor expression levels indicate a higher degree of receptor activation. When receptor expression equals the receptor half-saturation coefficient, the receptor activation percentage is 0.5%, corresponding to a 50% activation probability. When receptor expression is significantly higher than the receptor half-saturation coefficient, the receptor activation percentage approaches one, reaching a state of complete saturation. Conversely, when receptor expression is significantly lower than the receptor half-saturation coefficient, the receptor activation percentage approaches zero, indicating almost no activation response. It can be understood that the transcription factor activation percentage, similar to the receptor activation percentage, is used to quantitatively reflect the degree of functional activation of transcription factors.

[0082] Step 530: Obtain the physical activation probability of receptor-activated transcription factors based on the receptor activation ratio and the transcription factor activation ratio.

[0083] In one embodiment, the physical activation probability is obtained by multiplying the receptor activation percentage by the transcription factor activation percentage and taking the square root, which is the geometric mean of the two, and is expressed as:

[0084] This embodiment uses geometric mean for state coupling, which reflects the cascade and logic of signal transduction. It indicates that when this pathway is activated, the receptor expresses a signal, and the transcription factor is also active. In other words, effective signal transduction requires both receptor activation and transcription factor activation to be met simultaneously. Extremely low activation levels at either end will significantly reduce the overall physical activation probability. Compared to a direct arithmetic mean, it more closely reflects the actual cascade and transduction patterns of signals.

[0085] Step 540: Obtain the signal flow probability intensity based on the physical activation probability and the regulatory weights of transcription factors in the prior regulatory network data.

[0086] In one embodiment, the weight is adjusted By multiplying the physical activation probability by the prior regulatory network data derived from prior data by the prior regulatory weights of the corresponding pathways, the probability intensity of the receptor-to-transcription factor signal flow can be obtained. , is represented as:

[0087] Among them, signal flow probability intensity This indicates the probability that a receptor signal will successfully transduce to a transcription factor in a single cell, determined by both the receptor expression level and the transcription factor's activity level. This embodiment introduces prior regulatory weights to constrain data-driven physical activation probabilities using existing biological knowledge, reducing false positives caused by random expression fluctuations and improving the biological validity of signal inference results.

[0088] Step 430: For transcription factors and target genes, determine the intensity of signal regulation based on single-cell transcriptome expression matrix, prior regulatory network data, transcription factor half-saturation coefficient, target gene half-saturation coefficient, and transcription factor activity matrix.

[0089] In one embodiment, since transcription factors exhibit a synergistic effect in regulating target genes during normal biological processes—meaning a transcription factor needs to form a dimer or multiple molecules to bind together in the promoter region to effectively drive gene expression—this embodiment, based on this assumption, introduces a second-order kinetic equation to simulate the regulation of transcription factors and target genes. For example, the second-order kinetic equation could be a second-order Hill equation. (Refer to...) Figure 6 , Figure 6 This is a flowchart provided in this application embodiment for determining the intensity of signal regulation based on a single-cell transcriptome expression matrix, prior regulatory network data, transcription factor half-saturation coefficients, target gene half-saturation coefficients, and transcription factor activity matrices. The flowchart specifically includes the following steps: Step 610: Obtain the expression levels of target genes corresponding to target genes in the single-cell transcriptome expression matrix and the activity values ​​of transcription factors corresponding to transcription factors in the transcription factor activity matrix.

[0090] Step 620: Obtain the second-order activity value based on the transcription factor activity value and the second-order half-saturation coefficient based on the transcription factor half-saturation coefficient.

[0091] In one embodiment, target gene expression level It is the expression level of the target gene and its second-order activity value extracted from the single-cell transcriptome expression matrix. transcription factor activity value The square of the second-order half-saturation coefficient Transcription factor half-saturation coefficient The square of . A second-order kinetic model is used to form a hysteresis phase in the low signal region to filter background noise in single-cell sequencing, and to exhibit an on / off response in the threshold region to capture the core signal.

[0092] Step 630: Calculate the proportion of target gene expression level in the sum of target gene expression level and target gene half-saturation coefficient to obtain the target gene intensity proportion; and calculate the proportion of second-order activity value in the sum of second-order activity value and second-order half-saturation coefficient to obtain the transcription factor intensity proportion.

[0093] In one embodiment, the target gene intensity ratio is calculated using the first-order Hill equation to reflect the response saturation level corresponding to the target gene expression level, and the target gene expression level is calculated. and target gene half-saturation coefficient sum Then calculate the expression level of the target gene. exist The proportion of target gene activation is obtained by taking the percentage of the target gene activation, and is expressed as: The intensity ratio of transcription factors was calculated using a second-order Hill equation to effectively suppress false-positive activations caused by background noise in single-cell sequencing, improve the algorithm's robustness to sequencing noise, and calculate the second-order activity value. and second-order half-saturation coefficient sum Then calculate the second-order activity value. exist The proportion of transcription factor intensity is obtained by taking the percentage of the total number of transcription factors, and is expressed as: .

[0094] Step 640: Obtain the physical response intensity of transcription factor activation of target gene based on the proportion of target gene intensity and the proportion of transcription factor intensity.

[0095] In one embodiment, the physical response intensity is obtained by multiplying the proportion of target gene intensity and the proportion of transcription factor intensity and then taking the square root, which is the geometric mean of the two, and is expressed as:

[0096] Step 650: Obtain the signal regulation intensity based on the physical response intensity and the regulatory weights of transcription factors in the prior regulatory network data.

[0097] In one embodiment, the signal regulation intensity from the transcription factor to the target gene can be obtained by multiplying the physical response intensity by the regulatory weight of the transcription factor in the corresponding pathway. , is represented as:

[0098] Among them, signal modulation intensity This indicates the strength of regulation in which a transcription factor successfully drives the expression of a target gene in a single cell, determined by both the activity level of the transcription factor and the expression level of the target gene in that cell.

[0099] In one embodiment, reference is made to Figure 7 , Figure 7 This is a schematic diagram comparing the nonlinear response characteristics of the first-order Hill function and the second-order Hill function provided in the embodiments of this application. The horizontal axis represents the transcription factor activity value, and the vertical axis represents the activation probability output value, representing the degree of activation of the target gene by the transcription factor, with values ​​ranging from 0 to 1. The diagram contains two response curves. The dashed line corresponds to the activation process model of the first-order Hill equation (shown in the diagram as n=1), exhibiting a hyperbolic shape and high noise sensitivity (premature response). The solid line corresponds to the regulatory simulation model of the second-order Hill equation (shown in the diagram as n=2), exhibiting an S-shaped shape and noise filtering (hysteresis phase) in the low-activity region. The intersection of the horizontal and vertical dashed lines corresponds to the position of the transcription factor half-saturation coefficient. At this point, the transcription factor activity value equals the transcription factor half-saturation coefficient, and the output probability of both curves is 0.5, meaning the activation probability of the target gene is exactly fifty percent. Among them, the second-order model has a significant hysteresis phase in the low-activity region, which can be used to effectively suppress false positive activation caused by single-cell sequencing background noise and improve the robustness of the algorithm to sequencing noise; when the activity level is close to the half-saturation threshold, the response level rises rapidly, exhibiting switch-like activation characteristics, which can truly restore the transcription switch effect in biology and achieve high-sensitivity capture of core signal flow.

[0100] This application's embodiments achieve layer-by-layer verification of signal transduction by segmentally quantifying the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene. This restores the entire process of signal recognition from extracellular to intracellular transcriptional regulation, avoiding the analytical defects of inference based solely on receptor-ligand co-expression or direct ligand-target gene inference. It effectively eliminates false-positive communication events that only involve receptor-ligand interaction without triggering downstream transcriptional regulation, thereby improving the biological rationality and accuracy of communication inference.

[0101] Step 140: Based on the signal modulation intensity and signal flow probability intensity, perform significance verification on the biological pathway to be tested in order to determine the target signal flow pathway.

[0102] In one embodiment, to distinguish between real biological regulatory pathways and false positive signals generated by random fluctuations in data, it is necessary to construct a null hypothesis distribution of signal intensity through permutation tests, and then use multiple hypothesis testing to correct and screen statistically reliable biological pathways to be tested. Based on the verification results, the target signal pathway is determined to ensure the reliability of the results. (Refer to...) Figure 8 , Figure 8 This is a flowchart provided in this application embodiment for significant verification of the biological pathway under test based on signal modulation intensity and signal flow probability intensity, specifically including the following steps: Step 810: Obtain the pathway activation parameters of a single cell based on the product of the signal modulation intensity and the signal flow probability intensity.

[0103] In one embodiment, for each biological pathway to be tested, the probability intensity of the signal flow from the receptor to the transcription factor corresponding to the same pathway is first multiplied by the signal regulation intensity from the transcription factor to the target gene at single-cell resolution. The product of the two constitutes the pathway activation parameter of the entire pathway in that cell.

[0104] Step 820: Calculate the average value of the pathway activation parameters of the biological pathway under test in all cells within a specific cell population to obtain the true observed signal intensity of the biological pathway under test.

[0105] In one embodiment, the arithmetic mean of the pathway activation parameters of all cells within a specific cell population corresponding to the same cell type or cell state is taken to obtain the actual observed signal intensity of the biological pathway under test in the corresponding cell population. , is represented as:

[0106] Where N represents the total number of cells in the population. This represents the signal flow probability intensity of the i-th cell. This represents the signal regulation intensity of the i-th cell. The actual observed signal intensity represents the actual signal transduction level of the biological pathway under investigation in the corresponding cell population.

[0107] Step 830: Randomly shuffle the pairings of receptors and transcription factors or transcription factors and target genes in the biological pathways to be tested to obtain the corresponding randomized biological pathways.

[0108] In one embodiment, the purpose of the permutation test is to maintain the expression level and activity value distribution of each gene unchanged, while disrupting only the true regulatory pairing relationships between genes, thereby constructing a null hypothesis under a random background and calculating the signal intensity under the random background. Therefore, in this embodiment, for each biological pathway to be tested, regarding the receptor-transcription factor pairing, the expression levels of all receptors and the activity values ​​of all transcription factors are kept constant, and the pairing relationships between receptors and transcription factors are randomly shuffled to generate randomly paired receptor-transcription factor pairs. Similarly, regarding the transcription factor-target gene pairing, the activity values ​​of all transcription factors and the expression levels of all target genes are kept constant, and the pairing relationships between transcription factors and target genes are randomly shuffled to generate randomly paired transcription factor-target gene pairs. Each time a randomly paired receptor-transcription factor or transcription factor-target gene pairing is generated, replacing the original pairing relationship, a randomized biological pathway can be obtained.

[0109] Understandably, by repeating the above random shuffling operation multiple times according to a pre-set total number of permutations, for example, 1000 times, a sufficient set of random biological pathways can be obtained to construct the null distribution of signal intensity for the corresponding biological pathway to be tested. Alternatively, permutations can be performed only on the receptor-transcription factor or transcription factor-target gene pairings to verify the significance of these pairings.

[0110] Step 840: Calculate the corresponding virtual signal intensity for each random biological pathway, calculate the significance P-value of the biological pathway under test based on the actual observed signal intensity and the virtual signal intensity, and obtain the significance verification result based on the significance P-value.

[0111] In one embodiment, for each random biological pathway, a virtual signal intensity is calculated using the same computational process as the actual observed signal intensity of the real pathway. The virtual signal intensities of all random pathways collectively constitute a zero-distribution of signal intensity under a random background, representing the signal intensity level achievable solely by random fluctuations in the data. Next, the significance p-value of the biological pathway under test is calculated based on the zero-distribution, expressed as:

[0112] Where P represents the significance P-value, and n represents the total number of permutations. Indicates the actual observed signal strength. This represents the virtual signal intensity corresponding to the random biological pathway obtained by the i-th random permutation. I is an indicator function; it takes the value 1 when the condition in parentheses is true (i.e., the virtual signal intensity is greater than or equal to the actual observed signal intensity), and 0 when the condition is false. Increasing both the numerator and denominator by 1 performs pseudo-counting to ensure the significance P-value is never zero, avoiding absolute significance conclusions due to the limited number of permutations. Here, the smaller the significance P-value, the lower the probability that the signal of the tested biological pathway is generated by random fluctuations, and the higher the biological authenticity of the pathway.

[0113] In one embodiment, after calculating the significance P-values ​​for all biological pathways to be tested, the Benjamini-Hochberg method can be used to perform multiple hypothesis testing corrections on all pathways to prevent false positives, effectively eliminate random noise interference caused by insufficient single-cell sequencing depth and heterogeneity, and ensure the statistical significance of complex signal flow path inference. The specific correction formula is expressed as follows:

[0114] in, Let be the corrected significance P-value for the i-th biological pathway to be tested, m be the total number of biological pathways to be tested, and i be the rank of the significance P-values ​​in ascending order. The significance P-value is the value for the i-th position after sorting. This represents the p-value of the corrected significance of the pathway following the highest rank. It can be seen that the correction process is calculated sequentially forward from the highest rank, and the minimum value is used to ensure the monotonicity of the correction result.

[0115] Next, using the corrected significance P-value as the screening criterion, biological pathways to be tested with a corrected significance P-value lower than a preset threshold are selected as statistically significant target signaling pathways. The preset threshold can be determined according to the actual scenario.

[0116] Therefore, this embodiment constructs a random background distribution through permutation testing, eliminating the need for pre-setting the data distribution pattern and adapting to the characteristics of strong heterogeneity and complex distribution in single-cell data. After validation, false positive pathways caused by insufficient single-cell sequencing depth, cellular heterogeneity, and random expression fluctuations can be effectively eliminated, ensuring that the final output signal pathway has statistical significance, thus improving the accuracy and reproducibility of communication analysis results in complex biological samples.

[0117] In one embodiment, once the target signaling pathway is identified, an upset diagram can be used to view the differential pathway performance based on the pathway characteristics of different cell types, and to view information such as key differentially expressed transcription factors. The differential pathway characteristics are visualized in the form of a chord diagram. This embodiment does not specifically limit the visualization format of the target signaling pathway.

[0118] In one embodiment, under different fluid mixing intensities, the full-pathway signal flow analysis method of this embodiment is compared with pathway signal analysis methods in related technologies, such as cell communication inference methods that attenuate and weight pathway signals based on spatial physical distance, to verify the technical advantages of this embodiment in resolving target biological pathways in suspended cell samples such as blood. This embodiment designs and implements rigorous baseline pressure tests to resist fluid disturbances. In the experimental design, the core verification objective is to evaluate the accuracy of intracellular full-pathway signal flow recovery when receiving spatial microenvironment signal inputs of different qualities.

[0119] Reference Figure 9 , Figure 9This is a comparative diagram of verification results provided in this application embodiment. The horizontal axis represents the spatial perturbation ratio, i.e., the fluid mixing intensity. A higher value indicates a higher degree of disorder in the spatial location of cells. 0 represents solid tissue data such as lungs and heart, and 100 corresponds to a highly fluid state. The gray area, ranging from 80% to 100%, corresponds to a highly fluid state, such as the physiological environment of blood or suspended cell samples. The vertical axis represents the prediction accuracy, using the AUC-ROC value as a quantitative indicator. A higher value indicates better accuracy in the signal flow path analysis results. The dotted lines in the diagram represent the analysis results of analytical models in related technologies, the solid lines with squares represent the analysis results of this embodiment, and the horizontal dotted lines correspond to the baseline performance level of random guessing, with an AUC value of 0.5.

[0120] The simulation experiments employed a closed-loop control design to ensure the robust persuasiveness of the comparative conclusions. The first experimental design element was a gold standard setting biased towards the control scheme. When generating cell pairs labeled as real communication events, a mandatory spatial proximity constraint was added, requiring the spatial distance between paired cells to be less than a set threshold. In a static tissue scenario with zero disturbance, the relevant technical solutions were under optimal prediction conditions. This setting was entirely beneficial to the relevant technical solutions, and the performance comparison results obtained under this premise had sufficient evidentiary power. The second element was a high-fidelity single-cell sequencing noise simulation. The expression levels of ligand and receptor genes were generated using discrete sampling from a Poisson distribution, simulating the low expression and missed gene detection phenomena caused by insufficient sequencing depth in real single-cell sequencing. This restored the inherent noise characteristics of single-cell transcriptome data, ensuring the experimental environment closely resembled the real sequencing scenario. The third element was a precise mapping of fluid physics processes. By progressively increasing the scrambling ratio of cell spatial coordinates to simulate the fluid microenvironment, while maintaining the gene expression profiles of all cells completely unchanged, the degree of disorder in cell spatial positions was gradually increased. Physically, this corresponded to the physiological process of continuous dynamic changes in the spatial topology of cells in fluid microenvironments such as blood and lymph due to fluid scouring.

[0121] according to Figure 9Experimental results show that when the fluid mixing intensity is 0, corresponding to a static solid tissue scenario, the AUC-ROC value of the relevant technical solutions is close to 1, and the prediction accuracy is at its highest level, which is consistent with the optimal applicable scenario of its spatial proximity assumption. The AUC-ROC value of this embodiment is stable at around 0.78, maintaining stable and reliable prediction performance. As the fluid mixing intensity gradually increases, the prediction accuracy of the relevant technical solutions shows a continuous downward trend. The higher the degree of spatial disturbance, the greater the inference deviation based on the spatial distance assumption. When the fluid mixing intensity reaches 80% or more, entering a highly fluid state, the AUC-ROC value of the relevant technical solutions drops below 0.76 and continues to decline. Its spatial distance assumption fails in the suspended cell environment, and the prediction accuracy decreases significantly. The prediction accuracy of this embodiment remains stable throughout, unaffected by the degree of spatial disturbance, and its prediction accuracy in the highly fluid state is significantly higher than that of the relevant technical solutions.

[0122] Therefore, it can be seen that related technologies, such as cell communication models based on spatial distance, are only suitable for static solid tissue samples, and their performance is significantly reduced in highly fluid suspended cell samples such as blood. The whole-pathway signal flow analysis method in this embodiment does not rely on cell spatial location information. It completes the signal flow inference of target biological pathways based on the single-cell transcriptome expression matrix and prior regulatory network data. It can maintain stable prediction performance in scenarios with highly dynamic cell spatial locations. It has higher inference accuracy and scenario adaptability in suspended single-cell data such as blood. It is suitable for application scenarios such as disease target screening of large clinical cohort blood samples and has good clinical translational application value.

[0123] In one embodiment, a specific application of the present application in blood PBMC single-cell transcriptome data is described. Using BGIDNBelabC4 samples from the centenarian cohort and CIMA cohort in a public database, the signaling pathways received by peripheral blood mononuclear cells are analyzed at the cell type level. Ultimately, it is verified that the upregulated LTA protein in the blood of centenarians can regulate the upregulation of target genes such as ISG20 and LYAR by immune cells.

[0124] The first step is data loading and preprocessing. Single-cell transcriptome data files are imported, reading single-cell data objects in formats compatible with Seurat V5 and earlier versions. The single-cell transcriptome expression matrix is ​​extracted, and specific cell types are selected as target cell populations for signal reception for subsequent pathway analysis.

[0125] Alternatively, epigenetic regulation reference data can be imported to read the correspondence file between transcription factors, regulatory regions, and target genes. Epigenetic evidence can be introduced to strongly constrain the regulatory relationship between transcription factors and target genes. This step is suitable for analysis scenarios with paired ATAC sequencing data and can further reduce false positives in regulatory relationships.

[0126] Then, the prior NicheNet prior regulatory network data is loaded, including multiple prior pathway databases containing weighted regulatory networks, ligand-transcription factor association matrices, human ligand-receptor interaction networks, and ligand-target gene association matrices, as prior references for constructing target biological pathways. This dataset can be replaced with the corresponding reference dataset according to the research species and research direction.

[0127] Next, the target analysis ligand is defined, LTA is selected as the extracellular signaling protein of interest, the starting point for subsequent pathway analysis is limited, and the analysis scope is narrowed to focus on the downstream transduction pathway of the target signal.

[0128] Signal analysis is performed for each biological pathway to be tested. First, a transcription factor activity matrix is ​​calculated. Based on the regulatory weights in the prior regulatory network data, the expression levels of corresponding target genes in each cell are weighted and summed, converting the genome-wide expression data into transcription factor activity values, and outputting the transcription factor activity matrix. For example, in this process, the input data consisted of 33,559 genes and 33,202 cells. Gene matching yielded 12,613 intersection genes. After constructing a sparse weight matrix, activity calculations were performed using matrix multiplication, ultimately obtaining a single-cell resolution transcription factor activity matrix.

[0129] Then, the target biological pathway backbone is extracted. Based on prior regulatory network data and the defined target ligand, pathway splicing is completed through multi-level regulatory relationship mapping, resulting in multiple initial biological pathways. Expression filtering is then performed on all initial genes in the initial biological pathways, retaining only pathways where all genes meet the expression threshold, resulting in a set of experimental biological pathways containing ligands, receptors, transcription factors, and target genes, while retaining the regulatory weights corresponding to each segment. Example code is as follows: sankey_skeleton<-Sigflow_pp_sankey_skeleton(seurat_obj=target_obj,ligand_list=ligands_oi,lr_network=lr_network,weighted_networks=weighted_networks,regulon_df=regulon_df,min_pct=0.2,min_weight_gr=0.1) The sankey_skeleton function, when applied, ultimately splices together 3696 complete biological pathways to be tested, including 1 ligand, 3 receptors, 26 transcription factors, and 1543 target genes.

[0130] Next, the half-saturation coefficients are calculated. The non-zero expression levels of receptors, target genes, and transcription factors are extracted from the single-cell transcriptome expression matrix, and the arithmetic mean of the corresponding non-zero quantitative values ​​is calculated. This yields the half-saturation coefficients for receptors, transcription factors, and target genes. Example code is as follows: k_vals<-Sigflow_tl_calc_k(sankey_df=sankey_df,seurat_obj=target_obj,tf_activity_matrix=tf_act) The `k_vals` function extracts the non-zero expression levels of receptors, target genes, and transcription factors from the single-cell transcriptome expression matrix, calculates the arithmetic mean of the corresponding non-zero quantitative values, and obtains the receptor half-saturation coefficient, transcription factor half-saturation coefficient, and target gene half-saturation coefficient. After execution, a total of 3 receptor half-saturation coefficients, 26 transcription factor half-saturation coefficients, and 1543 target gene half-saturation coefficients were calculated.

[0131] Next, the probability intensity of receptor-to-transcription factor signal flow is calculated. For each receptor-transcription factor pair, the probability intensity of receptor-to-transcription factor signal flow is obtained by combining the corresponding receptor expression level, receptor half-saturation coefficient, transcription factor activity value, and transcription factor half-saturation coefficient. The calculation is performed grouped by cell type, and the probability intensity of signal flow at the cell type level is output. Example code is as follows: rtf_prob_df<-Sigflow_tl_rtf_prob(sankey_df=sankey_df,seurat_obj=target_obj,tf_activity_matrix=tf_act,k_values=k_vals,group_by="celltype") The `Sigflow_tl_rtf_prob` function is responsible for calculating the probability intensity of signal flow from the receptor to the transcription factor. It uses the Hill equation to simulate the response logic of biological signals and uses the geometric mean to represent biological signal transduction. During execution, all 3696 original pathways passed validity checks. After deduplication, 44 receptor-transcription factor pairs were obtained. Signal intensity calculations were performed by grouping by cell type, and the signal intensity results at the cell type level were output.

[0132] Next, the signal regulation strength from transcription factors to target genes is calculated. For each transcription factor-target gene pairing, the signal regulation strength is determined based on the single-cell transcriptome expression matrix, prior regulatory network data, transcription factor half-saturation coefficient, target gene half-saturation coefficient, and transcription factor activity matrix. The calculation is also performed grouped by cell type, outputting the regulation strength results at the cell type level. Example code is as follows: target_prob_df<-Sigflow_tl_target_prob(sankey_df=sankey_df,seurat_obj=target_obj,tf_activity_matrix=tf_act,k_values=k_vals,regulon_df=regulon_df,group_by="celltype") The Sigflow_tl_target_pro function is responsible for calculating the signal regulation intensity of the pathway from transcription factor to target gene. Its logic is similar to that of the Sigflow_tl_rtf_prob function, except that it uses a second-order Hill equation to simulate the cascade amplification effect in biology.

[0133] Finally, statistical significance verification and pathway screening were performed. Permutation tests were conducted on the receptor-to-transcription factor pathway, with a total of 1000 permutations. Randomized biological pathways were generated by randomly shuffling the receptor-to-transcription factor pairings. Virtual signal intensities were generated using computational logic identical to that of real pathways, and a null hypothesis distribution of the signal intensities was constructed. Based on this, the significance P-value for each receptor-to-transcription factor pathway was calculated, and the corrected significance P-value was obtained after multiple hypothesis testing.

[0134] For the permutation test of the transcription factor-to-target gene pathway, 1000 random permutation operations were performed. By randomly shuffling the pairing relationship between transcription factors and target genes to construct a random background, the corrected significance P-value of each transcription factor-target gene pathway was calculated.

[0135] Next, we perform saliency filtering, as shown in the example code below: final_results<-Sigflow_st_filter_pathways(rtf_res=rtf_res_final,target_res=target_res_final,p_adj_cutoff=0.05) The `Sigflow_st_filter_pathways` function is a filtering module that filters out pathways with significance levels higher than expected. Ultimately, 2506 statistically significant target signal flow paths are selected as the final signal flow path analysis results.

[0136] Finally, the signal flow path analysis results are visualized and analyzed. For example, an intersection analysis Upset diagram can be used for visualization, as shown in the following example code: upset_plot<-Sigflow_pl_pathway_upset(df=final_results,group_col="group") The Sigflow_pl_pathway_upset function takes as input the final signal flow path analysis result after saliency filtering. The function automatically counts all intersection combinations of paths between different groups and generates a visualization graph.

[0137] In one embodiment, reference is made to Figure 10 , Figure 10 This is a visualization of the results. The left side is a horizontal bar chart, with the vertical axis listing eight immune cell subsets: MK, DC, NK, monocytes, CD8T cells, B cells, plasma cells, and CD4T cells. The length of the bar corresponds to the number of target signaling pathways detected in each cell subset, visually reflecting the richness of pathways regulated by different immune cells. In the dot matrix chart, each column corresponds to a pathway intersection combination. Solid dots in a column indicate that the corresponding cell type is included in the intersection combination; no dots indicate that the corresponding cell type does not participate in the intersection. By arranging and combining the dots, the combination patterns of all pathway sets can be identified, including single-cell-specific pathways, pathways shared by two cell types, and pathways shared by multiple cell types. The vertical axis of the bar chart represents the gene intersection scale, i.e., the number of pathways under the corresponding intersection combination, with a specific value marked at the top of each bar. The higher the bar, the more pathways are in that intersection pattern. In this embodiment, the largest shared pathway contains 788 pathways, corresponding to a set of signaling pathways shared by multiple types of immune cells, reflecting the common regulatory characteristics of the immune regulatory network. The remaining bars, which decrease in size, correspond to shared pathways and cell-specific pathways of different sizes, reflecting the specificity of signal responses of different cell subpopulations.

[0138] In one embodiment, reference is made to Figure 11 , Figure 11 This is a schematic diagram of cell signal flow results. It can visually display all statistically significant target biological pathways starting with LTA as the initiating ligand, and generate a visual network by arranging all effective pathways hierarchically, presenting the hierarchical transduction network of signals from extracellular ligands through receptors and transcription factors to target genes.

[0139] The diagram illustrates a four-layered ring structure from the inside out. The central node is the ligand LTA, a key differentially expressed signaling protein identified through analysis. Upregulated in the blood of centenarians, it serves as the initiation point for extracellular signal transduction. The second layer contains receptor nodes, including TNFRSF1B, TNFRSF14, LTBR, ​​and TNFRSF1A. These are signal recognition nodes located on the cell membrane surface, specifically recognizing the ligand LTA and transducing the signal into the cell. The third layer contains transcription factor nodes, including multiple transcription factors such as STAT2, IRF3, FLI1, GATA1, MAFB, SPI1, ELF1, REL, and JUND. These are intermediate regulatory nodes for intracellular signal transduction, and their activity values ​​are derived from the transcription factor activity matrix. The outermost layer contains target gene nodes, including multiple downstream functional genes such as ISG20, LYAR, VAMP5, EMP3, and MIDN. These are functional output nodes for signal transduction, and their expression levels are co-regulated by upstream transcription factors.

[0140] The connections between different nodes correspond to effective transduction pathways that have been quantitatively calculated and statistically verified. The connections between ligands and receptors correspond to ligand-receptor interactions verified by prior networks; the connections between receptors and transcription factors correspond to the probability intensity of signal flow from receptors to transcription factors; and the connections between transcription factors and target genes correspond to the intensity of signal regulation from transcription factors to target genes. This visualization clearly reconstructs the entire regulatory chain of ligand LTA through multi-receptor and multi-transcription factor levels regulating downstream target genes, and intuitively verifies that the upregulated LTA protein in the blood of centenarians can drive the upregulation of immune-related target genes such as ISG20 and LYAR in immune cells through a complete intracellular signaling pathway.

[0141] The technical solution provided in this application constructs a complete signal transduction pathway including receptors, transcription factors, and target genes, restoring the biological process from extracellular recognition to intracellular transcriptional regulation, and reducing the false positive probability caused by relying solely on receptor-ligand co-expression or direct inference between ligand and target genes. It is adaptable to cell sample analysis in fluid environments such as blood, improving the completeness, accuracy, and scenario adaptability of cell communication analysis results.

[0142] This application also provides an omnichannel signal flow analysis device, which can implement the above-described omnichannel signal flow analysis method, as described above. Figure 12 The device includes: Activity matrix calculation module 1210: used to obtain the single-cell transcriptome expression matrix and prior regulatory network data of the sample, and to determine the transcription factor activity matrix based on the single-cell transcriptome expression matrix and prior regulatory network data.

[0143] Pathway initial filtering module 1220: Used to identify at least one experimental biological pathway from a single-cell transcriptome expression matrix based on a priori pathway database. The experimental biological pathway includes at least a ligand, receptor, transcription factor, and target gene.

[0144] Quantitative calculation module 1230: used to traverse the biological pathway to be tested, and determine the pathway activation parameters of the biological pathway to be tested in a single cell based on the single-cell transcriptome expression matrix, prior regulatory network data and transcription factor activity matrix; the pathway activation parameters include the signal flow probability intensity from the receptor to the transcription factor and the signal regulation intensity from the transcription factor to the target gene.

[0145] Verification output module 1240: used to perform significance verification on the biological pathway under test based on signal modulation intensity and signal flow probability intensity, so as to determine the target signal flow pathway.

[0146] The specific implementation of the omnichannel signal flow analysis device in this embodiment is basically the same as the specific implementation of the omnichannel signal flow analysis method described above, and will not be repeated here.

[0147] This application also provides an electronic device, including: At least one memory; At least one processor; At least one program; The program is stored in a memory, and the processor executes the at least one program to implement the full-path signal flow analysis method described above. The electronic device can be any smart terminal, including mobile phones, tablets, personal digital assistants (PDAs), and in-vehicle computers.

[0148] Please see Figure 13 , Figure 13 The hardware structure of an electronic device according to another embodiment is illustrated. The electronic device includes: The processor 1301 can be implemented using a general-purpose central processing unit (CPU), microprocessor, application specific integrated circuit (ASIC), or one or more integrated circuits, and is used to execute relevant programs to implement the technical solutions provided in the embodiments of this application. The memory 1302 can be implemented as a read-only memory (ROM), static storage device, dynamic storage device, or random access memory (RAM). The memory 1302 can store the operating system and other application programs. When the technical solutions provided in the embodiments of this specification are implemented through software or firmware, the relevant program code is stored in the memory 1302 and is called and executed by the processor 1301 using the full-path signal flow analysis method of the embodiments of this application. The input / output interface 1303 is used to implement information input and output; The communication interface 1304 is used to enable communication and interaction between this device and other devices. Communication can be achieved through wired means (such as USB, Ethernet cable, etc.) or wireless means (such as mobile network, WIFI, Bluetooth, etc.). Bus 1305 transmits information between various components of the device (e.g., processor 1301, memory 1302, input / output interface 1303, and communication interface 1304); The processor 1301, memory 1302, input / output interface 1303 and communication interface 1304 are connected to each other within the device via bus 1305.

[0149] This application also provides a storage medium storing a computer program that, when executed by a processor, implements the above-described omnichannel signal flow analysis method.

[0150] Memory, as a non-transitory storage medium, can be used to store non-transitory software programs and non-transitory computer-executable programs. Furthermore, memory may include high-speed random access memory, and may also include non-transitory memory, such as at least one disk storage device, flash memory device, or other non-transitory solid-state storage device. In some embodiments, memory may optionally include memory remotely located relative to the processor, and these remote memories can be connected to the processor via a network. Examples of such networks include, but are not limited to, the Internet, intranets, local area networks, mobile communication networks, and combinations thereof.

[0151] The method, apparatus, device, and storage medium for full-pathway signal flow analysis proposed in this application acquire the expression matrix and regulatory network data of a single-cell transcriptome, and determine the transcription factor activity matrix based on the expression matrix and regulatory network data; identify at least one experimental biological pathway based on the regulatory network data, the experimental biological pathway including at least receptors, transcription factors, and target genes; traverse the experimental biological pathway, and determine the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene based on the expression matrix, regulatory network data, and transcription factor activity matrix; calculate the significance verification result for the experimental biological pathway based on the signal regulation intensity and signal flow probability intensity, and determine the target signal flow pathway based at least on the significance verification result. This application embodiment constructs a complete signal transduction pathway including receptors, transcription factors, and target genes by quantitatively calculating transcription factor activity, supplementing the intermediate links of intracellular signal transduction, restoring the entire process from extracellular recognition to intracellular transcriptional regulation, and obtaining analysis results that more closely match real biological processes, improving the completeness of cell communication analysis pathways. Simultaneously, by quantitatively calculating the signal flow probability intensity from receptor to transcription factor and the signal regulation intensity from transcription factor to target gene, layer-by-layer verification of signal transduction is achieved. This effectively eliminates false-positive communication events where only receptor-ligand co-expression exists but downstream transcriptional regulation is not actually triggered, improving the accuracy of communication inference. Furthermore, quantitative calculations of signal intensity based on transcriptome expression characteristics do not rely on fixed cell spatial location information, breaking free from the constraints of static tissue spatial assumptions. This allows for adaptation to the physiological characteristics of dynamic cell flow and dispersed ligand distribution in fluid environments such as blood and lymph, enhancing the rationality and applicability of cell sample communication analysis results in fluid environments. Finally, after statistical significance verification, the output signal flow pathways possess statistical significance support, improving the completeness, accuracy, and scenario adaptability of cell communication analysis results.

[0152] The embodiments described in this application are for the purpose of more clearly illustrating the technical solutions of the embodiments of this application, and do not constitute a limitation on the technical solutions provided by the embodiments of this application. As those skilled in the art will know, with the evolution of technology and the emergence of new application scenarios, the technical solutions provided by the embodiments of this application are also applicable to similar technical problems.

[0153] Those skilled in the art will understand that the technical solutions shown in the figures do not constitute a limitation on the embodiments of this application, and may include more or fewer steps than shown, or combine certain steps, or different steps.

[0154] The device embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs.

[0155] Those skilled in the art will understand that all or some of the steps in the methods disclosed above, as well as the functional modules / units in the systems and devices, can be implemented as software, firmware, hardware, or suitable combinations thereof.

[0156] The terms “first,” “second,” “third,” “fourth,” etc. (if present) in the specification and accompanying drawings of this application are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be interchanged where appropriate so that the embodiments of this application described herein can be implemented in orders other than those illustrated or described herein. Furthermore, the terms “comprising” and “having,” and any variations thereof, are intended to cover non-exclusive inclusion; for example, a process, method, system, product, or apparatus that comprises a series of steps or units is not necessarily limited to those steps or units explicitly listed, but may include other steps or units not explicitly listed or inherent to such processes, methods, products, or apparatus.

[0157] It should be understood that in this application, "at least one (item)" means one or more, and "more than" means two or more. "Or" is used to describe the relationship between related objects, indicating that three relationships can exist. For example, "A or B" can mean: only A exists, only B exists, and both A and B exist simultaneously, where A and B can be singular or plural. The character " / " generally indicates that the preceding and following related objects are in an "or" relationship. "At least one (item) of the following" or similar expressions refer to any combination of these items, including any combination of single or plural items. For example, at least one (item) of a, b, or c can mean: a, b, c, "a and b", "a and c", "b and c", or "a and b and c", where a, b, and c can be single or multiple.

[0158] In the several embodiments provided in this application, it should be understood that the disclosed apparatus and methods can be implemented in other ways. For example, the apparatus embodiments described above are merely illustrative; for instance, the division of the units described above is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be through some interfaces; the indirect coupling or communication connection between apparatuses or units may be electrical, mechanical, or other forms.

[0159] The units described above as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.

[0160] Furthermore, the functional units in the various embodiments of this application can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software functional unit.

[0161] If the integrated unit is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, or all or part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes multiple instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods of the various embodiments of this application. The aforementioned storage medium includes various media capable of storing programs, such as USB flash drives, portable hard drives, read-only memory (ROM), random access memory (RAM), magnetic disks, or optical disks.

[0162] The preferred embodiments of the present application have been described above with reference to the accompanying drawings, but this does not limit the scope of the claims of the present application. Any modifications, equivalent substitutions, and improvements made by those skilled in the art without departing from the scope and substance of the embodiments of the present application shall be within the scope of the claims of the present application.

Claims

1. A method for analyzing signal flow across all pathways, characterized in that, include: The single-cell transcriptome expression matrix and prior regulatory network data of the sample are obtained, and the transcription factor activity matrix is ​​determined based on the single-cell transcriptome expression matrix and the prior regulatory network data. Based on a priori pathway database, at least one experimental biological pathway is identified from the single-cell transcriptome expression matrix, and the experimental biological pathway includes at least a ligand, a receptor, a transcription factor, and a target gene. By traversing the biological pathways to be tested, and based on the single-cell transcriptome expression matrix, the prior regulatory network data, and the transcription factor activity matrix, the pathway activation parameters of the biological pathways to be tested in single cells are determined. The pathway activation parameters include the signal flow probability intensity from the receptor to the transcription factor and the signal regulation intensity from the transcription factor to the target gene; Based on the signal modulation intensity and the signal flow probability intensity, the biological pathway to be tested is validated for significance in order to determine the target signal flow pathway.

2. The method for analyzing the entire signal flow path according to claim 1, characterized in that, The single-cell transcriptome expression matrix contains the expression levels of different genes in each cell; the prior regulatory network data contains the regulatory weights between different genes and each transcription factor; the determination of the transcription factor activity matrix based on the single-cell transcriptome expression matrix and the prior regulatory network data includes: Obtain the intersection genes between the single-cell transcriptome expression matrix and the prior regulatory network data, construct an intersection expression matrix based on the expression levels of the intersection genes, and construct an intersection weight matrix based on the regulatory weights of the intersection genes. The transcription factor activity matrix is ​​calculated based on the transpose of the intersection expression matrix and the intersection weight matrix; the transcription factor activity matrix contains the activity value of each transcription factor in different cells.

3. The method for analyzing the entire signal flow path according to claim 1, characterized in that, The method of identifying at least one experimental biological pathway from the single-cell transcriptome expression matrix based on a priori pathway database includes: Based on a priori pathway database, at least one initial biological pathway is obtained, wherein the initial biological pathway includes at least a ligand, a receptor, a transcription factor, and a target gene; From the single-cell transcriptome expression matrix, obtain the expression levels of the ligands, receptors, transcription factors, and target genes in the initial biological pathway in different cells; Determine whether the proportion of cells with non-zero expression levels in the ligand, receptor, transcription factor, and target gene is greater than or equal to a preset threshold. If so, the initial biological pathway is taken as the biological pathway to be tested.

4. The method for analyzing the entire signal flow path according to claim 1, characterized in that, The determination of pathway activation parameters for the biological pathway under investigation in single cells based on the single-cell transcriptome expression matrix, the prior regulatory network data, and the transcription factor activity matrix includes: Based on the single-cell transcriptome expression matrix, the mean expression levels of the receptor, transcription factor, and target gene in each cell were determined, and the receptor half-saturation coefficient, transcription factor half-saturation coefficient, and target gene half-saturation coefficient were obtained. For the receptor and the transcription factor, the probability intensity of the signal flow is determined based on the single-cell transcriptome expression matrix, the prior regulatory network data, the receptor half-saturation coefficient, the transcription factor half-saturation coefficient, and the transcription factor activity matrix; The signal regulation intensity is determined for the transcription factor and the target gene based on the single-cell transcriptome expression matrix, the prior regulatory network data, the half-saturation coefficient of the transcription factor, the half-saturation coefficient of the target gene, and the transcription factor activity matrix.

5. The method for analyzing the entire signal flow path according to claim 4, characterized in that, The determination of the signal flow probability intensity based on the single-cell transcriptome expression matrix, the prior regulatory network data, the receptor half-saturation coefficient, the transcription factor half-saturation coefficient, and the transcription factor activity matrix includes: Obtain the receptor expression level corresponding to the receptor in the single-cell transcriptome expression matrix, and the transcription factor activity value corresponding to the transcription factor in the transcription factor activity matrix; The receptor activation percentage is obtained by calculating the proportion of the receptor expression level in the sum of the receptor expression level and the receptor half-saturation coefficient. The transcription factor activation percentage is obtained by calculating the proportion of the transcription factor activity value in the sum of the transcription factor activity value and the transcription factor half-saturation coefficient. The physical activation probability of the receptor activating the transcription factor is obtained based on the receptor activation percentage and the transcription factor activation percentage. The signal flow probability intensity is obtained based on the physical activation probability and the regulatory weights of the transcription factors in the prior regulatory network data.

6. The method for analyzing the entire signal flow path according to claim 4, characterized in that, The determination of the signal regulation intensity based on the single-cell transcriptome expression matrix, the prior regulatory network data, the transcription factor half-saturation coefficient, the target gene half-saturation coefficient, and the transcription factor activity matrix includes: Obtain the expression levels of the target genes corresponding to the target genes in the single-cell transcriptome expression matrix and the activity values ​​of the transcription factors corresponding to the transcription factors in the transcription factor activity matrix; The second-order activity value is obtained from the transcription factor activity value, and the second-order half-saturation coefficient is obtained from the transcription factor half-saturation coefficient. The target gene intensity ratio is obtained by calculating the proportion of the target gene expression level in the sum of the target gene expression level and the target gene half-saturation coefficient; and the transcription factor intensity ratio is obtained by calculating the proportion of the second-order activity value in the sum of the second-order activity value and the second-order half-saturation coefficient. The physical response intensity of the transcription factor activating the target gene is obtained based on the target gene intensity ratio and the transcription factor intensity ratio. The signal regulation intensity is obtained based on the physical response intensity and the regulatory weights corresponding to the transcription factors in the prior regulatory network data.

7. The method for analyzing the entire signal flow path according to claim 1, characterized in that, The significance verification of the biological pathway under test based on the signal modulation intensity and the signal flow probability intensity includes: The pathway activation parameters of a single cell are obtained by multiplying the signal modulation intensity and the signal flow probability intensity. The true observed signal intensity of the biological pathway under test is obtained by calculating the average value of the pathway activation parameters in all cells within a specific cell population. Randomly shuffle the pairings of the receptor and the transcription factor or the transcription factor and the target gene in the biological pathway to be tested to obtain the corresponding randomized biological pathway; For each of the random biological pathways, a corresponding virtual signal intensity is calculated. The significance P-value of the biological pathway to be tested is calculated based on the actual observed signal intensity and the virtual signal intensity. The significance verification result is obtained based on the significance P-value.

8. A full-path signal flow analysis device, characterized in that, include: Activity matrix calculation module: used to obtain the single-cell transcriptome expression matrix and prior regulatory network data of the sample, and to determine the transcription factor activity matrix based on the single-cell transcriptome expression matrix and the prior regulatory network data; Initial Pathway Filtering Module: Used to identify at least one experimental biological pathway from the single-cell transcriptome expression matrix based on a priori pathway database, wherein the experimental biological pathway includes at least a ligand, receptor, transcription factor, and target gene; Quantitative calculation module: used to traverse the biological pathway to be tested, and determine the pathway activation parameters of the biological pathway to be tested in a single cell based on the single-cell transcriptome expression matrix, the prior regulatory network data and the transcription factor activity matrix; The pathway activation parameters include the signal flow probability intensity from the receptor to the transcription factor and the signal regulation intensity from the transcription factor to the target gene; Verification output module: used to perform saliency verification on the biological pathway to be tested based on the signal modulation intensity and the signal flow probability intensity, so as to determine the target signal flow pathway.

9. An electronic device, characterized in that, The electronic device includes a memory and a processor, the memory storing a computer program, and the processor executing the computer program to implement the full-path signal flow analysis method according to any one of claims 1 to 7.

10. A storage medium storing a computer program, characterized in that, When the computer program is executed by the processor, it implements the full-path signal flow analysis method according to any one of claims 1 to 7.

Images