A method and system for joint analysis of gene-level DNA methylation, transcriptome and proteome based on snakemake framework and application thereof

CN121641196BActive Publication Date: 2026-08-11SHANGHAI OE BIOTECH CO LTD
View PDF 2 Cites 0 Cited by

Patent Information

Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2024-08-28
Publication Date
2026-08-11

AI Technical Summary

Technical Problem

[0005](1)DNA甲基化分析局限于位点水平,缺乏对基因甲基化水平的整体评估,难以聚焦关键基因的调控模式

Benefits of technology

[0063] The analytical method in this invention provides a holistic assessment of gene methylation levels, allowing the analytical results to focus on the regulatory patterns of key genes.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN121641196B_ABST
    Figure CN121641196B_ABST
Patent Text Reader

Abstract

This invention discloses a method for joint analysis of DNA methylation, transcriptome, and proteome at the gene level based on the Snakemake framework, comprising the following steps: preprocessing of raw DNA methylation and / or transcriptome and / or proteome data; gene methylation analysis; correlation analysis and / or differential intersection analysis and / or pathway enrichment joint analysis based on gene methylation analysis data, preprocessed transcriptome data, and preprocessed proteome data; and finally, organizing the analysis results and generating a visualization report. This invention provides comprehensive results, involving quantification of gene methylation, differential analysis, pathway enrichment analysis, and correlation analysis, differential intersection analysis, and pathway enrichment joint analysis integrating DNA methylation and / or transcriptome and / or proteome data; it achieves automatic organization, verification, visualization, and report generation of analysis results; all operation steps are traceable and generate corresponding analysis log records. This invention also discloses related systems and applications.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention belongs to the field of bioinformatics technology and relates to a method, system and application of gene-level DNA methylation, transcriptome and proteome joint analysis based on the Snakemake framework. Background Technology

[0002] DNA methylation is the process of adding a methyl group (-CH3) to a DNA molecule. Specifically, this process is catalyzed by DNA methyltransferases and typically occurs at the 5th carbon atom of cytosine (C) in CpG dinucleotides (DNA sequences linked by C and G, known as CpG islands), ultimately forming 5-methylcytosine. As a crucial mechanism for gene expression regulation, high methylation in regions such as promoters usually inhibits gene expression, while low methylation or demethylation can activate it. DNA methylation has been shown to play a key role in various diseases, including genomic imprinting, X chromosome inactivation, cell differentiation and development, genome stability, and cancer. As a relatively inexpensive and convenient method for DNA methylation detection, methylation chips provide stable and accurate measurements of DNA methylation, laying the foundation for widespread research on DNA methylation regulation.

[0003] With the development of sequencing technology, the cost of next-generation sequencing has been continuously decreasing and the accuracy has been continuously improving, leading to its widespread application in scientific research. Meanwhile, advancements in protein quantification technology have made proteomics quantitative detection increasingly widely used in the field of transcriptional regulation. Integrating DNA methylation data with transcriptome sequencing and proteomics quantitative analysis results is an important approach to analyzing and elucidating the transcriptional regulatory role of DNA methylation.

[0004] Existing combined analysis of DNA methylation data with transcriptomic and proteomic data has the following limitations:

[0005] (1) DNA methylation analysis is limited to the site level and lacks an overall assessment of gene methylation levels, making it difficult to focus on the regulatory patterns of key genes.

[0006] (2) Incomplete results presentation: The analysis results are too simple, the data mining is not in-depth, and there is a lack of corresponding visualization content. Summary of the Invention

[0007] To address the shortcomings of existing technologies, the present invention aims to provide a method, system, and application for combined analysis of gene-level DNA methylation, transcriptomics, and proteomics based on the snakemake framework.

[0008] To achieve the objectives of this invention, the technical solution adopted is as follows:

[0009] This invention provides a method for combined analysis of gene-level DNA methylation, transcriptomics, and proteomics based on the snakemake framework, comprising the following steps:

[0010] Step 1: Document Preparation Steps

[0011] Prepare raw analysis reports for DNA methylation and / or transcriptomics and / or proteomics, and prepare config files for subsequent data processing and analysis.

[0012] The config file includes: the original report storage location and the corresponding sample analysis name and group name, group information for differential analysis, differential fold and significance parameters, and GO and KEGG database annotation information;

[0013] The original analysis report and the config file are used for subsequent data processing and analysis.

[0014] Step 2: Preprocess DNA methylation and / or transcriptome and / or proteome data:

[0015] Preprocessing steps for DNA methylation data:

[0016] For the raw data from the methylation chip, the site Beta methylation level, the positional relationship between the gene and the nearby CpG islands were extracted, and site functional regions were divided and site-gene mapping was performed;

[0017] The Beta methylation level refers to the quantitative methylation value of each site in the original DNA methylation microarray report. The positional relationship between the gene and nearby CpG islands includes Island region (CpG island region), Shore region (2000bp region near CpG island), Shelf region (2000bp region near Shore), and OpenSea region (region other than Island, Shore, and Shelf). Functional region division refers to dividing the site into different gene regions, including promoter regions, gene regions covering promoter regions, gene body regions, and other functional regions defined based on UCSC RefGene annotation information. Site-gene mapping refers to mapping the DNA methylation sites contained in different functional regions to their corresponding genes.

[0018] Preprocessing steps for transcriptome data:

[0019] For data from the transcriptome, the results of gene quantification, differential analysis, and enrichment analysis were extracted;

[0020] The gene quantification results include gene IDs and their expression levels in each sample; the differential analysis results include gene IDs and their P-values ​​for differential analysis between different experimental and control groups, the corrected P-values ​​(Q-values) obtained after multiple test correction of the differential P-values, the fold change (FC) value, and the log2 fold change (log2FC) value; the enrichment analysis results include pathway IDs, pathway names, pathway description information, enrichment significance P-values, enrichment scores, the gene ratio of the enriched pathway, and the specific differentially expressed genes contained in the pathway.

[0021] Preprocessing steps for proteomic data:

[0022] For data from the proteome, extract the results of protein quantification, differential analysis, and enrichment analysis;

[0023] The protein quantification results include the protein ID and its expression level in each sample; the differential analysis results include the protein ID and its P-values ​​for differential analysis between different experimental groups and control groups, the corrected P-value (Q-value) obtained after multiple test correction of the differential P-values, the fold change (FC) value, and the log2 fold change (log2FC) value; the enrichment analysis results include the pathway ID, pathway name, pathway description information, enrichment significance P-value, enrichment score, protein ratio of the enriched pathway, and the specific differentially expressed proteins contained in the pathway.

[0024] Step 3: Gene DNA methylation analysis steps:

[0025] For the methylation data after functional region division, methylation sites that are not CpG islands are filtered out based on site annotation information;

[0026] For each functional region, based on the site-gene mapping relationship, methylation sites were mapped to genes to obtain the site set of each gene in each functional region. The average methylation level of each gene on the site set of each functional region was calculated. For each gene, the t-test was used to calculate the significance of the difference in methylation level between different experimental groups and the control group in each functional region. The average methylation level of each experimental group and the control group in each functional region was calculated. The difference in the average methylation level between each experimental group and the control group in each functional region was calculated as the deltaBeta value. The ratio of the average methylation level between each experimental group and the control group in each functional region was calculated as the fold change (FC) value. The log2 transformed value of the fold change (FC) value between each experimental group and the control group in each functional region was calculated as the log2FC value.

[0027] In one specific implementation, the default screening criteria for the significantly differentially methylated genes are a P-value < 0.05 and |deltaBeta| > 0.1;

[0028] After screening, GO and KEGG pathway enrichment analysis was performed on methylated genes.

[0029] The CpG islands refer to DNA regions rich in CpG dinucleotides, which play an important role in gene regulation. In step three, methylation sites located in CpG islands are determined based on site annotation information. Methylation sites that are not CpG islands in different functional regions are filtered out, allowing for focused analysis of regions that have a more significant impact on gene expression.

[0030] Step 4: Perform correlation analysis and / or differential intersection analysis and / or enrichment pathway joint analysis based on gene methylation analysis data, preprocessed transcriptome data, and preprocessed proteome data;

[0031] Steps for correlation analysis between DNA methylation levels and transcription and protein expression levels:

[0032] For the DNA methylation level and transcription and protein expression level of genes, the average methylation, transcription and / or protein expression level of each gene and protein is calculated, and a scatter plot of the overall gene correlation is drawn; the correlation coefficient of each gene expression is calculated; the distribution density plot of the correlation coefficient in the range of 0 to 1 is drawn; significantly correlated genes are screened, and the correlation scatter plot of single genes is drawn; finally, pathway enrichment analysis is performed on significantly correlated genes.

[0033] The correlation coefficients include two calculation methods: Spearman and Pearson. Spearman correlation coefficients are calculated based on the order of sample expression levels and have no specific requirements for the data distribution, making them more widely applicable and the default method. Pearson correlation coefficients require the data to follow a normal distribution, and are therefore typically used when the expression values ​​of genes / proteins in each sample conform to a normal distribution. Significant correlation is defined as a correlation coefficient with a p-value < 0.05. The pathway enrichment analysis performed on significantly correlated genes includes GO and KEGG enrichment analyses.

[0034] Differential intersection analysis of DNA methylation with transcriptome and proteome:

[0035] For the DNA methylation and transcriptional / protein differential analysis results, a violin plot of the degree of transcriptional and protein expression difference was drawn after grouping according to the significance of DNA methylation differences to show the significant changes in the degree of transcriptional and protein difference with the degree of methylation difference. A quadrant plot of the degree of DNA methylation difference with transcriptome and proteome difference was also drawn to further demonstrate the impact of changes in gene methylation difference on transcriptional and protein expression differences, and to visually reflect the scale of the intersection between different methylation difference screening results and transcriptome and proteome difference screening results. Intersection analysis was performed on genes with significant differences in each omics, and the resulting genes reflected the direct impact of significant methylation differences on significant differences in transcriptional and protein expression. Pathway enrichment analysis was performed on the intersecting genes to obtain the significant impact of coordinated changes in gene methylation and transcriptome and / or proteome expression on biological functional pathways. Furthermore, a multi-omics expression profile of the intersecting genes was drawn. If multiple comparison groups of gene methylation and transcriptome and / or proteome were involved, a bubble plot of the correlation of the intersecting genes in each comparison group was drawn to specifically show the coordinated changes in the expression levels of the intersecting genes in each omics.

[0036] In the differential intersection analysis, the DNA methylation difference significance grouping includes hypomethylation, no significant change, and hypermethylation; hypomethylation refers to the group with P value < 0.05 and deltaBeta < -0.1, no significant change refers to the group with P ≥ 0.05 or -0.1 ≤ deltaBeta ≤ 0.1, and hypermethylation refers to the group with P value < 0.05 and deltaBeta > 0.1; the pathway enrichment analysis of significantly related genes includes GO and KEGG enrichment analysis;

[0037] Joint analysis of DNA methylation and differential gene enrichment pathways in the transcriptome and proteome:

[0038] For the enrichment pathways of differentially expressed genes and / or proteins in DNA methylation, transcriptomics, and proteomics, scatter plots and bar charts of pathway significance were plotted to show the significance P-values ​​of pathways included in the enrichment results of different omics pathways (P-values ​​in the figures were converted to -log10). Intersection analysis was performed on the significantly enriched pathways in different omics data. For the intersection pathways, pathway bubble charts were plotted to show the common characteristics of pathway significance and enrichment scores. For the intersection pathways, pathway stacking plots were plotted to show the ranking results after the gene ratios of each omics pathway were superimposed. Finally, the expression profiles of significant genes and proteins in each omics pathway were plotted to reveal the coordinated changes in the expression of key genes contained in the significantly consistent pathways in each omics.

[0039] In the combined enrichment pathway analysis, the enriched pathways are distinct from the enrichment analysis results performed in correlation analysis and differential intersection analysis. They refer to the pathways and their enrichment results obtained from the KEGG enrichment analysis of the individual analysis results of gene methylation, transcriptomics, and proteomics. The enrichment results include, but are not limited to, pathway ID, pathway name, pathway description information, enrichment significance P-value, enrichment score, gene / protein ratio of the enriched pathway, and specific gene / protein differentials contained in the pathway.

[0040] The gene / protein ratio is obtained by calculating "number of differentially expressed genes or proteins in the pathway" / "total number of differentially expressed genes or proteins"; the enrichment score is obtained by calculating the ratio of the two parts, where the numerator is the gene / protein ratio and the denominator is obtained by calculating "number of genes or proteins expressed in the pathway" / "total number of genes or proteins expressed".

[0041] Step 5: Web-based report preparation steps:

[0042] Based on the analysis results, a web-based report can be generated with one click, providing a combined analysis of DNA methylation, transcriptomics, and proteomics. The web-based report summarizes the entire analysis results, describes each analysis step, displays corresponding charts, and includes pop-up help documentation. The web-based report also includes internal shortcuts and links to analysis method introductions / external websites, enabling quick navigation within the web-based system and rapid access to online resources.

[0043] In this invention, step one is a document preparation step, step two is a data preprocessing step, step three is a gene DNA methylation analysis step, step four is a combined analysis step, and step five is a report preparation step.

[0044] In a preferred embodiment of the present invention, the config file in the file preparation step includes: the original report storage location and the corresponding sample analysis name and group name, grouping information for differential analysis, differential fold and significance parameters, GO and KEGG database annotation information, etc.

[0045] In another preferred embodiment of the present invention, in the preprocessing step of the DNA methylation data, the functional region division includes, but is not limited to, promoter regions, gene regions covering promoter regions and gene body regions, and functional regions can be customized based on UCSC RefGene annotation information; the site-gene mapping is functional region specific, therefore it is necessary to map the site-gene correspondence according to different functional region division standards.

[0046] The promoter region refers to a DNA sequence located at the 5' end of a gene, which is the binding site of RNA polymerase and other transcription factors, initiating the transcription process of the gene; the gene body region refers to all sequences from the transcription start site (TSS) to the transcription termination site (TTS), including exons, introns, and some regulatory sequences; the custom functional regions include the UTR5 region (the transcription and post-transcriptional regulatory region between the TSS and the first exon), the UTR3 region (the transcription and post-transcriptional regulatory region upstream of the TTS), the first exon region (the region where the first exon is located), the TSS1500 region (the region from the TSS to its transcription upstream of 1500 bp), the TSS200 region (the region from the TSS to its transcription upstream of 200 bp), etc.

[0047] The site-gene mapping refers to mapping the DNA methylation sites contained in different functional regions to their corresponding genes.

[0048] In another preferred embodiment of the present invention, in the gene DNA methylation analysis step, the t-test difference test results for calculating gene methylation levels in each comparison group include a significance P-value, the difference between the mean gene methylation values ​​of the two groups (deltaBeta), the fold change (FC), and the log2 logarithmic transformation value of the fold change (log2FC). The significance Q-value is obtained by correcting the significance P-value using the Benjamini & Hochberg multiple test. Differential genes are screened according to the threshold filled in by the config, with the default P-value < 0.05 and deltaBeta > 0.1. The GO and KEGG pathway enrichment analysis results for differentially methylated genes include an enrichment significance result table, a GO enrichment significant pathway bar chart, and a KEGG enrichment significant pathway bubble chart.

[0049] The deltaBeta, FC, and log2FC represent the degree of difference between the experimental and control groups in gene methylation level Beta values ​​from different dimensions. The significance P-value reflects the reliability of the difference test results; the smaller the P-value, the higher the reliability of the difference test results. The commonly used significance index is P<0.05. The Benjamini & Hochberg multiple test correction for the significance P-value refers to correcting the false positive P-value <0.05 results caused by differential analysis of a large number of genes based on the P-value results. The resulting Q-value can achieve more stringent significance filtering.

[0050] In another preferred embodiment of the present invention, in the correlation analysis step between DNA methylation level and transcription and protein expression levels, the overall gene correlation scatter plot includes linear fitting of gene points using the lm method in R language, and provides a goodness-of-fit correction value, a significance P-value for the linear fitting F-test, and a Pearson correlation coefficient, illustrating the relationship between gene methylation level and transcriptional expression level and / or protein expression level; the calculation of the correlation coefficient of each gene expression can be performed using the Pearson or Spearman method, calculating and illustrating the correlation between the methylation level of each gene and the transcriptional or protein expression level; significant correlation refers to a significance P-value of the correlation coefficient < 0.05; the pathway enrichment analysis of significantly correlated genes includes a KEGG and GO enrichment significance result table, a GO enrichment significant pathway bar chart, and a KEGG enrichment significant pathway bubble chart;

[0051] The Spearman or Pearson method is used for correlation calculation. Spearman calculates correlation based on the order of sample expression levels and has no requirements on the specific data distribution, so it is more widely applicable and is the default method for calculating correlation coefficients. Pearson correlation coefficient requires the data to meet a normal distribution when calculating it. Therefore, as a preferred method, it is usually used when the expression values ​​of each sample of a gene / protein meet a normal distribution.

[0052] The 1m method refers to a method for linear regression analysis in R language. It can achieve a linear fit of y = ax + b for given independent variables (x, such as the methylation beta value of a single gene in each sample) and dependent variables (y, such as the transcriptome expression value of a single gene in each sample), and return analytical results such as the goodness-of-fit correction value and the significance p-value of the linear fit F test.

[0053] In another preferred embodiment of the present invention, in the differential intersection analysis of DNA methylation with the transcriptome and proteome, the significant differences in DNA methylation are grouped into hypomethylation (P < 0.05 and deltaBeta < -0.1), no significant change (P ≥ 0.05 or -0.1 ≤ deltaBeta ≤ 0.1), and hypermethylation (P < 0.05 and deltaBeta > 0.1); the differential degree violin plot includes a violin plot displaying the log2FC values ​​of the transcriptome and proteome of the above three groups of genes; the differential degree quadrant plot includes DNA... Scatter plots and quadrant divisions were performed for methylation deltaBeta and log2FC of the transcriptome and proteome. Significantly differentially expressed genes were also marked with color based on the omics differential screening criteria in the config. The first quadrant represents upregulated methylation and upregulated transcription / protein expression; the second quadrant represents downregulated methylation and upregulated transcription / protein expression; the third quadrant represents downregulated methylation and downregulated transcription / protein expression; and the fourth quadrant represents upregulated methylation and downregulated transcription / protein expression. All plots in this step were created using the R language packages ggplot2, pheatmap, and ggpubr.

[0054] In another preferred embodiment of the present invention, in the joint analysis of DNA methylation and enrichment pathways of transcriptomics and proteomics, the intersection plot, scatter plot, bar chart, bubble chart and stacking plot are drawn using the ggplot2 package of R language; the multi-omics expression profile of genes in the intersection pathway is drawn using the ComplexHeatmap package of R language.

[0055] This invention also proposes a multi-omics joint analysis system, which can be used to implement the above-mentioned joint analysis method, including: a file data input module, a data preprocessing module, a methylation analysis module, a multi-omics joint analysis module, and a visualization generation and report output module;

[0056] The file data input module is used to input prepared raw analysis reports of DNA methylation, transcriptomics, and proteomics, as well as config configuration files;

[0057] The data preprocessing module is used to preprocess DNA methylation, transcriptome, and proteome data;

[0058] The methylation analysis module is used to analyze the DNA methylation level of genes and perform differential analysis and pathway enrichment analysis.

[0059] The multi-omics joint analysis module is used to perform correlation analysis, differential intersection analysis and / or enrichment pathway joint analysis on DNA methylation, transcriptome and proteome data.

[0060] The visualization generation and report output module is used to visualize the analysis results and generate the final report.

[0061] This invention also discloses the application of the above-mentioned combined analysis method, or multi-omics combined analysis system, in gene-level DNA methylation and transcriptome and / or proteome combined analysis, site-level DNA methylation and transcriptome and / or proteome combined analysis, screening of genes associated with methylation and transcriptome and / or proteome expression profiles, screening of genes with differential intersection of methylation and transcriptome and / or proteome, and screening of key regulatory pathways of methylation and transcriptome and / or proteome.

[0062] The beneficial effects of this invention include:

[0063] The analytical method in this invention provides a holistic assessment of gene methylation levels, allowing the analytical results to focus on the regulatory patterns of key genes.

[0064] The results are comprehensive, including analysis of gene DNA methylation levels, as well as correlation analysis, differential intersection analysis, and pathway analysis of combined transcriptome and proteome data. The combined analysis modules provide detailed data analysis and rich plotting displays.

[0065] The highly modular nature of the analysis content greatly reduces process maintenance costs and enhances process scalability, which is beneficial for subsequent updates and iterations.

[0066] It boasts high operational efficiency, and its modular analysis process allows for detailed breakdown of the analysis code. This facilitates the separate delivery and execution of each module's sub-contents based on the snakemake framework, significantly improving analysis efficiency. It is more than 10 times faster than single-threaded processing, and this efficiency improvement is even more pronounced when there are many combinations of experimental and control groups.

[0067] It reduces manual time consumption, has strong compatibility in the data preprocessing section, can adapt to various data types and formats, and supports personalized adjustment of analysis thresholds, plotting parameters, etc., meeting more than 90% of analysis scenarios; it can complete a complex whole analysis process with only simple file configuration; compared with the previous analysis process, it reduces manual time consumption by at least 50%.

[0068] It automatically organizes all analysis results, and automatically summarizes, statistically analyzes, visualizes, and logically categorizes the results after each analysis step is completed. The result files can be directly used to generate web-based reports.

[0069] All operational steps are traceable, facilitating error lookup. If an error is analyzed, corresponding error log information will be available. Attached Figure Description

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

[0071] Figure 1 This is a schematic diagram of the process of the present invention.

[0072] Figure 2 This is a schematic diagram illustrating the DNA methylation data analysis content of the present invention.

[0073] Figure 3 This is a schematic diagram illustrating the multi-omics joint analysis content of the present invention.

[0074] Figure 4 This is an example of a scatter plot showing the overall correlation of the present invention.

[0075] Figure 5 This is an example diagram showing the correlation coefficient distribution of the present invention.

[0076] Figure 6 This is an example diagram of the joint distribution of correlation coefficients in this invention.

[0077] Figure 7 This is an example of a scatter plot illustrating the single-gene expression association of the present invention.

[0078] Figure 8 This is an example of a bar chart illustrating the significant pathway enrichment of GO according to the present invention.

[0079] Figure 9 This is an example bubble diagram of the significant pathway enrichment of KEGG according to the present invention.

[0080] Figure 10 This is an example diagram of the difference degree violin diagram of the present invention.

[0081] Figure 11 This is an example diagram of the difference quadrant diagram of the present invention.

[0082] Figure 12 This is an example diagram of the multi-omics differentially significant gene intersection analysis of the present invention.

[0083] Figure 13 This is an example diagram of the intersection gene multi-omics expression profile (single gene) of the present invention.

[0084] Figure 14 This is an example diagram of the intersection gene multi-omics expression profile of the present invention.

[0085] Figure 15 This is an example bubble chart illustrating the correlation between the intersection genes of the present invention in each comparison group.

[0086] Figure 16 This is an example diagram of the intersection analysis of multi-omics differential base enrichment pathways in this invention.

[0087] Figure 17 This is an example of a scatter plot of the path P-values ​​for the present invention.

[0088] Figure 18 This is an example of a bar chart showing the pathway P-value of the present invention.

[0089] Figure 19 This is an example diagram of the pathway stacking diagram of the present invention.

[0090] Figure 20 This is an example diagram of the pathway bubble diagram of the present invention.

[0091] Figure 21 This is an example diagram of the expression profiles of genes and proteins with significant differences in each omics contained in the pathway of the present invention. Detailed Implementation

[0092] The present invention will be further described in detail below with reference to the specific embodiments and accompanying drawings. Except for the contents specifically mentioned below, the processes, conditions, and experimental methods for implementing the present invention are all common knowledge and general knowledge in the art, and the present invention does not have any particular limitations.

[0093] This invention provides a method for joint analysis of DNA methylation, transcriptome, and proteome at the gene level based on the Snakemake framework, comprising the following steps: preprocessing of raw DNA methylation and / or transcriptome and / or proteome data; gene methylation analysis; correlation analysis and / or differential intersection analysis and / or pathway enrichment joint analysis based on gene methylation analysis data, preprocessed transcriptome data, and preprocessed proteome data; and finally, organizing the analysis results and generating a visualization report. This invention provides comprehensive results, involving quantification of gene methylation, differential analysis, pathway enrichment analysis, as well as correlation analysis, differential intersection analysis, and pathway enrichment joint analysis integrating DNA methylation and / or transcriptome and / or proteome data. The workflow based on the Snakemake framework enables automatic organization, verification, visualization, and report generation of analysis results. Furthermore, all operational steps are traceable and generate corresponding analysis log records to facilitate querying and analyzing error information caused by uncontrollable factors.

[0094] This invention provides a method for combined analysis of gene-level DNA methylation, transcriptomics, and proteomics based on the snakemake framework, such as... Figure 1 As shown, it includes the following steps:

[0095] Step 1: Document Preparation Steps

[0096] Prepare raw analysis reports for DNA methylation and / or transcriptomics and / or proteomics; prepare config files, which will be used for automated data quality control and subsequent data analysis.

[0097] Step 2: Preprocess DNA methylation and / or transcriptome and / or proteome data:

[0098] Preprocessing steps for DNA methylation data:

[0099] For the raw data from the methylation chip, the Beta methylation level of the site and the positional relationship between the gene and the nearby CpG islands are extracted, and site functional regions are divided and site-gene mapping is performed to achieve differential site pathway enrichment.

[0100] Preprocessing steps for transcriptome data:

[0101] For data from the transcriptome, the results of gene quantification, differential analysis, and enrichment analysis were extracted;

[0102] Preprocessing steps for proteomic data:

[0103] For data from the proteome, extract the results of protein quantification, differential analysis, and enrichment analysis;

[0104] Step 3: Gene DNA methylation analysis steps, such as... Figure 2 As shown:

[0105] For the methylation data after functional region division, methylation sites other than CpG islands were filtered out based on site annotation information; the average methylation level of each gene at available sites in each functional region was calculated based on the site-gene correspondence; the t-test difference test of gene methylation levels was performed for each comparison group; genes with significant differential methylation were screened; and GO and KEGG pathway enrichment analysis was performed on differentially methylated genes.

[0106] Step 4: Based on gene methylation analysis data, preprocessed transcriptome data, and preprocessed proteome data, perform correlation analysis and / or differential intersection analysis and / or combined enrichment pathway analysis, such as... Figure 3 As shown;

[0107] Steps for correlation analysis between DNA methylation levels and transcription and protein expression levels:

[0108] For the DNA methylation level, transcription, and protein expression level of genes, a scatter plot of overall gene correlation is drawn; the correlation coefficient of each gene expression is calculated; then a distribution plot of the correlation coefficient is drawn; significantly correlated genes are screened, and a scatter plot of single gene correlation is drawn; pathway enrichment analysis is performed on significantly correlated genes.

[0109] Differential intersection analysis of DNA methylation with transcriptome and proteome:

[0110] For the DNA methylation and transcriptional and protein differential analysis results of genes, a violin plot of the degree of difference in transcriptional and protein expression is drawn after grouping according to the significance of DNA methylation differences; a quadrant plot of the degree of difference between DNA methylation and transcriptome and proteome is drawn; intersection analysis is performed on genes with significant differences in each omics, and pathway enrichment analysis is performed on the intersection genes; multi-omics expression profiles of intersection genes are drawn; if multiple comparison groups are involved, a bubble plot of the correlation of intersection genes in each comparison group is drawn.

[0111] Steps for combined analysis of DNA methylation and differential gene enrichment pathways in transcriptomics and proteomics:

[0112] For DNA methylation and differentially enriched gene pathways in the transcriptome and proteome, pathway intersection analysis was performed, and pathway significance scatter plots and bar charts were drawn; for the intersection pathways, pathway bubble charts, pathway stacking charts, and expression profiles of significant genes and proteins in each omics contained in the pathways were drawn.

[0113] Step 5: Web-based report preparation steps:

[0114] Based on the results, a web-based report on the combined analysis of DNA methylation, transcriptomics, and proteomics can be generated with a single click. The web-based report summarizes the entire analysis results, describes each analysis step, displays corresponding charts, and includes pop-up help documentation. The web-based report also includes internal shortcuts and links to analysis method introductions / external websites, enabling quick navigation within the web-based system and rapid access to online resources.

[0115] Example 1

[0116] The following Example 1 further illustrates the application of the present invention in the combined analysis of DNA methylation and transcriptomics and / or proteomics at the gene level, but these descriptions are not intended to limit the present invention.

[0117] Accept user's DNA methylation chip data, transcriptome sequencing data, proteome detection data, and GO and KEGG database information; at the same time, fill in the config file, including: the original report storage location, sample analysis name and corresponding group name, group information used for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information.

[0118] First, the data from each omics dataset were preprocessed using Python. For methylation microarray data, the methylation sites were first divided according to the UCSC RefGene annotation information, based on functional region definitions such as promoter regions and gene regions (user-defined functional regions). Based on the region division results, the site-gene correspondences were mapped. For transcriptome and proteome data, gene (or protein) quantification, differential analysis, and enrichment analysis results were directly extracted from their standard reports.

[0119] Based on the region division and site-gene mapping results in the DNA methylation pretreatment, non-CpG island methylation sites were first filtered out. Then, according to the site-gene mapping relationship, the mean methylation value of each gene in each functional region was calculated as the gene methylation level. Next, based on the comparison group information provided by the config, the significance P-value of gene methylation t-test for each comparison group was tested. The difference in mean gene methylation value between the two groups (deltaBeta), fold change (FC), and log2 logarithmic transformation of the fold change (log2FC) were also provided. The significance Q value was obtained by correcting the significance P-value using the Benjamini & Hochberg multiple test. Furthermore, based on the difference test threshold provided in the config (default P-value < 0.05 and |deltaBeta| > 0.1), significantly differentially methylated genes were screened, and pathway enrichment analysis was performed. A table of enrichment significance results and a GO enrichment significant pathway bar chart drawn using the R language ggplot2 plotting package are provided (see reference). Figure 8 KEGG enrichment significant pathway bubble chart (reference) Figure 9 This is used to demonstrate the significance of pathway enrichment, such as the P-value.

[0120] Further correlation analysis was performed based on gene DNA methylation levels and transcription and protein expression levels. First, the mean values ​​of genes in each sample group were calculated from the quantitative results of each omics study. A scatter plot of overall gene correlation was then plotted, and the linear fitting results obtained from the lm method using R language were displayed (e.g., ...). Figure 4 This includes the goodness-of-fit correction value, the significance p-value of the F-test for linear fit, and the Pearson correlation coefficient. Further, the correlation coefficients between gene expression values ​​in each omics are calculated (Pearson or Spearman methods can be used, but considering that the sample expression value distribution of a single gene usually does not meet the normality requirement, the Spearman method is used by default), and a correlation coefficient distribution plot is drawn (see reference). Figure 5 To compare the correlation between single-gene methylation-transcriptome correlation and methylation-proteome correlation, a joint distribution plot of correlation coefficients was plotted (see reference). Figure 6The figure shows that the methylation level of genes in the third quadrant is negatively correlated with both transcriptomic and proteomic expression values, consistently demonstrating the inhibitory effect of methylation on gene expression. Further, significantly correlated genes were screened, and scatter plots of single-gene expression values ​​were drawn (see reference). Figure 7 This demonstrates the specific expression patterns of significantly related genes across omics. Pathway enrichment analysis is performed on the gene sets, providing a table of enrichment significance results and a GO enrichment bar chart (see [link to relevant documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9 This was used to display indicators such as the significance p-value of pathway enrichment. All plotting in this step was performed using R language, with the expression spectrum heatmap created using the pheatmap plotting package and the other plots created using the ggplot2 plotting package.

[0121] Next, based on the differential analysis results of gene DNA methylation with the transcriptome and proteome, a differential intersection analysis was performed. First, based on the DNA methylation differential analysis results, all genes were divided into three groups: hypomethylation (P < 0.05 and deltaBeta < -0.1), no significant change, and hypermethylation (P < 0.05 and deltaBeta > 0.1). Based on this grouping, a violin plot of the log2FC values ​​of the transcriptome and proteome was then plotted (see reference). Figure 10 This study aimed to demonstrate the impact of DNA methylation differences on downstream gene transcription and protein expression. Furthermore, differential quadrant diagrams were constructed based on DNA methylation deltaBeta and the log2FC values ​​of the transcriptome and proteome (see reference). Figure 11 Based on the differential screening criteria for each omics in the config, significantly differentially expressed genes were marked with color to demonstrate the association between gene methylation differences and transcriptomic / proteomic expression differences. Next, intersection analysis was performed on differentially expressed genes from each omics, and an intersection plot was generated (see reference). Figure 12 ) and expression profile heatmap of intersecting genes (reference) Figure 13 , 14 This displays the patterns of change in methylation, transcription, and protein expression levels of the intersecting genes. Pathway enrichment analysis is performed on the intersecting gene set, providing a table of enrichment significance results and a GO enrichment bar chart (see [link to documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9 This is to display indicators such as the significance p-value of pathway enrichment. If the analysis involves multiple comparison groups, the 10 genes with the most significant differences (mean of each omics significance p-value after -log10 transformation) are selected from each comparison group, and a bubble chart of inter-omics correlation coefficients is plotted for all comparison groups (see reference). Figure 15This was used to visualize the similarities and differences in changes in methylation, transcription, and protein expression levels across comparison groups. All plotting in this step was performed using R language, with the expression profile heatmap created using the pheatmap package and the other plots created using the ggplot2 package.

[0122] Next, based on the DNA methylation and differential gene pathway enrichment results from the transcriptome and proteome, a joint analysis of the enriched pathways was performed. First, an intersection analysis was conducted on the significant pathways, and an intersection diagram was drawn (see reference). Figure 16 Further, a scatter plot of pathway significance was drawn (see reference). Figure 17 ), bar chart (reference) Figure 18 This was used to demonstrate the differences in p-values ​​of enrichment analyses for the same pathway across different omics systems. For intersecting pathways, a pathway bubble chart was created (see reference). Figure 19 Stacked diagram (reference) Figure 20 This section provides a detailed overview of the enrichment significance (P-value), enrichment score, and gene ratio of the intersecting pathways. Furthermore, for each intersecting pathway, significantly differentially expressed genes from each omics group are extracted, and a tri-omics expression profile heatmap is generated for each pathway (see reference). Figure 21 This study aimed to reveal the coordinated changes in the expression of key genes within significantly enriched pathways across different omics. All plotting in this step was performed using R language, with the expression profile heatmap created using the ComplexHeatmap package and the remaining plots created using the ggplot2 package.

[0123] Finally, all the analysis results were compiled and arranged into different directories according to categories: the data preprocessing and gene methylation analysis results of each omics were placed in the data preprocessing directory; the gene expression correlation analysis results of each omics were placed in the correlation analysis directory; the intersection analysis results of differentially expressed genes of each omics were placed in the differential intersection analysis directory; and the joint analysis results of enriched pathways of each omics were placed in the enriched pathway analysis directory.

[0124] Based on the above directory, use a Python script to generate the corresponding web-based analysis report with a single click.

[0125] Example 2

[0126] The following Example 2 further illustrates the application of the present invention in the combined analysis of site-level DNA methylation and transcriptomics and / or proteomics, but these descriptions are not intended to limit the present invention.

[0127] Accept user's DNA methylation chip data, transcriptome sequencing data, proteome detection data, and GO and KEGG database information; at the same time, fill in the config file, including: the original report storage location, sample analysis name and corresponding group name, group information used for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information.

[0128] First, the data from each omics dataset were preprocessed using Python. For the methylation microarray data, methylation sites were first divided according to the UCSC RefGene annotation information, based on functional region definitions such as promoter regions and gene regions (user-defined functional regions). Based on the region division results, the site-gene correspondences were mapped. Significantly differentially methylated sites were screened based on the methylation site difference analysis results in the original report. Sites were then mapped back to functional genes based on the site-gene correspondences, and pathway enrichment analysis was performed based on the mapped genes. A table of enrichment significance results and a GO enrichment significant pathway bar chart drawn using the ggplot2 plotting package in R (see reference). Figure 8 KEGG enrichment significant pathway bubble chart (reference) Figure 9 This is used to demonstrate indicators such as the significance of pathway enrichment. For transcriptomic and proteomic data, gene (or protein) quantification, differential analysis, and enrichment analysis results are directly extracted according to their standard reports.

[0129] Further, based on the locus DNA methylation level and gene transcription and protein expression levels, correlation analysis was performed according to the locus-gene mapping relationship. First, the mean values ​​of loci / genes in each sample group were calculated in the quantitative results of each omics study. A scatter plot of overall gene correlation was then plotted, and the linear fitting results obtained from the lm method using R language (e.g., locus methylation level and gene expression value) were displayed. Figure 4 This includes the linear fit goodness-of-fit correction value, the significance p-value of the linear fit F-test, and the Pearson correlation coefficient. Further, the correlation coefficient between methylation values ​​and gene expression values ​​at each locus in each omics is calculated (either the Pearson or Spearman method can be used, but considering that the sample expression value distribution of a single locus / gene usually does not meet the normality requirement, the Spearman method is used by default), and a correlation coefficient distribution plot is drawn (see reference). Figure 5 To compare the correlation between single-gene methylation-transcriptome correlation and methylation-proteome correlation, a joint distribution plot of correlation coefficients was plotted (see reference). Figure 6 The figure shows that the methylation level of genes in the third quadrant is negatively correlated with both transcriptomic and proteomic expression values, consistently demonstrating the inhibitory effect of methylation on gene expression. Furthermore, significantly correlated site-gene pairs were screened, and a scatter plot of the correlation between site methylation values ​​and gene expression values ​​was plotted (see reference). Figure 7 This demonstrates the specific expression patterns of significantly related loci and genes across omics. Pathway enrichment analysis is performed on the gene set, providing a table of enrichment significance results and a GO enrichment bar chart (see [link to relevant documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9This was used to display indicators such as the significance p-value of pathway enrichment. All plotting in this step was performed using R language, with the expression spectrum heatmap created using the pheatmap plotting package and the other plots created using the ggplot2 plotting package.

[0130] Next, based on the differential analysis results of site-specific DNA methylation with the transcriptome and proteome, differential intersection analysis was performed according to the site-gene mapping relationship. First, a differential quadrant diagram was drawn based on the site-specific DNA methylation deltaBeta and the log2FC of the transcriptome and proteome (see reference). Figure 11 Based on the differential screening criteria for each omics in the config, significantly differentially expressed genes were marked with color to demonstrate the association between differential methylation sites and differential transcriptome / proteome expression. Next, differentially methylated sites were mapped to genes, and intersection analysis was performed with differentially expressed genes / proteins in the transcriptome and / or proteome, resulting in an intersection plot (see reference). Figure 12 ) and expression profile heatmap of intersection sites-gene-protein (reference) Figure 13 , 14 This demonstrates the changing patterns of methylation levels, transcription levels, and protein expression levels at the intersection sites. Pathway enrichment analysis was performed on the intersection gene set, providing a table of enrichment significance and a GO enrichment histogram of significant pathways (see [link to documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9 This is used to display indicators such as the significance p-value of pathway enrichment. If the analysis involves multiple comparison groups, the 10 most significantly different locus-gene pairs (mean of the p-values ​​after -log10 transformation) are selected from each comparison group, and a bubble chart of inter-omics correlation coefficients is plotted for all comparison groups (see reference). Figure 15 This was used to visualize the similarities and differences in changes in methylation, transcription, and protein expression levels across comparison groups. All plotting in this step was performed using R language, with the expression profile heatmap created using the pheatmap package and the other plots created using the ggplot2 package.

[0131] Next, based on the site-specific DNA methylation and differential gene pathway enrichment results from the transcriptome and proteome, a joint analysis of the enriched pathways was performed. First, an intersection analysis was conducted on the significant pathways, and an intersection diagram was drawn (see reference). Figure 16 Further, a scatter plot of pathway significance was drawn (see reference). Figure 17 ), bar chart (reference) Figure 18 This was used to demonstrate the differences in p-values ​​of enrichment analyses for the same pathway across different omics systems. For intersecting pathways, a pathway bubble chart was created (see reference). Figure 19 Stacked diagram (reference) Figure 20This step uses the ggplot2 plotting package in R to display detailed indicators such as the significance p-value, enrichment score, and gene ratio of the intersection pathways.

[0132] Finally, all the analysis results were compiled and arranged into different directories according to categories: the data preprocessing results of each omics were placed in the data preprocessing directory; the results of the site-gene / protein expression correlation analysis were placed in the correlation analysis directory; the results of the intersection analysis of differential sites / genes / proteins of each omics were placed in the differential intersection analysis directory; and the results of the joint analysis of enriched pathways of each omics were placed in the enriched pathway analysis directory.

[0133] Based on the above directory, use a Python script to generate the corresponding web-based analysis report with a single click.

[0134] Example 3

[0135] The following Example 3 further illustrates the application of the present invention in screening genes associated with methylation and transcriptome and / or proteome expression profiles, but these descriptions are not intended to limit the present invention.

[0136] Accept user's DNA methylation chip data, transcriptome sequencing data, proteome detection data, and GO and KEGG database information; at the same time, fill in the config file, including: the original report storage location, sample analysis name and corresponding group name, group information used for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information.

[0137] First, the data from each omics dataset were preprocessed using Python. For methylation microarray data, the methylation sites were first divided according to the UCSC RefGene annotation information, based on functional region definitions such as promoter regions and gene regions (user-defined functional regions). Then, based on the region division results, the site-gene correspondences were mapped. For transcriptome and proteome data, gene (or protein) quantification results were directly extracted from their standard reports.

[0138] Based on the region division and site-gene mapping results in DNA methylation pretreatment, non-CpG island methylation sites were first filtered out. Then, according to the site-gene mapping relationship, the average methylation value of each site contained in each gene in each functional region was sorted out as the gene methylation level.

[0139] Further correlation analysis was performed based on gene DNA methylation levels and transcription and protein expression levels. First, the mean values ​​of genes in each sample group were calculated from the quantitative results of each omics study. A scatter plot of overall gene correlation was then plotted, and the linear fitting results obtained from the lm method using R language were displayed (e.g., ...). Figure 4This includes the goodness-of-fit correction value, the significance p-value of the F-test for linear fit, and the Pearson correlation coefficient. Further, the correlation coefficients between gene expression values ​​in each omics are calculated (Pearson or Spearman methods can be used, but considering that the sample expression value distribution of a single gene usually does not meet the normality requirement, the Spearman method is used by default), and a correlation coefficient distribution plot is drawn (see reference). Figure 5 To compare the correlation between single-gene methylation-transcriptome correlation and methylation-proteome correlation, a joint distribution plot of correlation coefficients was plotted (see reference). Figure 6 The figure shows that the methylation level of genes in the third quadrant is negatively correlated with both transcriptomic and proteomic expression values, consistently demonstrating the inhibitory effect of methylation on gene expression. Further, significantly correlated genes were screened, and scatter plots of single-gene expression values ​​were drawn (see reference). Figure 7 This demonstrates the specific expression patterns of significantly related genes across omics. Pathway enrichment analysis is performed on the gene sets, providing a table of enrichment significance results and a GO enrichment bar chart (see [link to relevant documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9 This was used to display indicators such as the significance p-value of pathway enrichment. All plotting in this step was performed using R language, with the expression spectrum heatmap created using the pheatmap plotting package and the other plots created using the ggplot2 plotting package.

[0140] Finally, all the analysis results were compiled and arranged into different directories according to categories: the data preprocessing and gene methylation analysis results of each omics were placed in the data preprocessing directory; and the gene expression correlation analysis results of each omics were placed in the correlation analysis directory.

[0141] Based on the above directory, use a Python script to generate the corresponding web-based analysis report with a single click.

[0142] Example 4

[0143] The following Example 4 further illustrates the application of the present invention in screening genes with differential intersections between methylation and transcriptome and / or proteome, but these descriptions are not intended to limit the present invention.

[0144] Accept user's DNA methylation chip data, transcriptome sequencing data, proteome detection data, and GO and KEGG database information; at the same time, fill in the config file, including: the original report storage location, sample analysis name and corresponding group name, group information used for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information.

[0145] First, the data from each omics dataset were preprocessed using Python. For methylation microarray data, the methylation sites were first divided according to the UCSC RefGene annotation information, based on functional region definitions such as promoter regions and gene regions (user-defined functional regions). Then, based on the region division results, the site-gene correspondences were mapped. For transcriptome and proteome data, gene (or protein) quantification results were directly extracted from their standard reports.

[0146] Based on the region division and site-gene mapping results in the DNA methylation pretreatment, non-CpG island methylation sites were first filtered out. Then, according to the site-gene mapping relationship, the mean methylation value of each gene in each functional region was compiled as the gene methylation level. Next, based on the comparison group information provided by the config, the significance P-value of gene methylation t-test for each comparison group sample was tested. At the same time, the difference in mean gene methylation value between the two groups (deltaBeta), the fold change (FC), and the log2 logarithmic transformation value of the fold change (log2FC) were provided. The significance Q value was obtained by correcting the significance P value with Benjamini & Hochberg multiple test. Furthermore, based on the difference test threshold provided in the config (the default is P value < 0.05 and |deltaBeta| > 0.1), genes with significant differential methylation were screened.

[0147] Next, based on the differential analysis results of gene DNA methylation with the transcriptome and proteome, a differential intersection analysis was performed. First, based on the DNA methylation differential analysis results, all genes were divided into three groups: hypomethylation (P < 0.05 and deltaBeta < -0.1), no significant change, and hypermethylation (P < 0.05 and deltaBeta > 0.1). Based on this grouping, a violin plot of the log2FC values ​​of the transcriptome and proteome was then plotted (see reference). Figure 10 This study aimed to demonstrate the impact of DNA methylation differences on downstream gene transcription and protein expression. Furthermore, differential quadrant diagrams were constructed based on DNA methylation deltaBeta and the log2FC values ​​of the transcriptome and proteome (see reference). Figure 11 Based on the differential screening criteria for each omics in the config, significantly differentially expressed genes were marked with color to demonstrate the association between gene methylation differences and transcriptomic / proteomic expression differences. Next, intersection analysis was performed on differentially expressed genes from each omics, and an intersection plot was generated (see reference). Figure 12 ) and expression profile heatmap of intersecting genes (reference) Figure 13 , 14 This displays the patterns of change in methylation, transcription, and protein expression levels of the intersecting genes. Pathway enrichment analysis is performed on the intersecting gene set, providing a table of enrichment significance results and a GO enrichment bar chart (see [link to documentation]). Figure 8 ), bubble chart of significant KEGG enrichment pathways (see) Figure 9 This is to display indicators such as the significance p-value of pathway enrichment. If the analysis involves multiple comparison groups, the 10 genes with the most significant differences (mean of each omics significance p-value after -log10 transformation) are selected from each comparison group, and a bubble chart of inter-omics correlation coefficients is plotted for all comparison groups (see reference). Figure 15 This was used to visualize the similarities and differences in changes in methylation, transcription, and protein expression levels across comparison groups. All plotting in this step was performed using R language, with the expression profile heatmap created using the pheatmap package and the other plots created using the ggplot2 package.

[0148] Finally, all the analysis results were compiled and arranged into different directories according to categories: the data preprocessing and gene methylation analysis results of each omics were placed in the data preprocessing directory; and the intersection analysis results of differentially expressed genes of each omics were placed in the differential intersection analysis directory.

[0149] Based on the above directory, use a Python script to generate the corresponding web-based analysis report with a single click.

[0150] Example 5

[0151] The following Example 5 further illustrates the application of the present invention in screening key regulatory pathways of methylation and transcriptomics and / or proteomics, but these descriptions are not intended to limit the present invention.

[0152] Accept user's DNA methylation chip data, transcriptome sequencing data, proteome detection data, and GO and KEGG database information; at the same time, fill in the config file, including: the original report storage location, sample analysis name and corresponding group name, group information used for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information.

[0153] First, the data from each omics dataset were preprocessed using Python. For methylation microarray data, the methylation sites were first divided according to the UCSC RefGene annotation information, based on functional region definitions such as promoter regions and gene regions (user-defined functional regions). Then, based on the region division results, the site-gene correspondences were mapped. For transcriptome and proteome data, gene (or protein) quantification results were directly extracted from their standard reports.

[0154] Based on the region division and site-gene mapping results in the DNA methylation pretreatment, non-CpG island methylation sites were first filtered out. Then, according to the site-gene mapping relationship, the mean methylation value of each gene in each functional region was calculated as the gene methylation level. Next, based on the comparison group information provided by the config, the significance P-value of gene methylation t-test for each comparison group was tested. The difference in mean gene methylation value between the two groups (deltaBeta), fold change (FC), and log2 logarithmic transformation of the fold change (log2FC) were also provided. The significance Q value was obtained by correcting the significance P-value using the Benjamini & Hochberg multiple test. Furthermore, based on the difference test threshold provided in the config (default P-value < 0.05 and |deltaBeta| > 0.1), significantly differentially methylated genes were screened, and pathway enrichment analysis was performed. A table of enrichment significance results and a GO enrichment significant pathway bar chart drawn using the R language ggplot2 plotting package are provided (see reference). Figure 8 KEGG enrichment significant pathway bubble chart (reference) Figure 9 ( ), to demonstrate indicators such as the significance of pathway enrichment.

[0155] Next, based on the DNA methylation and differential gene pathway enrichment results from the transcriptome and proteome, a joint analysis of the enriched pathways was performed. First, an intersection analysis was conducted on the significant pathways, and an intersection diagram was drawn (see reference). Figure 16 Further, a scatter plot of pathway significance was drawn (see reference). Figure 17 ), bar chart (reference) Figure 18 This was used to demonstrate the differences in p-values ​​of enrichment analyses for the same pathway across different omics systems. For intersecting pathways, a pathway bubble chart was created (see reference). Figure 19 Stacked diagram (reference) Figure 20 This section provides a detailed overview of the enrichment significance (P-value), enrichment score, and gene ratio of the intersecting pathways. Furthermore, for each intersecting pathway, significantly differentially expressed genes from each omics group are extracted, and a tri-omics expression profile heatmap is generated for each pathway (see reference). Figure 21 This study aimed to reveal the coordinated changes in the expression of key genes within significantly enriched pathways across different omics. All plotting in this step was performed using R language, with the expression profile heatmap created using the ComplexHeatmap package and the remaining plots created using the ggplot2 package.

[0156] Finally, all the analysis results were compiled and arranged into different directories according to categories: the data preprocessing and gene methylation analysis results of each omics were placed in the data preprocessing directory; and the joint analysis results of enriched pathways of each omics were placed in the enriched pathway analysis directory.

[0157] Based on the above directory, use a Python script to generate the corresponding web-based analysis report with a single click.

[0158] The scope of protection of this invention is not limited to the above embodiments. Any variations and advantages that can be conceived by those skilled in the art without departing from the spirit and scope of this invention are included in this invention and are protected by the appended claims.

Claims

1. A method for combined analysis of gene-level DNA methylation, transcriptomics, and proteomics based on the snakemake framework, characterized in that, The analytical method includes the following steps: Step 1: Prepare raw analysis reports and config files for DNA methylation, transcriptomics, and proteomics; Step 2: Preprocess DNA methylation, transcriptome, and proteome data; In step two, the preprocessing of the DNA methylation data includes extracting the site Beta methylation level, the positional relationship between the gene and nearby CpG islands from the raw data from the methylation chip, and performing site functional region division and site-gene mapping. The Beta methylation level refers to the quantitative methylation value of each site in the original DNA methylation microarray report; the positional relationship between the gene and nearby CpG islands includes Island region, Shore region, Shelf region, and OpenSea region; the functional regions include promoter regions, gene regions covering promoter regions and gene body regions, and other functional regions defined based on UCSC RefGene annotation information; the site-gene mapping includes the correspondence between DNA methylation sites contained in different functional regions and their corresponding genes. The preprocessing of the transcriptome data includes extracting gene quantification, differential analysis, and enrichment analysis results from the transcriptome data. The preprocessing of the proteomic data includes extracting protein quantification, differential analysis, and enrichment analysis results from the proteomic data. Step 3: Perform gene methylation analysis on the data after DNA methylation pretreatment; In step three, the gene methylation analysis includes: Based on site annotation information, methylation sites located in CpG islands are identified, and methylation sites that are not CpG islands in different functional regions are filtered out. For each functional region, methylation sites are mapped to genes based on the site-gene mapping relationship; for each gene, the average methylation level of all methylation sites in the current functional region is calculated. For each gene, the t-test was used to calculate the significance of the difference in methylation levels between different experimental groups and the control group in each functional region. The average methylation level of each experimental group and the control group in each functional region was calculated. The difference in the average methylation level between each experimental group and the control group in each functional region was calculated as the deltaBeta value. The ratio of the average methylation level between each experimental group and the control group in each functional region was calculated as the fold change FC value. The log2 transformed value of the fold change FC value between each experimental group and the control group in each functional region was calculated as the log2FC value. Genes with differential methylation levels were screened based on a preset threshold for methylation level differences, and GO and KEGG pathway enrichment analyses were performed on the screened genes. The threshold for methylation level differences was defined as P-value < 0.05 and |deltaBeta| > 0.

1. Step 4: Based on gene methylation analysis data, preprocessed transcriptome data, and preprocessed proteome data, perform correlation analysis, differential intersection analysis, and enrichment pathway joint analysis. In step four, the correlation analysis involves plotting an overall gene correlation scatter plot based on the DNA methylation level, transcription, and protein expression level of the genes; calculating the correlation coefficient of each gene expression; plotting the distribution density of the coefficient in the range of 0 to 1; screening for significantly correlated genes and plotting the correlation scatter plot of single genes; and performing pathway enrichment analysis on significantly correlated genes. In the differential intersection analysis, for the DNA methylation and transcriptional / protein differential analysis results of genes, a violin plot of the degree of difference in transcriptional and protein expression is drawn after grouping according to the significance of DNA methylation differences; a quadrant plot of the degree of difference between DNA methylation and transcriptome and proteome is drawn; an intersection analysis is performed on genes with significant differences in each omics, and pathway enrichment analysis is performed on the intersection genes; a multi-omics expression profile of the intersection genes is drawn; if multiple comparison groups of gene methylation and transcriptome and proteome are involved, a bubble plot of the correlation of intersection genes in each comparison group is drawn. In the enrichment pathway joint analysis, for gene methylation and differentially enriched genes and proteins in the transcriptome and proteome, pathway intersection analysis is performed, and pathway significance scatter plots and bar charts are drawn; for the intersection pathways, pathway bubble charts, pathway stacking charts, and expression profiles of significant genes and proteins in each omics contained in the pathway are drawn. Step 5: Based on the analysis results in Step 4, generate a visualization report.

2. The analytical method as described in claim 1, characterized in that, In step one, the config file includes: the storage location of the original report and the corresponding sample analysis name and group name, grouping information for differential analysis, fold change and significance parameters, and GO and KEGG database annotation information; the original analysis report and the config file are used for subsequent data processing and analysis.

3. The analytical method as described in claim 1, characterized in that, In the correlation analysis, the correlation coefficient is calculated using methods including the Spearman method and the Pearson method; significant correlation refers to a correlation coefficient with a significance P-value < 0.05; the pathway enrichment analysis of significantly correlated genes includes GO and KEGG enrichment analysis. In the differential intersection analysis, the DNA methylation difference significance grouping includes hypomethylation, no significant change, and hypermethylation; hypomethylation refers to the group with P value < 0.05 and deltaBeta < -0.1, no significant change refers to the group with P ≥ 0.05 or -0.1 ≤ deltaBeta ≤ 0.1, and hypermethylation refers to the group with P value < 0.05 and deltaBeta > 0.1; the pathway enrichment analysis of significantly related genes includes GO and KEGG enrichment analysis; In the combined enrichment pathway analysis, the enriched pathway refers to the pathways and their enrichment results obtained from the KEGG enrichment analysis of the individual analysis results of gene methylation, transcriptomics, and proteomics; the enrichment results include pathway ID, pathway name, pathway description information, enrichment significance P-value, enrichment score, gene ratio or protein ratio of the enriched pathway, and specific gene or protein differentials contained in the pathway.

4. The analytical method as described in claim 1, characterized in that, In step five, when generating the visualization report, a web-based report on the combined analysis of DNA methylation, transcriptomics, and proteomics is generated with a single click based on the analysis results. The entire analysis results are summarized, and each analysis step is described, accompanied by corresponding charts and pop-up help documents. The web-based report includes internal shortcuts, descriptions of analysis methods, and links to external websites, enabling quick navigation within the web-based version and rapid access to online resources.

5. A multi-omics joint analysis system for implementing the joint analysis method as described in any one of claims 1-4, characterized in that, The multi-omics joint analysis system includes: a file data input module, a data preprocessing module, a methylation analysis module, a multi-omics joint analysis module, and a visualization generation and report output module; The file data input module is used to input prepared raw analysis reports of DNA methylation, transcriptomics, and proteomics, as well as config configuration files; The data preprocessing module is used to preprocess DNA methylation, transcriptome, and proteome data; The methylation analysis module is used to analyze the DNA methylation level of genes and perform differential analysis and pathway enrichment analysis. The multi-omics joint analysis module is used to perform correlation analysis, differential intersection analysis, and enrichment pathway joint analysis on DNA methylation, transcriptomics, and proteomics data. The visualization generation and report output module is used to visualize the analysis results and generate the final report.

6. The combined analysis method as described in any one of claims 1-4, or the multi-omics combined analysis system as described in claim 5, in the combined analysis of DNA methylation with transcriptome and proteome at the gene level, the combined analysis of DNA methylation with transcriptome and proteome at the site level, the screening of genes associated with methylation expression profiles with transcriptome and proteome, the screening of genes with differential intersections between methylation and transcriptome and proteome, and the screening of key regulatory pathways between methylation and transcriptome and proteome.

Citation Information

Patent Citations

  • Transcriptome and DNA methylation data correlation analysis method and system

    CN112201302A

  • Analysis of genomic word frameworks on genomic methylation data

    WO2023183907A2