Computer device and method for mining stress response-related variant sites or molecular QTL on basis of cross-condition comparison analysis

Abstract

Provided are a computer device and method for mining stress response-related variant sites or molecular QTL on the basis of cross-condition comparison analysis. Firstly, molecular traits are transformed by means of the Bray-Curtis distance, so that the transformed molecular traits can represent the influence of natural variation on molecular traits, including the genetic effects and stress response effects of natural variation. Then, comparative PCA is used to compare and reduce the dimensionality of molecular traits under a stress treatment and molecular traits under control conditions, thereby accurately obtaining the stress response effects of the molecular traits and using them as the stress response index of the molecular traits. Finally, the molecular QTL of the stress response is obtained by means of a genome-wide association analysis of the stress response index, which is used for identifying cross-condition molecular response QTL, such as those for stress stimulation.

Description

Computer devices and methods for mining stress-related variant sites or molecular QTLs based on cross-condition comparative analysis Technical Field

[0001] This invention belongs to the field of bioinformatics technology, specifically relating to computer devices and methods for mining stress response-related variant sites or molecular QTLs based on cross-condition comparative analysis. Background Technology

[0002] Molecular traits refer to characteristics of an organism determined by genetic material at the molecular level. These traits typically involve the molecular structure and function of the organism, including karyotype, genetic material configuration, genotype, base sequence, messenger RNA, and specific surface antigens.

[0003] Under stimulus conditions, the influence of natural variation on molecular traits involves two aspects. First, different alleles lead to differences in genetic effects among individuals; second, different alleles lead to differences in stimulus response effects among individuals. However, when identifying molecular traits associated with stimulus responses, assessment is usually based solely on the mean differences of molecular traits between populations under control and treatment conditions. This approach ignores the differences in stress responses among different alleles and cannot accurately identify population-level molecular traits associated with stimulus responses.

[0004] Genetic variations identified through genome-wide association studies (GWAS) may simultaneously influence both molecular traits and complex traits. This suggests that some molecular QTLs affect complex traits by influencing molecular traits (Gusev et al., 2018; Liu et al., 2020a; Liu et al., 2020b; Watanabe et al., 2017). Utilizing the functional information conferred by the omics characteristics of molecular QTLs, association analyses of multi-omics QTL data will help elucidate the regulatory mechanisms of complex traits (Gamazon et al., 2015; Wu et al., 2023; Wu et al., 2018b; Zhu et al., 2016). Due to the food yield crisis caused by drastic global climate change, researchers are increasingly focusing on improving crop stress tolerance (Gao and Cui, 2023; Yang et al., 2023). Subsequently, multi-omics association analysis has also been used to analyze the regulatory mechanisms of crop stress responses (Jiang et al., 2022; Liang et al., 2021; Liu et al., 2020b; Sun et al., 2023; Wu et al., 2021). However, existing multi-omics association analysis studies mainly focus on omics data under single conditions, and it remains challenging to analyze stress multi-omics data that require comparison between controls and treatments.

[0005] A key step in stress multi-omics association analysis is the identification of stress-responding molecules. This identification process requires comparing molecular QTL loci under control and stress treatment conditions to identify naturally occurring variants that show differences in response between stress treatment and control. Existing studies have included methods such as comparing the effect sizes of QTLs under two conditions to identify stress-responding molecular QTLs (Wu et al., 2018a), and statistical tests to screen molecular features before identifying QTLs, using features that differ between the two conditions as stress-responding features. These features are then used to identify stress-responding molecular QTLs through GWAS (Jianget et al., 2022; Wu et al., 2021; Zhang et al., 2021). These two approaches consider the genetic effects of natural variation on molecular traits and the stress responses of molecular traits, respectively. However, when analyzing stress-related population omics data, the genetic effects of natural variation on molecular traits and the stress response effects of natural variation on molecular traits are intertwined (Lea et al., 2019; Napier et al., 2023). This heterogeneity prevents either approach from accurately identifying stress-responsive molecular QTLs. Invention Overview

[0006] The technical problem to be solved by this invention is how to identify stress response-related molecular QTLs and / or how to identify stress response-related genes and / or how to identify natural variant sites that lead to stress response.

[0007] To address the aforementioned technical problems, the present invention first provides a computer device, including a memory, a processor, and a computer program stored in the memory, characterized in that the processor executes the computer program to perform the following steps:

[0008] S1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; Receive SNP genotyping data of the genome of each individual in the genetic population;

[0009] The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition.

[0010] S2) Data Processing: The expression level of each gene in the genetic population under the stimulation condition is transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression level of each gene in the genetic population under the control condition is transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA.

[0011] Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association (GWAS) analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators.

[0012] S3) Data output: Output the significant association SNPs.

[0013] In the aforementioned computer device, the control condition may be a non-stimulatory condition, specifically the normal growth conditions of the genetic population.

[0014] S1) further includes the step of obtaining SNP genotyping data for the genome of each individual in the genetic population.

[0015] The significance analysis in S2) may include the following steps:

[0016] The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition are calculated to obtain the first regression coefficient (also known as regression coefficient a). The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition are calculated to obtain the second regression coefficient (also known as regression coefficient b). The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

[0017] In the aforementioned computer device, the genetic population may be a plant genetic population, and the stimulus may be a stress treatment.

[0018] The stress treatment may be salt stress treatment, or other abiotic stress treatments, such as drought stress, flooding stress, low temperature stress, high temperature stress, hormone stress, mechanical stress, etc.

[0019] The PCA algorithm described above can be a comparative PCA algorithm based on sparse PCA.

[0020] To address the aforementioned technical problems, the present invention also provides an apparatus for identifying stimulus-response related QTLs, the apparatus comprising the following modules:

[0021] A1) Data receiving module: used to receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; and to receive SNP genotyping data of the genome of each individual in the genetic population;

[0022] The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition.

[0023] A2) Data Processing Module: This module is used to convert the expression level of each gene in the genetic population under the stimulation condition using Bray-Curtis distance to obtain a Bray-Curtis distance matrix under the stimulation condition; convert the expression level of each gene in the genetic population under the control condition using Bray-Curtis distance to obtain a Bray-Curtis distance matrix under the control condition; perform dimensionality reduction on the difference between the Bray-Curtis distance matrix of each gene under the stimulation condition and the Bray-Curtis distance matrix under the control condition using the PCA algorithm to obtain the principal components of candidate contrastive PCA; and perform significance analysis on the principal components of candidate contrastive PCA to obtain the principal components of contrastive PCA.

[0024] Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association (GWAS) analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators.

[0025] A3) Haplotype Analysis Module: This module uses the expression level of each gene in the genetic population as a molecular trait. Based on the transcriptome data under the control conditions and the SNP genotyping data of each individual's genome in the genetic population, it performs genome-wide association (GWAS) analysis on the molecular trait to obtain the significantly associated SNPs of the molecular trait. Then, it performs cluster analysis on the significantly associated SNPs of the molecular trait to obtain the haplotypes associated with the molecular trait.

[0026] A4) QTL Analysis Module: Used to obtain stimulus-response-related candidate QTLs based on the significant association SNPs of the stimulus response index, and to perform a two-way ANOVA on the stimulus-response-related candidate QTLs using the molecular trait-related haplotypes and the stimulus conditions as factors, to obtain the stimulus-response-related QTLs.

[0027] In the above-mentioned device, the control condition can be a non-stimulatory condition, and the control condition can specifically be the normal growth conditions of the genetic population.

[0028] The PCA algorithm described above can be a comparative PCA algorithm based on sparse PCA.

[0029] In the above-described apparatus, the significance analysis in A2) may include the following steps:

[0030] The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition are calculated to obtain the first regression coefficient (also known as regression coefficient a). The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition are calculated to obtain the second regression coefficient (also known as regression coefficient b). The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

[0031] In the above-mentioned device, the genetic population can be a plant genetic population, and the stimulus can be a stress treatment. The stress treatment can be salt stress treatment, or it can be other abiotic stress treatments, such as drought stress, flooding stress, low temperature stress, high temperature stress, hormone stress, mechanical stress, etc.

[0032] To address the aforementioned technical problems, the present invention also provides a method for identifying causal (mutated) genes related to stimulus phenotypes, the method comprising the following steps:

[0033] B1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; receive SNP genotyping data of the genome of each individual in the genetic population;

[0034] The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition.

[0035] B2) Data Processing: The expression levels of each gene in the genetic population under the stimulation condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression levels of each gene in the genetic population under the control condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA.

[0036] Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association (GWAS) analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators.

[0037] B3) Haplotype analysis: The expression level of each gene in the genetic population is taken as a molecular trait. Based on the transcriptome data under the control conditions and the SNP genotyping data of each individual genome in the genetic population, genome-wide association (GWAS) analysis is performed on the molecular trait to obtain the significantly associated SNPs of the molecular trait. Cluster analysis is performed on the significantly associated SNPs of the molecular trait to obtain the haplotypes associated with the molecular trait.

[0038] B4) QTL analysis: Based on the significant association of the stimulus response index with SNP, candidate QTLs related to stimulus response are obtained. Using the haplotypes related to the molecular traits and the stimulus conditions as factors, a two-way ANOVA is performed on the candidate QTLs related to stimulus response to obtain the stimulus response-related QTLs.

[0039] B5) Mendelian randomization analysis: Receive stimulus-phenotype related data of the genetic population, and based on the stimulus-response related QTLs and the stimulus-phenotype related data, obtain the contribution of genes containing the stimulus-response related QTLs to the stimulus phenotype through Mendelian randomization analysis, and obtain causal genes related to the stimulus phenotype based on the contribution.

[0040] In the above method, the control condition can be a non-stimulatory condition, and specifically, the control condition can be the normal growth condition of the genetic population.

[0041] The PCA algorithm can be a contrastive PCA algorithm based on sparse PCA.

[0042] In the above method, the significance analysis in B2) may include the following steps:

[0043] The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition are calculated to obtain the first regression coefficient (also known as regression coefficient a). The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition are calculated to obtain the second regression coefficient (also known as regression coefficient b). The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

[0044] The genetic population mentioned in the above method can be a plant genetic population, and the stimulus can be a stress treatment. The stress treatment can be salt stress treatment, or other abiotic stress treatments, such as drought stress, flooding stress, low temperature stress, high temperature stress, hormone stress, mechanical stress, etc.

[0045] To address the aforementioned technical problems, the present invention also provides a method for identifying stimulus-response related QTLs, the method comprising the following steps:

[0046] C1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; receive SNP genotyping data of the genome of each individual in the genetic population;

[0047] The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition.

[0048] C2) Data Processing: The expression levels of each gene in the genetic population under the stimulation condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression levels of each gene in the genetic population under the control condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA.

[0049] Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association (GWAS) analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators.

[0050] C3) Haplotype analysis: The expression level of each gene in the genetic population is taken as a molecular trait. Based on the transcriptome data under the control conditions and the SNP genotyping data of each individual genome in the genetic population, genome-wide association (GWAS) analysis is performed on the molecular trait to obtain the significantly associated SNPs of the molecular trait. Cluster analysis is performed on the significantly associated SNPs of the molecular trait to obtain the haplotypes associated with the molecular trait.

[0051] C4) QTL analysis: Based on the significant association SNPs of the stimulus response index, candidate QTLs related to stimulus response are obtained. Using the haplotypes related to the molecular traits and the stimulus conditions as factors, a two-way ANOVA is performed on the candidate QTLs related to stimulus response to obtain the stimulus response-related QTLs.

[0052] In the above method, the control condition can be a non-stimulatory condition, and specifically, the control condition can be the normal growth condition of the genetic population.

[0053] The PCA algorithm can be a contrastive PCA algorithm based on sparse PCA.

[0054] In the above method, the significance analysis in C2) may include the following steps:

[0055] The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition are calculated to obtain the first regression coefficient (also known as regression coefficient a). The regression coefficients of the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition are calculated to obtain the second regression coefficient (also known as regression coefficient b). The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

[0056] The genetic population mentioned in the above method can be a plant genetic population, and the stimulus can be a stress treatment. The stress treatment can be salt stress treatment, or other abiotic stress treatments, such as drought stress, flooding stress, low temperature stress, high temperature stress, hormone stress, mechanical stress, etc.

[0057] The causal genes mentioned above can be causal mutant genes. A causal mutant gene refers to a mutation that occurs in the genome and has a causal relationship with a specific disease, phenotype, or physiological phenomenon. This mutation may lead to abnormal function or inactivation of certain genes, thereby affecting the production or function of related proteins, and consequently affecting cellular physiological processes.

[0058] To address the aforementioned technical problems, the present invention also provides a computer-readable storage medium storing a computer program that enables a computer to perform steps B1)-B2) of the method described above, or steps B1)-B4) of the method described above, or steps B1-B5) of the method described above.

[0059] The application of the computer device, the device, and / or the computer-readable storage medium described above in plant stress breeding is also within the scope of protection of this invention.

[0060] This invention first transforms molecular traits using Bray-Curtis distance, enabling the transformed molecular traits to represent the influence of natural variation on the traits, including the genetic effects of natural variation and stress response effects. Then, comparative PCA is used to perform dimensionality reduction on molecular traits under stress treatment and those under control, thereby accurately obtaining the stress response effect of the molecular traits and using it as a stress response index. Finally, genome-wide association analysis of the stress response index is used to obtain molecular QTLs of stress response, which are then used to identify cross-conditional molecular response QTLs such as those related to stress stimuli. Attached Figure Description

[0061] Figure 1 is a flowchart illustrating the computer device of the present invention.

[0062] Figure 2 shows the effect of natural variation on gene expression using Bray-Curtis distance.

[0063] Figure 3 shows the dynamic regulatory effect of ZmRH3 by comparative PCA analysis.

[0064] Figure 4 shows a comparison of the regression coefficients of relative gene expression and cPCA principal components under control and salt stress treatment.

[0065] Figure 5 shows the functional analysis of QTLs expressed in response to salt stress. Embodiments of the present invention

[0066] The present invention will now be described in further detail with reference to specific embodiments. The given embodiments are merely illustrative of the invention and not intended to limit its scope. The embodiments provided below can serve as a guide for further improvements by those skilled in the art and do not constitute a limitation on the invention in any way.

[0067] Unless otherwise specified, the experimental methods used in the following examples are conventional methods, performed according to the techniques or conditions described in the literature in this field or according to the product instructions. Unless otherwise specified, the materials and reagents used in the following examples are commercially available.

[0068] Example 1. Identification of stimulus-response related natural variation sites based on cross-conditional comparative analysis

[0069] 1. The effect of natural variation on gene expression using the Bray-Curtis distance.

[0070] To accurately identify molecular traits associated with stress-related molecular responses, it is necessary to visually represent the impact of natural variation on these traits. This would allow comparisons of molecular traits under control and stress conditions to reflect the differences in stimulus responses between different genotypes caused by natural variation. A valuable metric is the Bray-Curtis distance, which measures differences in species composition across different plots through species abundance. By treating genotypes as "species" and transforming molecular traits using the Bray-Curtis distance, the differences in the impact of natural variation (different genotypes) on molecular traits among different maize inbred lines can be measured.

[0071] 1.1 GWAS analysis of molecular traits to determine haplotypes associated with molecular traits

[0072] To demonstrate whether the Bray-Curtis distance can effectively reflect the differences in the influence of natural variations on molecular traits, this invention performed a GWAS analysis on the expression level of the gene ZmPUM4 in maize inbred lines under salt stress (transcriptome data of 162 maize inbred lines under salt stress and control treatments). The GWAS results of ZmPUM4 showed that it had a significant signal at approximately 192 Mb on chromosome 1 (Figure 2a). Cluster analysis of significant SNPs (P<=1e-7) in the GWAS results revealed that, based on the cluster analysis results, the SNPs associated with the molecular trait of ZmPUM4 gene expression in the maize inbred line population were divided into two haplotypes (SNP combinations) hap1 and hap2. Among them, there were 44 maize inbred lines containing haplotype hap1 and 118 maize inbred lines containing haplotype hap2 (Figure 2b); and there were significant differences in the molecular trait (ZmPUM4 expression level) between the two haplotypes (Figure 2c).

[0073] 1.2 The Bray-Curtis distance matrix can effectively reflect the influence of natural variation on molecular traits.

[0074] By transforming molecular traits using Bray-Curtis distance, the expression levels of ZmPUM4 were transformed using Bray-Curtis distance. The Bray-Curtis distance matrix obtained by transforming ZmPUM4 expression using Bray-Curtis distance was compared with haplotypes related to ZmPUM4 gene expression. It was found that there were large differences in Bray-Curtis distances between inbred lines of different haplotypes, while the differences in Bray-Curtis distances between inbred lines of the same haplotype were not significant (d in Figure 1). This indicates that transforming molecular traits using Bray-Curtis distance can reflect the differences in molecular traits (ZmPUM4 gene expression levels) among different haplotypes under salt stress treatment.

[0075] 2. By comparing the stress response indicators of genes constructed using PCA

[0076] The molecular traits (gene expression levels) transformed using Bray-Curtis distance in step 1 can reflect the influence of natural variation (different haplotypes) on molecular traits. By comparing molecular traits under stimulus and control conditions, this invention can assess the stimulus-response effect of molecular traits. However, these molecular traits transformed using Bray-Curtis distance form a distance matrix, making it difficult to extract differential components from the two matrices under stress and control conditions. This involves feature extraction from the two matrices and comparison of differences between features. Therefore, this invention uses a contrastive PCA algorithm (a variant of contrastive PCA: a contrastive PCA algorithm based on sparse PCA), which can process two datasets simultaneously and automatically extract the most significant differences between the two datasets.

[0077] Contrastive PCA (related literature: Abid, A., Zhang, MJ, Bagaria, VK et al. Exploring patterns enriched in a dataset with contrastive principal component analysis. Nat Commun 9, 2134 (2018). https: / / doi.org / 10.1038 / s41467-018-04608-8) is an algorithm that uses PCA to reduce the dimensionality of the contrastive covariance matrices of a target and background dataset, extracting the differential components between the datasets. Specifically, it is a novel machine learning algorithm that uses a reference dataset as a background to remove noise from the reference data in the target dataset, thereby obtaining the differences between the target and reference datasets. It provides a new analytical approach for studies involving comparisons between datasets. This method first calculates the covariance matrices C_X and C_Y of the target and reference datasets, then filters out noise in the reference dataset by subtracting the two matrices, and finally performs PCA analysis on the difference between the two matrices to obtain the differential components between the two datasets. To control the noise ratio in the reference dataset, a contrast coefficient α is added during the difference calculation of the two matrices. The specific optimization formula for contrastive PCA is shown in equation (1) below:

[0078] ;

[0079] Where V is the principal component of contrastive PCA, used to represent the difference components between datasets.

[0080] The contrastive PCA algorithm based on sparse PCA (related literature: Philippe Boileau, Nima S Hejazi, Sandrine Dudoit, Exploring high-dimensional biological data with sparse contrastive principal component analysis, Bioinformatics, Volume 36, Issue 11, June 2020, Pages 3422–3430, https: / / doi.org / 10.1093 / bioinformatics / btaa176) modifies the PCA algorithm in the dimensionality reduction process, using sparse PCA to reduce the dimensionality of the contrastive covariance matrix. The optimization problem of this process is a PCA problem that penalizes components.

[0081] ;

[0082] in Represents the Frobenius norm. denoted by , it represents the term-by-term matrix norm, which is the sum of the absolute values ​​of all entries in the matrix. U represents decomposition coefficients, and V represents dictionary element. C represents the contrast covariance matrix between the stress treatment and control datasets, and its calculation formula is as follows (4):

[0083] ;

[0084] Where C treat C represents the covariance matrix of molecular traits after transformation under stress treatment. control The covariance matrix of molecular traits after transformation under the control treatment is shown. For comparison parameters.

[0085] In ecology and biology, the Bray-curtis distance is a measure of species diversity to assess differences in species composition across different plots. Since differences in molecular phenotypic values ​​between different inbred lines are caused by allelic differences, natural variations influencing molecular phenotypic values ​​are considered "species." The Bray-curtis distance can quantify genotypic differences between inbred lines with differing molecular phenotypic values, thus reflecting the impact of natural variation on molecular traits. We assume that u and v represent the genotypic data of SNPs influencing molecular phenotypic values ​​in any two inbred lines, which can be represented as an m x 1 vector, where m represents the number of SNPs influencing the molecular phenotypic value. This represents the genotype of the i-th SNP in the first inbred line. Let represent the genotype of the i-th SNP in the second inbred line. Then the genotypic difference between any two inbred lines can be expressed as:

[0086] ;

[0087] Since the phenotypic values ​​of molecular traits reflect the influence of natural variation on phenotypic values, the genotypes of SNPs represented by u and v in the above formula can be replaced by the phenotypic values ​​of molecular traits, and the formula is simplified as follows:

[0088] ;

[0089] Here, u and v represent the molecular phenotypic values ​​between any two self-crossing lines.

[0090] To verify whether PCA can be used to extract molecular trait stimulus response indicators, this invention uses the gene ZmRH3 of the DEAD-box RNA Helicases family in maize, which is associated with abiotic stress, for analysis.

[0091] 2.1 GWAS analysis of molecular traits to determine haplotypes associated with molecular traits

[0092] GWAS analysis was performed on the expression level of ZmRH3 gene in a maize inbred line population. Cluster analysis was conducted on the 36 SNPs (P<=1e-7) that were significantly associated with ZmRH3 gene expression. The clustering results showed that the maize inbred line population could be divided into two haplotypes, ZmRH3-hap1 and ZmRH3-hap2. ZmRH3-hap1 contained 39 inbred lines, and ZmRH3-hap2 contained 123 inbred lines. Furthermore, there was a significant difference in the expression level of ZmRH3 between the two haplotypes.

[0093] 2.2 Transformation of molecular traits via Bray-Curtis distance

[0094] First, the Bray-curtis distance was used to transform the expression levels of ZmRH3 under salt stress (maize seedlings germinating in nutrient soil fully soaked in NaCl solution) and control conditions (using transcriptome data from a population of 162 maize inbred lines measured in our laboratory) (steps as in 1.1), resulting in the Bray-curtis distance matrices of ZmRH3 expression levels under salt treatment and control conditions.

[0095] The results after conversion revealed differences in the effects of different haplotypes (natural variants) on ZmRH3 expression levels (molecular traits) between salt treatment and control conditions. Furthermore, the effect of natural variants on molecular traits was higher under salt treatment than under control conditions, indicating that the effects of different haplotypes on ZmRH3 expression levels differed between salt stress and control conditions. This difference was attributed to the response of natural variants (different haplotypes) to salt stress (Figure 3a).

[0096] 2.3 Comparative dimensionality reduction of molecular traits after Bray-Curtis distance transformation under different stress conditions by comparing PCA

[0097] Subsequently, to more accurately extract the stimulus-response effect of natural variation, this invention used a variant of the contrastive PCA algorithm, sparse contrastive PCA, to perform dimensionality reduction analysis on ZmRH3 expression levels transformed using Bray-Curtis distance under salt stress treatment and control conditions. The analysis revealed that in inbred lines with principal components greater than 0 in the contrastive PCA, there was a difference in ZmRH3 expression levels between the control and treatment; however, in inbred lines with principal components less than 0, there was no significant difference in ZmRH3 expression levels between the control and treatment. This indicates that the dimensionality reduction analysis results (principal components of contrastive PCA) of the Bray-Curtis distance matrix transformed from ZmRH3 expression levels under salt stress treatment and control conditions reflect the differences in ZmRH3 expression levels (molecular traits) in response to salt stress among different inbred lines (different haplotypes caused by natural variation) (Figure 3b).

[0098] 2.4 Perform GWAS analysis on the principal components compared to PCA

[0099] Next, SNP genotyping data obtained from maize inbred line population resequencing data were received. GWAS analysis was performed on the principal components of the contrastive PCA (using the principal components of the contrastive PCA as stimulus-response indices for molecular traits). This revealed that it could identify significantly associated sites for salt stress, and the peak SNP (significantly associated SNP) was located at 1,991,212 bp on chromosome 5 of the maize genome, precisely on the intron of ZmRH3 (Figure 3c). Further haplotype analysis using the peak SNP (the inbred lines were reclassified based on the major and minor alleles of the peak SNP) revealed different haplotypes affecting ZmRH3 expression levels with varying responses under the two conditions (Figure 3d). This indicates the presence of cis-variant sites near the ZmRH3 gene that regulate the salt stress response of ZmRH3, and also demonstrates that the method of this invention using the contrastive PCA algorithm can mine naturally occurring variant sites related to stimulus responses (such as stress responses).

[0100] 2.5 Comparison of the applicability of PCA algorithm for dimensionality reduction across the entire genome.

[0101] To verify the applicability of contrastive PCA across the entire genome, this invention further first performed Bray-Curtis distance transformation on the expression levels of 29,552 genes (29,552 molecular traits) in the transcriptome data of maize inbred lines under salt stress and control conditions (including population transcriptome data of 162 maize inbred lines), obtaining 29,552 pairs of Bray-Curtis distance matrices under stress and control conditions. Then, contrastive PCA dimensionality reduction analysis was performed on each pair of Bray-Curtis distance matrices, yielding 59,104 contrastive PCA principal components related to salt stress response.

[0102] 2.6 Significance analysis of principal components compared with PCA

[0103] Since contrastive PCA does not provide a p-value for statistical significance, this invention further develops a significance analysis method for evaluating the correlation between contrastive PCA principal components and stimulus response differences.

[0104] Specifically, this invention uses least-squares linear regression to calculate the regression coefficients between the principal components of the comparative PCA and gene expression levels under control and stimulus conditions (such as salt stress treatment). By comparing the differences in their effects on gene expression levels under the two conditions, the significance of the differences in stimulus responses reflected by the principal components of the comparative PCA is assessed.

[0105] By analyzing the differences in regression coefficients between salt stress treatment and control conditions, this invention found 32,494 comparative PCA principal components, whose regression coefficients with gene expression levels under salt stress treatment and control conditions showed significant differences. This indicates that comparative PCA principal components can effectively represent the natural variation in transcriptome data that affects the differences in gene salt stress response (Figure 4a).

[0106] However, this invention also notes that there was no significant difference between the regression coefficients of 26,610 principal components based on contrastive PCA and gene expression levels under salt stress and control conditions, which is understandable under normal circumstances (Figure 4b). Since two principal components were used for analysis of each gene, some noise may have been included. The purpose of using two principal components is to more comprehensively explore the natural variations influencing gene responses to salt stress. Therefore, this invention defines the principal components of contrastive PCA that show significant differences in regression coefficients between gene expression levels under stress and control conditions (significant contrastive PCA principal components) as stress response indicators for subsequent identification of natural variation sites associated with stress response.

[0107] 3. Use stress response indicators to perform GWAS analysis to screen QTLs related to stress response.

[0108] To verify whether salt stress response indicators can be used for the identification of salt-tolerant QTLs, this invention performed GWAS analysis on 32,494 salt stress response indicators (compared to PCA principal components) obtained in the maize transcriptome in step 2.6 as traits. The results identified 4,399 QTL sites associated with salt stress response, belonging to 3,531 genes.

[0109] Because some salt stress response indicators show weak salt stress responses, the differences in salt stress responses among different QTL alleles may not be significant. Therefore, this invention further uses haplotypes represented by QTL peak SNPs and salt stress treatment conditions as factors to conduct a two-way ANOVA on the genes containing identified QTLs to assess the strength of QTL salt stress responses. The two-way ANOVA results identified significant differences in the responses of different alleles of 2666 QTLs to salt stress. Analysis of the expression levels of genes containing these QTLs revealed that these genes significantly respond to salt stress under one of the allele types, and the detailed trends can be divided into 5 categories (Figure 5a).

[0110] Gene ontology (GO) enrichment analysis was performed on genes containing QTLs with significant two-way ANOVA p-values. The analysis revealed that these genes were significantly enriched with terms related to signal transduction and abiotic stress responses, including terms related to osmotic stress responses, salt stress responses, and intracellular potassium ion transport (Figure 5b). These results indicate that the salt stress response indices obtained by the method of this invention can be used for the identification of salt-tolerant QTLs, but screening using two-way ANOVA is required.

[0111] 4. Identification of causal genes associated with salt tolerance phenotype based on Mendelian randomization analysis

[0112] Based on the 2,666 salt tolerance QTLs identified and previously measured salt tolerance phenotype data, this invention uses Mendelian randomization analysis to assess the contribution of genes containing salt tolerance QTLs to the salt tolerance phenotype, thereby identifying causal (mutated) genes associated with the salt tolerance phenotype.

[0113] In Na +In the Mendelian randomization results of Na+ content, ZmHAK4 had the most significant contribution. Haplotype analysis revealed that under minor alleles, ZmHAK4 was upregulated under salt stress, and maize inbred lines had lower Na+ content (Figure 5c and d). This is consistent with the results reported in the literature: Zhang et al. found that ZmHAK4 increases the salt tolerance of maize by excluding Na+ from the shoot (related literature: Zhang, M., Liang, X., Wang, L. et al. A HAK family Na+ transporter confers natural variation of salt tolerance in maize. Nat. Plants 5, 1297–1308 (2019). https: / / doi.org / 10.1038 / s41477-019-0565-y).

[0114] Therefore, the method of the present invention can effectively identify stimulus-response related natural variant sites, stimulus-response related QTLs or genes, and apply them to crop stress tolerance breeding to improve crop stress tolerance.

[0115] The present invention also provides a computer device for identifying stimulus-response related natural variant sites, the computer device including a memory, a processor, and a computer program stored in the memory, the processor executing the computer program to perform the following steps:

[0116] B1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; receive SNP genotyping data of the genome of each individual in the genetic population;

[0117] The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition.

[0118] B2) Data Processing: The expression levels of each gene in the genetic population under the stimulation condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression levels of each gene in the genetic population under the control condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA.

[0119] Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association (GWAS) analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators.

[0120] B3) Haplotype analysis: The expression level of each gene in the genetic population is taken as a molecular trait. Based on the transcriptome data under control conditions and the SNP genotyping data of each individual genome in the genetic population, genome-wide association (GWAS) analysis is performed on the molecular trait to obtain the significantly associated SNPs of the molecular trait. Cluster analysis is performed on the significantly associated SNPs of the molecular trait to obtain the haplotypes associated with the molecular trait.

[0121] B4) QTL analysis: Based on the significant association between stimulus response indicators and SNPs, candidate QTLs related to stimulus response are obtained. Using molecular trait-related haplotypes and the stimulus conditions as factors, a two-way ANOVA is performed on the candidate QTLs related to stimulus response to obtain the stimulus response-related QTLs.

[0122] B5) Mendelian randomization analysis: Receive stimulus-phenotype related data from a genetic population, and based on stimulus-response related QTLs and stimulus-phenotype related data, use Mendelian randomization analysis to obtain the contribution of genes containing stimulus-response related QTLs to the stimulus phenotype, and obtain causal genes related to the stimulus phenotype based on the contribution.

[0123] The present invention has been described in detail above. Those skilled in the art will recognize that the invention can be practiced in a wide range of ways with equivalent parameters, concentrations, and conditions without departing from its spirit and scope, and without requiring unnecessary experiments. While specific embodiments have been provided, it should be understood that further modifications can be made to the invention. In summary, according to the principles of the invention, this application is intended to include any changes, uses, or improvements to the invention, including changes made using conventional techniques known in the art that depart from the scope disclosed herein. Industrial applicability

[0124] This invention transforms molecular traits using Bray-Curtis distance, enabling the transformed molecular traits to represent the influence of natural variation on molecular traits, including the genetic effects of natural variation and stress response effects. Subsequently, comparative PCA is used to perform dimensionality reduction on molecular traits under stress treatment and control, accurately obtaining the stress response effect of molecular traits, which is then used as a stress response index. Finally, genome-wide association analysis of the stress response index yields molecular QTLs of stress response, used for the identification of cross-conditional molecular response QTLs such as those related to stress stimuli. This invention, when analyzing stress-related population omics data, can consider both the genetic effects of natural variation on molecular traits and the stress response effects of natural variation on molecular traits, effectively analyzing stress response differences between different alleles, thereby accurately identifying population-level molecular QTLs related to stimulus responses. This can be applied to genome-wide selection breeding and quality improvement of crops.

[0125] Cross-reference to related applications

[0126] This application claims priority to Chinese patent application No. 202411982410.0, filed on December 31, 2024, the entire contents of which are incorporated herein by reference.

Claims

1. A computer device comprising a memory, a processor, and a computer program stored in the memory, characterized in that, The processor executes the computer program to perform the following steps: S1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; Receive SNP genotyping data of the genome of each individual in the genetic population; The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition. S2) Data Processing: The expression level of each gene in the genetic population under the stimulation condition is transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression level of each gene in the genetic population under the control condition is transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA. Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators. S3) Data output: Output the significant association SNPs.

2. The computer device according to claim 1, characterized in that: The significance analysis in S2) includes the following steps: The first regression coefficient is obtained by calculating the regression coefficient between the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition. The second regression coefficient is obtained by calculating the regression coefficient between the principal components of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition. The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

3. The computer device according to claim 2, characterized in that: The PCA algorithm is a contrastive PCA algorithm based on sparse PCA.

4. The computer device according to claim 3, characterized in that: The genetic population is a plant genetic population, and the stimulus is a stress treatment.

5. The computer device according to claim 1, characterized in that: The PCA algorithm is a contrastive PCA algorithm based on sparse PCA.

6. The computer device according to claim 1, characterized in that: The genetic population is a plant genetic population, and the stimulus is a stress treatment.

7. A device for identifying stimulus-response related QTLs, characterized in that: The device includes the following modules: A1) Data receiving module: used to receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; and to receive SNP genotyping data of the genome of each individual in the genetic population; The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition. A2) Data Processing Module: This module is used to convert the expression level of each gene in the genetic population under the stimulation condition using Bray-Curtis distance to obtain a Bray-Curtis distance matrix under the stimulation condition; convert the expression level of each gene in the genetic population under the control condition using Bray-Curtis distance to obtain a Bray-Curtis distance matrix under the control condition; perform dimensionality reduction on the difference between the Bray-Curtis distance matrix of each gene under the stimulation condition and the Bray-Curtis distance matrix under the control condition using the PCA algorithm to obtain the principal components of candidate contrastive PCA; and perform significance analysis on the principal components of candidate contrastive PCA to obtain the principal components of contrastive PCA. Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators. A3) Haplotype Analysis Module: This module uses the expression level of each gene in the genetic population as a molecular trait. Based on the transcriptome data under the control conditions and the SNP genotyping data of each individual's genome in the genetic population, it performs genome-wide association analysis on the molecular trait to obtain the significantly associated SNPs of the molecular trait. Then, it performs cluster analysis on the significantly associated SNPs of the molecular trait to obtain the haplotypes associated with the molecular trait. A4) QTL Analysis Module: Used to obtain stimulus-response-related candidate QTLs based on the significant association SNPs of the stimulus response index, and to perform a two-way ANOVA on the stimulus-response-related candidate QTLs using the molecular trait-related haplotypes and the stimulus conditions as factors, to obtain the stimulus-response-related QTLs.

8. The apparatus according to claim 7, characterized in that: The significance analysis described in A2) includes the following steps: The first regression coefficient is obtained by calculating the regression coefficient between the principal component of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition. The second regression coefficient is obtained by calculating the regression coefficient between the principal component of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition. The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

9. A method for identifying causal genes related to stimulus phenotypes, characterized in that: The method includes the following steps: B1) Data reception: Receive transcriptome data of the genetic population under stimulation conditions and transcriptome data under control conditions; receive SNP genotyping data of the genome of each individual in the genetic population; The transcriptome data under the stimulation condition includes the expression levels of all genes in the individuals of the genetic population under the stimulation condition; the transcriptome data under the control condition includes the expression levels of all genes in the individuals of the genetic population under the control condition. B2) Data Processing: The expression levels of each gene in the genetic population under the stimulation condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the stimulation condition; the expression levels of each gene in the genetic population under the control condition are transformed using Bray-Curtis distance to obtain the Bray-Curtis distance matrix under the control condition; the difference between the Bray-Curtis distance matrix under the stimulation condition and the Bray-Curtis distance matrix under the control condition for each gene is reduced in dimensionality using the PCA algorithm to obtain the principal components of candidate contrastive PCA; the principal components of candidate contrastive PCA are subjected to significance analysis to obtain the principal components of contrastive PCA. Using the principal components of the comparative PCA as stimulus response indicators, and based on the SNP genotyping data of each individual's genome in the genetic population, genome-wide association analysis was performed on the stimulus response indicators to obtain the significantly associated SNPs of the stimulus response indicators. B3) Haplotype analysis: The expression level of each gene in the genetic population is taken as a molecular trait. Based on the transcriptome data under the control conditions and the SNP genotyping data of each individual genome in the genetic population, a genome-wide association analysis is performed on the molecular trait to obtain the significantly associated SNPs of the molecular trait. The significantly associated SNPs of the molecular trait are then clustered to obtain the haplotypes associated with the molecular trait. B4) QTL analysis: Based on the significant association of the stimulus response index with SNP, candidate QTLs related to stimulus response are obtained. Using the haplotypes related to the molecular traits and the stimulus conditions as factors, a two-way ANOVA is performed on the candidate QTLs related to stimulus response to obtain the stimulus response-related QTLs. B5) Mendelian randomization analysis: Receive stimulus-phenotype related data of the genetic population, and based on the stimulus-response related QTLs and the stimulus-phenotype related data, obtain the contribution of genes containing the stimulus-response related QTLs to the stimulus phenotype through Mendelian randomization analysis, and obtain causal genes related to the stimulus phenotype based on the contribution.

10. The method according to claim 9, characterized in that: The significance analysis described in B2) includes the following steps: The first regression coefficient is obtained by calculating the regression coefficient between the principal component of the candidate contrastive PCA for each gene and the expression level of each gene under the stimulation condition. The second regression coefficient is obtained by calculating the regression coefficient between the principal component of the candidate contrastive PCA for each gene and the expression level of each gene under the control condition. The first regression coefficient and the second regression coefficient are subjected to a significance analysis, and the principal components of the candidate contrastive PCA that show significant differences are retained to obtain the principal components of the contrastive PCA.

11. A computer-readable storage medium storing a computer program, characterized in that: The computer program causes the computer to perform steps B1)-B2) of the method as described in claim 9 or 10.

12. A computer-readable storage medium storing a computer program, characterized in that: The computer program causes the computer to perform steps B1)-B4) of the method as described in claim 9 or 10.

13. A computer-readable storage medium storing a computer program, characterized in that: The computer program causes the computer to perform steps B1)-B5) of the method as described in claim 9 or 10.

14. The use of the computer device according to any one of claims 1-6, or the device according to claim 7 or 8, in plant breeding.

15. The use of the computer-readable storage medium according to any one of claims 11-13 in plant breeding.

Images