Anxiety disorders are common mental disorders (Huang et al., 2019) and a leading cause of disability (GBD 2017 Disease and Injury Incidence and Prevalence Collaborators, 2018). As an adaptive behavioral response, anxiety is usually a normal emotional experience in daily life. However, individuals affected by anxiety disorders frequently have feelings of excessive, intense, and persistent worry and fear about the actual or anticipated event. If the degree of anxiety exceeds a certain range and the abnormal feeling of fear and worry lasts for a long time, normal anxiety experiences may turn into anxiety disorders. Anxiety disorders are highly prevalent mental disorders (with a lifetime prevalence rate of over 20% (Kessler et al., 2012)) that bring huge economic burdens to the world (Baxter et al., 2014). The etiology of anxiety disorders remains poorly understood. However, it is found that genetic factors play pivotal roles. The heritability of anxiety disorders were estimated to range from 20 to 60% (Polderman et al., 2015; Craske et al., 2017), indicating substantial genetic components of anxiety disorders. To uncover the genetic basis of anxiety disorders/traits, several genetic studies have been performed (Sipilä et al., 2010; Stein et al., 2011; Nivard et al., 2014). Nevertheless, large-scale GWASs of anxiety disorders/traits have not emerged until recently (Otowa et al., 2016; Meier et al., 2019; Levey et al., 2020; Purves et al., 2020). For example, by using the UK biobank sample, Purves et al. identified five genome-wide significant (GWS) risk loci for anxiety disorder (Purves et al., 2020). In addition, Meier et al. reported that PDE4B is associated with anxiety and stress-related diagnoses in the Danish population (Meier et al., 2019). Otowa et al. also found that a noncoding RNA region located in 3q12.3 is associated with anxiety disorder (Otowa et al., 2016). Notably, a recent large-scale GWAS (N = 199,611) conducted by Levey et al. (i.e., the Million Veteran Program, MVP) reported six GWS loci for Generalized Anxiety Disorder 2-item scale (GAD-2) GWAS (Levey et al., 2020). There are also other GWASs for specific anxiety disorders/traits, such as Generalised Anxiety Disorder (Davies et al., 2015; Dunn et al., 2017)^, panic disorder (Otowa et al., 2009; Erhardt et al., 2011; Otowa et al., 2012), phobic anxiety (Walter et al., 2013), and social anxiety (Stein et al., 2017). Although these GWASs have identified several risk variants and loci for anxiety disorders/traits, the number of reported anxiety disorders/traits risk loci remains small compared with other psychiatric disorders (including schizophrenia (Schizophrenia Working Group of the Psychiatric Genomics Consortium, 2014; Pardiñas et al., 2018; Lam et al., 2019), depression (Wray et al., 2018), and bipolar disorder (Stahl et al., 2019)).

To date, most of the risk variants reported by anxiety disorders/traits GWASs are located in noncoding regions (Otowa et al., 2016; Meier et al., 2019; Levey et al., 2020; Purves et al., 2020), suggesting that these genetic variations confer risk of anxiety disorders/traits by modulating gene expression. However, we currently know little about the target genes regulated by these reported risk variants. In addition, the joint effects of the reported anxiety disorders risk loci on gene expression across human brain tissues remain elusive. The emergence of TWAS (Gusev et al., 2016) provides an opportunity to infer genes whose cis-regulated expression may be associated with diseases and traits. TWAS first leverages genotype and expression data from an external reference panel to identify the associations between genetic variations and gene expression. Then it identifies genes whose cis-regulated expression are associated with diseases or phenotypes by integrating gene expression data (i.e., external reference panel) with genome-wide associations from large-scale GWASs. TWASs have been widely used to identify risk genes for psychiatric disorders and neurodegenerative disease, including schizophrenia (Gandal et al., 2018), depression (Dall'Aglio et al., 2020), Parkinson’s disease (Li et al., 2019), and attention deficit hyperactivity disorder (Liao et al., 2019). Moreover, TWASs have identified multiple novel risk genes that were not implicated by the original GWASs (Liao et al., 2019).

In this study, we performed a TWAS (using FUSION software (Gusev et al., 2016)) on anxiety traits by integrating genome-wide associations from a recent large-scale GAD-2 score GWAS of anxiety disorder (Levey et al., 2020) (i.e., the Million Veteran Program) and gene expression data from multiple brain tissues. We identified 26 TWS genes, of which 19 were from the PsychENCODE gene expression reference and 17 were from GTEx. It is noteworthy that 10 genes (KANSL1-AS1, CRHR1, CRHR1-IT1, SPPL2C, RP11-707O23.5, RP11-259G18.1, MAPT-AS1, LRRC37A4P, PLEKHM1, and DND1P1) showed TWS in both datasets (i.e., PsychENCODE and GTEx), implying these genes are promising candidate genes for anxiety. Moreover, we noticed that some novel risk genes identified by our TWAS were not implicated in the original GWAS. Our study identifies genes whose expression change may confer risk of anxiety disorders/traits, shading lights on the genetic component of gene expression in anxiety disorders/traits, and providing a start point for mechanistic investigation and future drug development.

In this study, we reported the first TWAS of anxiety to identify genes whose genetically regulated expression may confer risk of anxiety. Our study identified 26 TWS genes, of which 25 were not nominated in the original GWAS, indicating that TWAS is a powerful approach to nominate novel risk genes that are not implicated by GWAS. And 10 TWS genes were identified in both PsychENCODE and GTEx datasets. Among the 10 overlapping genes, we found that LRRC37A4P is widely expressed in several neuronal cell types. None of these genes showed cell-type specific expression (Polioudakis et al., 2019) (Supplementary Figure S15–S20). Furthermore, we also found that genes located in 3p21.31 and 17q21.31 are significantly associated with anxiety. We noticed that these two loci also showed significant associations with other psychiatric disorders. For example, MAPT gene (located in 17q21.31) is a well-known candidate risk gene for Parkinson’s disease (Zabetian et al., 2007; Simón-Sánchez et al., 2009; Edwards et al., 2010). In addition, the 3p21.31 and 17q21.31 are also significantly associated with Post-traumatic stress disorder (PTSD) in a recent GWAS study by the MVP (Gelernter et al., 2019). Intriguingly, a recent TWAS on PTSD also indicated that genes located in these two regions are related to PTSD (Girgenti et al., 2021). Notably, RNF123 was one of the top genes that were significantly associated with PTSD at the transcriptome-wide significance level (Girgenti et al., 2021). The shared TWAS findings may be due to the high genetic correlation (approximately 50–60%) between PTSD and anxiety (Stein et al., 2021). In addition, RNF123 was also reported to be dysregulated in depression and may serve as a clinical biomarker for depression (Glahn et al., 2012; Teyssier et al., 2013). These results collectively indicate that RNF123 may play an important role in stress-related disorders, including anxiety, PTSD, and depression. Though anxiety disorder, PTSD, and depression are diagnosed as different mental disorders (as they represent potentially fundamentally different dimensional endophenotypes (Insel and Cuthbert, 2009)), the genetic risk factors that these diseases shared may represent general genetic risk factors of these mental diseases (Smoller, 2016). Taken together, these lines of evidence indicate that these two genomic regions may harbor authentic risk genes for anxiety and genetic variants likely confer risk of anxiety by modulating the expression of genes located in 3p21.31 and 17q21.31.

RNF123 (Ring Finger Protein 123) encodes E3 ubiquitin-protein ligase which is a subunit of the Kip1 ubiquitination-promoting complex (KPC). Previous studies revealed the important function of RNF123, including regulation of cell cycle process (Kamura et al., 2004) and innate antiviral signaling (Wang et al., 2016). Besides, several studies also showed that RNF123 may have a role in cancer, including aggressive glioblastoma and melanoma (Iida et al., 2017; Wang et al., 2020). However, the function of RNF123 in the human central nervous system and the corresponding mechanisms about how it confers the risk of mental disorders remain to be investigated.

Two major reasons might account for the observations that some of the loci were significant in original GWAS but not in TWAS. Firstly, TWAS was designed to detect the gene-trait association by performing transcriptome imputation of GWAS data with a reference trained gene expression predictive model (constructed from a dataset with both genotype and gene expression) (Gusev et al., 2016; Wainberg et al., 2019). If the genetic variations in the GWS loci are not functionally as cis-regulation variations (which affect gene expression), TWAS is not able to nominate candidate genes in these loci. Secondly, the power of the TWAS may be limited by the sample size of the training set we used, e.g., the psychENCODE dataset we used in this study can pinpoint more TWS genes than a single GTEx panel. We also identified loci that were not significant in original GWAS but significant in TWAS. In TWAS analysis, we can identify TWS loci that are not nominated by GWS, which is the power and advantages of TWAS (Gusev et al., 2016; Liao et al., 2019). For TWS loci not reported by GWS, a possible reason is that the statistic power is not enough in the current GWAS. For example, the leading GWAS association signal in TWS locus (17q21.31) is rs142925250 (5.95E-08, Table 1), which is close to GWS level (5.00E-08). With the increase of sample size, rs142925250 may reach the genome-wide significance level. Finally, we also identified loci that reached GWS and TWS simultaneously. The GWAS result is not easy to interpret because of complex LD (linkage disequilibrium) (several linked genes may show associations with a GWAS trait). TWAS leverages a reference expression panel (samples with gene expression and genotype data simultaneously) to perform transcriptome imputation on a target GWAS dataset, and the predicted gene expression of the GWAS dataset was used to identify the gene-trait relationship. So TWAS can help to nominate gene-trait relationships by integrating eQTL data and GWAS data, which may help to explain the functions of the GWS loci (Wainberg et al., 2019).

As the predicted expression change of the TWS genes may contribute to disease risk, we also explored the expression of TWS genes in individuals with anxiety disorder and controls, using the published differential gene expression study by Wingo et al. (Wingo and Gibson, 2015) Gene expression in the blood of 157 anxiety disorder cases and 179 controls (GSE61672) were measured by Wingo et al (Wingo and Gibson, 2015). We performed differential expression analysis using NCBI GEO2R (Barrett et al., 2013), and found no overlap between the differentially expressed genes (DEGs) and our TWAS findings. We hypothesized that this might be due to the tissue difference between DEG analysis (blood) and TWAS (brain).

We performed MAGMA analysis for two purposes. First, the original GAD-2 GWAS papers did not perform MAGMA gene-set enrichment analysis to interpret the GWAS results (Levey et al., 2020). We therefore performed gene-set enrichment analysis by MAGMA (de Leeuw et al., 2015) to explore if genome-wide associations were enriched in specific tissues or cell types (using expression data from GTEx and single-cell RNA sequencing dataset). Second, we aimed at finding out the most-relevant tissues for TWAS analysis when using GTEx reference panels. Our MAGMA tissue enrichment analysis indicated that the brain tissues in GTEx were suitable for performing TWAS. As expected, we found that the anxiety associations were enriched in brain tissues and neurons. Of note, significant enrichment in the testis was also observed. A possible reason for this observation is that only limited genome-wide associations were reported by the original GWAS. Though the sample size of the MVP is relatively large (N = 175,163), we noticed that only five genome-wide significant (GWS) loci have been identified in the original GWAS of anxiety. The number of GWS risk loci for anxiety is quite small compared with other psychiatric disorders such as schizophrenia and depression. In addition to sample size, we found that the heritability were different among anxiety disorders (20–60%) (Polderman et al., 2015; Craske et al., 2017), SCZ (approximately 80%) (Sullivan et al., 2003), and depression (30–40%) (Sullivan et al., 2000). The heritability estimated from the GWAS is lower than the true heritability as GWAS can only capture the common variations (Manolio et al., 2009). The heritability of anxiety traits estimated from GAD-2 based GWAS is about 5.58% (n = 175,163) (Levey et al., 2020), the heritability of schizophrenia explained by GWAS is approximately 23–24% (n = 135,236) (Lam et al., 2019), and the heritability of depression estimated from GWAS is about 8.7% (n = 480,359) (Wray et al., 2018). These results suggest that the polygenetic nature varies among these diseases. Second, gene expression patterns are highly similar between human brain and testis (Guo et al., 2005), which may result in the observation of significant enrichment of anxiety associations in testis. Finally, it is possible that the anxiety risk genes also play a role in testis. More work is needed to explain this interesting observation.

Our study highlighted the importance of TWAS (which integrates gene expression data with GWAS summary) in nominating risk genes for anxiety. However, genetics can only explain a small proportion of anxiety disorders/traits, which indicates that environmental factors are also involved in anxiety disorders/traits. Among the various environmental factors, psychosocial stress plays a vital role in anxiety (Bartlett et al., 2017). For example, early life adversity (Nemeroff CB, 2004; Maccari et al., 2014; Wiegand et al., 2021) is an important environmental risk factor for social anxiety disorder. The stress environmental factors like early life adversity would trigger the DNA methylation change (Provençal and Binder, 2015), and epigenome-wide association study (EWAS) have reported many genes that showed epigenetic changes in anxiety disorders, including CFAP46 (Ziegler et al., 2019), SLC43A2, and TNXB (Wiegand et al., 2021). Cognitive-behavioral therapy (CBT) (Barlow et al., 2000) is an useful psychosocial intervention for anxiety disorders and previous study has reported that the methylation level of IL1R1 was changed in subjects of panic disorder (one of the anxiety disorders) treated with CBT (Ziegler et al., 2019). These findings indicate that epigenetic regulation also plays an import role in anxiety disorders/traits.

Our study has several limitations. Firstly, we only included the GAD-2 score based GWAS summary statistics from the MVP project in our current study. We compared the phenotyping methods of this GWAS and other published anxiety GWASs (Otowa et al., 2009; Erhardt et al., 2011; Otowa et al., 2012; Walter et al., 2013; Davies et al., 2015; Otowa et al., 2016; Dunn et al., 2017; Stein et al., 2017; Meier et al., 2019; Levey et al., 2020; Purves et al., 2020), and found that the phenotyping methods were not identical across studies. We believed that it would be best to perform an independent replication study based on the same phenotyping method. Secondly, though several TWS genes are identified, more efforts are needed to uncover the roles of these identified risk genes in anxiety etiology. Thirdly, although TWAS could identify candidate risk genes for anxiety, the causal variants in the GWAS/TWAS significant loci could not be pinpointed by TWAS. More efforts are needed to reveal the genetic mechanisms and pathogenesis of the genes in the GWAS/TWAS significant loci in anxiety.

In summary, our study uncovered the gene-trait relationships of anxiety traits for the first time. We identified 26 TWS genes for anxiety. Our results provide novel insight into anxiety disease etiology and facilitate further mechanism studies and future drug development.