Method and system for detecting genes associated with complex traits based on summary statistical data

Abstract

The present invention discloses a method and system for detecting genes associated with complex traits based on summary statistical data. The method can parse the summary statistical data and gene expression weight data of complex traits according to custom settings, and extract effective genetic information therefrom, including a set of genetic variation sites within the gene range and their corresponding gene expression weights and significance P values. Based on these weights and P values, a joint Cauchy association test is applied to identify genes that are significantly associated with complex traits. The present invention innovatively proposes a method for detecting genes associated with traits that is suitable for summary statistical data, which significantly expands the analysis scope of the data and thus improves the statistical effectiveness of the association. At the same time, based on the analytical algorithm, the computational efficiency of the association test is effectively improved, the computational burden is reduced, and the rapid and accurate identification of genes associated with complex traits is achieved.

Description

Technical Field

[0001] The present invention relates to a method and system for detecting genes associated with complex traits based on summary statistical data, and belongs to the application field of bioinformatics and statistical genetics technology. Background Art

[0002] With the rapid development of sequencing technology, research on rare variants has received widespread attention, especially in the association analysis of complex traits. Rare variants can explain some previously unknown heritability, thereby more accurately finding genes associated with traits and helping people to gain a deeper understanding of their underlying biological mechanisms. However, due to the low probability of rare variant sites in the population, the univariate analysis method commonly used in genome-wide association studies (GWAS) is not effective in practical applications when the research sample size is not large enough. Therefore, multivariate-based association analysis methods have emerged. These methods aggregate the effectiveness of multiple rare variant sites and use a set of variant sites for association analysis with complex traits, which can effectively make up for the shortcomings of univariate analysis.

[0003] Currently, multivariate association analysis methods, such as the sequence kernel association test (SKAT) and burden test, have been widely used. However, further studies have shown that in sparse cases where only a small number of variants have an effect on a trait, multivariate association analysis methods lose some power. In such cases, exponential combination tests are needed to integrate multiple studies to assess overall significance. However, due to its high computational requirements, this method imposes a significant computational burden on large-scale sequences, significantly hindering the study of the mechanisms underlying the development of complex traits. Furthermore, another key issue is that previous association methods have focused on using individual genotype data, which inherently limits sample size due to privacy concerns. However, summary statistics, which do not involve personal privacy, significantly increase the sample size available for analysis and can significantly improve the reliability of genetic analyses. Consequently, summary statistics have attracted widespread attention from researchers. Currently, a vast amount of summary statistical datasets awaits exploration. Therefore, a method and system for detecting genes associated with complex traits based on summary statistics is urgently needed. Summary of the Invention

[0004] The primary technical problem to be solved by the present invention is to be able to use summary statistical data to perform gene trait association analysis.

[0005] The present invention also aims to solve another technical problem, namely, to achieve fast and efficient processing of gene trait association analysis and reduce the computational burden.

[0006] The purpose of the present invention is achieved through the following technical solutions.

[0007] The present invention first discloses a method for detecting genes associated with complex traits based on summary statistical data, which realizes efficient and reliable identification of genes associated with complex traits based on summary statistical data, comprising the following steps:

[0008] Get the specified association analysis parameters. Set parameters to specify the association analysis mode. Parameters that can be set include: the location of the summary statistics file, the chromosome of the variant site in the summary statistics, the column number of the physical location and significant P value, the location of the gene expression weight file, the gene name of the association analysis, the upstream and downstream extension range of the gene region, the R language environment directory, and the result output file directory.

[0009] The input data is obtained in text format with space delimiters. The present invention extracts information from user-provided summary statistics and gene expression weight data. The summary statistics include the chromosome number of each genetic variation site, the physical location of the site on the chromosome, and the significance P value of the site's association with the target trait. The gene expression weight data includes the gene name, the chromosome number of the gene, the starting and ending positions of the gene region, and the genetic variation site information contained in the gene region. The variation site information specifically includes the variation site name, physical location, allele type, and the effect size of each variation site on the gene expression.

[0010] In addition to the proactively provided weight data, the present invention also provides standard gene expression weight data. This weight file uses an elastic net model to train the genotype data of the GTEx dataset and the expression data on whole blood tissue to obtain the expression weights of genetic variant sites in whole blood tissue. This is then combined with the corresponding genetic information to generate the required standard gene expression weight file.

[0011] Parse summary statistics and gene expression weight data. Read the aforementioned input data and, for each gene, segment and filter the summary statistics based on the gene region range in the gene expression weight data. Specifically, obtain the target gene's chromosome, the start and end positions of the gene region, and the set of weighted variant sites within the gene range from the expression weight data. Users can customize the upstream and downstream extension ranges of the gene region by specifying a window size: downstream position = lower bound of gene range - window size, upstream position = upper bound of gene range + window size.

[0012] For each variant site in the summary statistics, we need to determine: whether the chromosome of the variant site is consistent with the gene chromosome; whether the physical location of the variant site is between the start and end positions of the gene region; and whether the variant site is in the weighted site set of the target gene. Then, we obtain all eligible genetic variant sites within the gene region, that is, variant sites located in the gene region and with effect size data.

[0013] After obtaining the set of genetic variation sites corresponding to each gene and its associated expression weight and significance P value, a gene-site data subset for subsequent analysis is generated.

[0014] Gene-trait association significance analysis. Based on the above process, the expression weights and P values ​​of the variable site set were obtained to conduct the joint Cauchy association test (ACAT): For each gene, the formula tan{(1-p i )π}The original significance P value p i The distribution of the subject is transformed into the Cauchy distribution, and then the significance P value p′ after the transformation is i The joint Cauchy correlation test statistic is obtained by weighting with expression weights:

[0015]

[0016] where w i represents the expression weight of the i-th variant site, p i represents the significant P value of the variant site, and k represents the number of variant sites contained in the gene region. Based on the joint statistic T ACAT and scaling parameters A two-sided test was used to calculate the joint P value for each gene:

[0017]

[0018] Generate gene-trait association results. Since only P value and weight value are required in this calculation process, only summary statistics of traits and gene expression weight data are needed, and no genotype data is needed. In addition, since the statistic T ACAT The distribution of is approximately Cauchy, and is therefore insensitive to interactions between variables, does not require consideration of covariates between variables, and can be used to determine correlation significance through analytical methods. Regarding the above analysis process, the present invention also possesses exception handling capabilities, providing a certain degree of fault tolerance for missing parameters and abnormal instructions.

[0019] After outputting the gene names and corresponding P values ​​tested for association with the trait, the results can be subjected to a Bonferroni test. For each gene, if the combined Cauchy association test P value is less than or equal to 0.05 divided by the number of genes tested, the gene is considered to be significantly associated with the trait being studied. After this screening process, genes associated with the trait of interest are identified.

[0020] The present invention also discloses a complex trait-associated gene detection system based on summary statistical data, the system comprising at least one processor and a memory connected to the processor, wherein: the memory is used to store computer programs or instructions; the processor is used to obtain and execute the computer programs or instructions in the memory so that the system performs the detection method disclosed in the present invention.

[0021] Beneficial effects of the present invention

[0022] The present invention provides a method and system for detecting genes associated with complex traits based on summary statistical data. By using a joint Cauchy association test, it innovatively realizes gene-trait association analysis of trait summary statistical data, significantly expands the scope of application of the data, and indirectly increases the sample size that can be used for analysis, thereby effectively improving the statistical efficiency of the association test. Compared with traditional methods, the present invention can quickly and efficiently parse summary data and directly calculate the significance of gene associations through analytical methods, which not only greatly reduces the computational burden, but also improves the overall efficiency and accuracy of association analysis. Therefore, the present invention provides strong support for trait association research, helps researchers to efficiently carry out association analysis, promotes an in-depth understanding of the biological mechanisms behind complex traits, and also provides new ideas and methods for scientific research and application in genetics and related fields. BRIEF DESCRIPTION OF THE DRAWINGS

[0023] Figure 1 Schematic diagram of the process of identifying genes associated with complex traits based on summary statistical data provided by the present invention

[0024] Figure 2 The judgment process required for traversing the variant sites in the analysis of gene expression weight data in the present invention

[0025] Figure 3 Generate a simulation process for testing correlation efficiency and false positive rate in the embodiment of the present invention

[0026] Figure 4 A line graph showing simulation results under different assumptions in an embodiment of the present invention DETAILED DESCRIPTION

[0027] The embodiments of the present application provide a method and system for detecting genes associated with complex traits based on summary statistical data. The technical content of the present invention is described in detail below with reference to the accompanying drawings and specific examples.

[0028] In view of the problems of insufficient statistical power and huge consumption of computing resources in the current methods for identifying genes associated with complex traits, this paper combines Figure 1 The specific implementation of the method for identifying genes associated with complex traits based on summary statistical data proposed in the present invention is as follows:

[0029] According to the analysis requirements of the specific scenario, set parameters to specify the mode for performing association analysis. For the type 1 diabetes (T1D) data from the World Health Organization (WTCCC) database used in this embodiment, the settings are as follows: the summary statistical data parameter summary is "T1D.summary", the chromosome of the variant site in the summary statistical data, the physical position and the column number parameter info_columns where the significant P value is located is "1,3,9", the gene expression weight file parameter weights is "weight.txt", the gene name parameter gene of the association analysis is "H3C6", the upstream and downstream extension range windows of the gene region is "5000000", the R language environment address parameter Rscript is "Rscript", and the result output address parameter output_folder is "ACAT_result". The gene parameter can be set in turn to all genes in the file pointed to by the weights parameter to find all genes significantly associated with type 1 diabetes.

[0030] Next, extract information from the summary statistics file T1D.summary and the gene expression weight data file weight.txt. For summary statistics and gene expression weight data for complex traits, the data format must be text, with a space delimiter.

[0031] The summary statistics must include the chromosome number of each genetic variation site, the physical location of the site on the chromosome, and the significant P value associated with the target trait. In this embodiment, the required information is extracted from the 1st, 3rd, and 9th columns of the T1D.summary file by specifying the parameter info_columns to "1,3,9". The gene expression weight data weight.txt includes the gene name, the chromosome number where the gene is located, the starting and ending positions of the gene region, and the genetic variation site information contained in the gene region; and the variation site information specifically includes the variation site name, physical location, allele type, and the effect size of each variation site on the gene expression. Different information is separated by underscores, for example, "rs12565286_785910_G_C_0.0".

[0032] The weight.txt used in this embodiment is a standardized gene expression weight data file. The weight file is based on the elastic network model, and the genotype and expression data of the whole blood tissue in the GTEx data set are trained to generate the expression weights of the genetic variation sites in the whole blood tissue, and combined with the corresponding genetic information to generate a standardized gene expression weight file. Specifically, the genotype data of GTEx is screened and segmented according to the gene information appearing in the expression data, and the sample genotype under each gene range is obtained as the independent variable X. At the same time, the tissue expression value of the gene in each sample is used as the dependent variable Y, and the environmental factor is ε to build a model:

[0033] Y=WX+ε,

[0034] The weight vector W for each variant site is calculated using an elastic network. This is then combined with other information from the variant site file to generate the required standard gene expression weights file. In addition to the provided weight.txt file, this method supports user-defined weights files. Users simply provide a gene expression weights file in the standard format described above to adapt to different analysis needs.

[0035] Parse summary statistics and gene expression weight data. Read the input data above and, for each gene, segment and filter the summary statistics based on the gene region range in the gene expression weight data. Specifically, obtain the chromosome of the target gene, the start and end positions of the gene region, and the set of weighted variant sites within the gene range from the expression weight data. The specified window size defines the upstream and downstream extension range of the gene region: downstream position = lower bound of gene range - window size, upstream position = upper bound of gene range + window size. The window size setting determines the number of genetic variant sites within the effective gene range. The larger the range, the more genetic variant sites are assumed to affect gene expression.

[0036] Traverse each variant site in the summary statistics, such as Figure 2 As shown, the following conditions must be met: the chromosome of the variant site is consistent with the gene chromosome; the physical location of the variant site is between upstream and downstream locations in the gene region; and the variant site is included in the weighted locus set of the target gene. Next, all eligible genetic variant sites within the gene region, i.e., variant sites located within the gene region with effect size data, are obtained. After obtaining the set of genetic variant sites corresponding to each gene and its associated expression weight and significance P value, a gene-site data subset is generated for subsequent analysis.

[0037] Perform gene-trait association significance analysis. Based on the above process, the expression weight and P value of the genetic variation set corresponding to the gene are obtained to perform a joint Cauchy association test: For each gene, the formula tan{(1-p i )π} the original significance P value p of the variant sites in the set i Transformed into a Cauchy distribution, and then the significance P value p′ after transformation i The weighted sum of the expression weights corresponding to each variant site is used to obtain the joint Cauchy association test statistic:

[0038]

[0039] w i represents the expression weight of the i-th genetic variation site;

[0040] p i Represents the significant P value of the variant site;

[0041] k represents the number of variant sites contained in the gene region.

[0042] Based on the joint statistic T ACAT and scaling parameters A two-sided test was used to calculate the joint P value for each gene:

[0043]

[0044] Generate gene-trait association results. Since only P values ​​and weight values ​​are required in this calculation process, only summary statistics and gene expression weight data of type 1 diabetes are needed, and no genotype data are required. In addition, since the statistic T ACA The distribution of is approximately Cauchy distribution, so it is insensitive to the interaction between variables, does not need to consider the covariates between variables, and can obtain the significance of association through analytical methods.

[0045] Finally, after outputting the gene names and corresponding test P values ​​that have passed the association test, the results can also be subjected to a Bonferroni test. For each gene, if the combined Cauchy association test P value is less than or equal to 0.05 divided by the number of genes involved in the test, the gene is considered to be significantly associated with type 1 diabetes. After the Bonferroni test, 69 genes were finally obtained that were significantly associated with type 1 diabetes. The names of the significantly associated genes and their corresponding P values ​​are shown in Table 1:

[0046] Table 1: Genes significantly associated with type 1 diabetes and their corresponding P values ​​found in this example

[0047]

[0048]

[0049] In addition to performing gene association tests based on real type 1 diabetes data, in this embodiment, analysis can also be performed based on simulated data to test the association effectiveness of the method. Specifically, the simulation data is generated by the following steps: for each gene in the selected gene set, genetic variation sites within the gene range are randomly selected at a ratio of (1%, 2%, 5%, 10%), and the selected sites are randomly assigned effect sizes that obey the standard normal distribution and are regarded as effect sites. Secondly, based on the additive model and the genetic assumption of causal pleiotropy, simulated gene expression and phenotype are generated for each gene through the simulated effect sites, where the heritability from genetic variation sites to gene expression is Set to (0.05, 0.10, 0.15, 0.20, 0.25), the heritability from gene expression to phenotype Set to (0.05, 0.10, 0.15, 0.20, 0.25).

[0050] After generating 100 (4×5×5) sets of simulated data for training and test sets in the above manner, similar to the standard weight file, the elastic network model is trained on the genotype data and expression data of the training set to obtain the simulated gene expression weights. In the test set, the summary statistics and gene expression weights generated by the simulation are input into the system, and significant genes are identified as described above. The entire simulation process is as follows Figure 3 As shown, where X is the genotype; β is the simulation effect, β pred is the elastic network prediction effect; E is the simulation expression, where E ε Expression affected by environmental noise, E G The expression of E is affected by genetic factors. pre is the elastic network predicted expression; P is the simulated phenotype, where P ε Phenotypes affected by environmental noise, P G P is the phenotype affected by genetic factors;value is the P value of the association between genes and traits. Finally, the number of significantly associated genes in a total of 100 cases is shown in Table 2 below. The statistical effect results are as follows Figure 4 As shown in the figure, the results show that the statistical power of the method increases with the increase of heritability, and decreases with the increase of the number of effective variant sites, which is consistent with the hypothesis. This shows that the method can effectively and accurately identify associated gene sites under simulation conditions and has good practicality and reliability.

[0051] Table 2 Number of significant genes under multiple hypothesis simulation data

[0052]

[0053] To test the false positive rate of the present invention's test results, gene expression and phenotype data were randomly generated according to a standard normal distribution. Training and trait association based on the randomly generated gene expression and phenotype data ultimately yielded an association test false positive rate of 0.048.

[0054] The above embodiments are only some specific implementations of the present invention. Although the description is relatively detailed, it does not limit the scope of the present invention. It should be pointed out that those skilled in the art can still make various modifications and improvements without departing from the concept of the present invention. These modifications and improvements all fall within the scope of protection of the present invention. Therefore, the scope of protection of the patent of the present invention shall be based on the appended claims.

Claims

1. A method for detecting genes associated with complex traits based on summary statistical data, characterized in that The steps include: Step 1: Get the specified association analysis parameters and set the parameters to specify the mode for association analysis; Step 2: Read the necessary input data for association analysis, read the summary statistics of complex traits and gene expression weight data, and extract information from them; Step 3: Parsing summary statistics and gene expression weight data, specifically including: reading the input data described in step 2 and performing the following steps: Step 3.1: For each gene, segment and filter the summary statistics based on the gene region range in the gene expression weight data. Obtain the target gene's chromosome, the start and end positions of the gene region, and the set of weighted variant sites within the gene range from the expression weight data. Next, traverse each variant site in the summary statistics to determine whether it meets the criteria. Finally, obtain all eligible genetic variant sites within the gene region, that is, variant sites located in the gene region and with effect size data. Step 3.2: Obtain the set of genetic variation sites corresponding to each gene and its associated expression weight and significance P value to generate a gene-site data subset for subsequent analysis; Step 4: Gene association significance analysis, based on the expression weights and P values ​​obtained in step 3, conduct a joint Cauchy association test: For each gene, the formula tan{(1-p i )π}The original significance P value p i The distribution of the subject is transformed into the Cauchy distribution, and then the significance P value p after the transformation is i ′ The expression weights are used for weighting to obtain the joint Cauchy association test statistic: where w i represents the expression weight of the i-th variant site, p i represents the significant P value of the variant site, and k represents the number of variant sites contained in the gene region; Based on the joint statistic T, that is, T ACAT , scaling parameters A two-sided test was used to calculate the joint P value for each gene: And generate gene-trait association results.

2. The method according to claim 1, wherein The parameters that need to be set include: the location of the summary statistical data file, the chromosome of the variant site in the summary statistical data, the column number of the physical position and significant P value, the location of the gene expression weight file, the gene name of the association analysis, the upstream and downstream extension range of the gene region, the R language environment address, and the result output file address.

3. The method according to claim 1, wherein Summary statistics and gene expression weight data require the following genetic information: Summary statistics include the chromosome number of each genetic variation site, the physical position of the site on the chromosome, and the significant P value of the association between the site and the target trait; Gene expression weight data includes the gene name, the chromosome number where the gene is located, the starting and ending positions of the gene region, and the genetic variation site information contained in the gene region; the variation site information specifically includes the variation site name, physical location, allele type and the effect size of each variation site on the gene expression.

4. The method according to claim 1, wherein The conditions specifically include: the chromosome of the variant site is consistent with the gene chromosome; the physical position of the variant site is between the starting position and the ending position of the gene region; and the variant site is in the weighted site set of the target gene.

5. The method according to claim 1, wherein The gene expression weight data uses elastic network training of GTEx dataset genotype data and whole blood expression data to obtain the expression weight information of each genetic variation site for the gene; the weight model is constructed by screening and segmenting the GTEx genotype data according to the gene information appearing in the expression data, obtaining the sample genotype under each gene range as the independent variable X, and taking the tissue expression value of the gene in each sample as the dependent variable Y, and the environmental factor as ε, to construct the model: Y=WX+ε The weight vector W of the variant site is obtained based on elastic network calculation, and then combined with other information of the variant site in the file to generate the required standard gene expression weight file.

6. The method according to claim 1, wherein By actively capturing various abnormal situations that may occur during the analysis process within the method, it has a certain fault tolerance capability.

7. The method according to claim 1, wherein The process of parsing the summary statistics and gene expression weight data described in step 3 can adjust the set of variant sites used for the joint Cauchy association test according to the upstream and downstream extension ranges of the user-defined gene region.

8. The method according to claim 1, wherein Step 4 outputs the name of the gene that has been tested for association with the trait and the test P value corresponding to the gene.

9. The method according to claim 8, wherein The P value was subjected to a Bonferroni test. For each gene, if the P value of the joint Cauchy association test was ≤ 0.05 divided by the number of genes involved in the test, the gene was considered to be significantly associated with the studied trait.

10. A complex trait-associated gene detection system based on summary statistical data, characterized in that The system comprises at least one processor and a memory connected to the processor, wherein: The memory is used to store computer programs or instructions; The processor is configured to obtain and execute a computer program or instruction in a memory, so as to enable the system to perform the method according to any one of claims 1 to 9.

Images