Identifying pathways and networks associated with the SARS-CoV-2 cell receptor ACE2 based on gene expression proles in human tissues

The angiotensin-converting enzyme 2 (ACE2) is a host cell receptor of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) that has infected more than six million people worldwide and has caused more than 370,000 deaths as of May 31, 2020. An investigation of ACE2 expression in human tissues may provide insights into the mechanism of SARS-CoV-2 infection. We identied pathways associated with ACE2 expression and gene co-expression networks of ACE2 in pan-tissue based on the gene expression proles in human tissues. We found that the pathways signicantly associated with ACE2 upregulation were mainly involved in immune, stromal signature, metabolism, cell growth and proliferation, and cancer and other diseases. The number of genes having a signicant positive expression correlation with ACE2 in females far exceeded that in males. The estrogen receptors (ESR1 and ESR2) and androgen receptor (AR) genes had a signicant positive expression correlation with ACE2 in pan-tissue. Meanwhile, the enrichment levels of immune cells were positively associated with the expression levels of ESR1 and ESR2, while they were inversely associated with the expression levels of AR in pan-tissue and in multiple individual tissues. It suggests that females are likely to have a more robust immune defense system against SARS-CoV-2 than males, partially explaining why females have better clinical outcomes of SARS-CoV-2 infections than males. Our data warrant further investigation for understanding the mechanism of SARS-CoV-2 infection. hematopoietic cell lineage, Jak-STAT signaling, adipocytokine signaling, chemokine signaling, viral myocarditis, intestinal immune network for IgA production, systemic lupus erythematosus, pathogenic Escherichia coli infection, epithelial cell signaling in Helicobacter pylori infection, and NOD-like receptor signaling. The stromal signature-related pathways included focal adhesion, ECM-receptor interaction, cell adhesion molecules, tight junction, adherens junction, regulation of actin cytoskeleton, axon guidance, and gap junction. The metabolism-related pathways included PPAR signaling, insulin signaling, arachidonic acid metabolism, drug metabolism-cytochrome P450, glutathione metabolism, retinol metabolism, fatty acid metabolism, tyrosine metabolism, metabolism of xenobiotics by cytochrome P450, glycolysis/gluconeogenesis, glycerophospholipid metabolism, phenylalanine metabolism, butanoate metabolism, and glycerolipid metabolism. The cell growth and proliferation-related pathways included MAPK, ErbB, p53, Wnt, VEGF, Notch, and mTOR signaling. The cancer-related pathways included pathways in cancer, small cell lung


Background
The outbreak of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) has caused a global pandemic (1). This virus has infected nearly 6 million people and has caused more than 350,000 deaths as of May 27, 2020 (2). The angiotensin-converting enzyme 2 (ACE2) is a host cell receptor of SARS-CoV-2 that plays a crucial role in regulating the SARS-CoV-2 invasion (3)(4)(5)(6), which in turn can result in ACE2 up regulation (7). Many studies have investigated ACE2 expression in human tissues (8)(9)(10)(11)(12). Our recent study showed that ACE2 is expressed in a wide variety of human tissues (8). It suggests that SARS-CoV-2 may attack various human organs, although patients infected with SARS-CoV-2 primarily displayed pneumonia-associated symptoms (13,14). In this study, we identi ed pathways associated with ACE2 expression and gene co-expression networks of ACE2 in pan-tissue using the datasets from the Genotype-Tissue Expression (GTEx) project (15).

Datasets
We downloaded the gene expression pro ling datasets (RNA-Seq, TPM normalized) for human tissues from the GTEx Project (https://www.gtexportal.org/home/datasets/). All gene expression values were added to 1 and then log2-transformed before subsequent analysis.

Gene-set enrichment analysis
We performed pathway analyses of the differentially expressed genes between the high-ACE2-expressionlevel (upper third) and the low-HIF1A-expression-level (bottom third) samples in pan-tissue. The differentially expressed genes were identi ed by Student's t-test using a threshold of FDR < 0.05 and fold change > 2. The FDR represented the adjusted P-value calculated by the Benjamini and Hochberg method (23). We identi ed the KEGG (24) pathways highly enriched in both groups using GSEA (16) with a threshold of FDR < 0.05. We used WGCNA (25) to identify the gene modules (GO) that were highly enriched in the high-ACE2-expression-level and the low-ACE2-expression-level samples in pan-tissue.

Evaluation of the immune cell enrichment levels in tissue
We determined the enrichment level of an immune signature in a sample as the mean expression level of marker genes of the immune signature in the sample. A total of three immune signatures were analyzed, including B cells, CD8+ T cells, and NK cells. We presented the marker genes of these immune cells in Supplementary Table S5.

Statistical analysis
We used Pearson's correlation test to calculate the expression correlations between ACE2 and other genes and the correlations between the expression levels of sex hormone receptor genes (ESR1, ESR2, and AR) and the enrichment levels of immune cells (B cells, CD8+ T cells, and NK cells) in pan-tissue.

Results
Pathways and gene ontology associated with ACE2 expression We identi ed highly enriched KEGG pathways in pan-tissue having high ACE2 expression levels (upper third) versus low ACE2 expression levels (bottom third) (Student's t-test, adjusted P-value FDR < 0.05, fold change > 2) by GSEA (16) with a threshold of FDR < 0.05. These pathways were mainly involved in immune, stromal signature, metabolism, cell growth and proliferation, cancer, and other diseases (Fig.  1A). The immune-related pathways included cytokine-cytokine receptor interaction, leukocyte transendothelial migration, complement and coagulation cascades, TGF-β signaling, hematopoietic cell lineage, Jak-STAT signaling, adipocytokine signaling, chemokine signaling, viral myocarditis, intestinal immune network for IgA production, systemic lupus erythematosus, pathogenic Escherichia coli infection, epithelial cell signaling in Helicobacter pylori infection, and NOD-like receptor signaling. The stromal signature-related pathways included focal adhesion, ECM-receptor interaction, cell adhesion molecules, tight junction, adherens junction, regulation of actin cytoskeleton, axon guidance, and gap junction. The metabolism-related pathways included PPAR signaling, insulin signaling, arachidonic acid metabolism, drug metabolism-cytochrome P450, glutathione metabolism, retinol metabolism, fatty acid metabolism, tyrosine metabolism, metabolism of xenobiotics by cytochrome P450, glycolysis/gluconeogenesis, glycerophospholipid metabolism, phenylalanine metabolism, butanoate metabolism, and glycerolipid metabolism. The cell growth and proliferation-related pathways included MAPK, ErbB, p53, Wnt, VEGF, Notch, and mTOR signaling. The cancer-related pathways included pathways in cancer, small cell lung cancer, bladder cancer, melanoma, prostate cancer, glioma, acute and chronic myeloid leukemia, basal cell carcinoma, thyroid cancer, pancreatic cancer, renal cell carcinoma, and endometrial cancer, and the other diseases-related pathways included dilated cardiomyopathy, hypertrophic cardiomyopathy, arrhythmogenic right ventricular cardiomyopathy, and prion diseases.
WGCNA generated 11 gene modules (indicated in green-yellow, cyan, yellow, light green, magenta, red, salmon, green, grey, black, and turquoise color, respectively) that were more highly enriched in the high-ACE2-expression-level than in the low-ACE2-expression-level pan-tissue (Fig. 1B). In contrast, two gene modules (indicated in blue and light-yellow color, respectively) were more highly enriched in the low-ACE2expression-level pan-tissue (Fig. 1B). The gene ontology (GO) terms in the highly enriched gene modules in the high-ACE2-expression-level pan-tissue were mainly associated with immune response, cell cycle, reproductive process, brush border assembly, skin development, thyroid hormone metabolic process, muscle system process, cell projection organization, muscle contraction, blood vessel development, and cell growth and proliferation. The GO terms in the highly enriched gene modules in the low-ACE2expression-level pan-tissue were mainly associated with nervous system development and muscle system process. The positive associations of ACE2 expression with immune response, cell cycle, and cell growth and proliferation in pan-tissue were consistent with the pathway analysis results.
Gene co-expression networks of ACE2 We found 2,983 and 74 genes having a signi cant positive and a signi cant negative expression correlation with ACE2 in pan-tissue, respectively (Pearson correlation coe cient |r| > 0.3) (Supplementary  Table S1). Interestingly, when we analyzed female and male pan-tissue individually, we found that 3, 940 (or 587) and 87 (or 93) genes having a signi cant positive and a signi cant negative expression correlation with ACE2 in female (or male) pan-tissue, respectively (|r| > 0.3) (Supplementary Tables S2&S3). It indicates that more genes have a signi cant positive expression correlation with ACE2 in females than in male pan-tissue. Strikingly, we found 77 genes showing a signi cant positive expression correlation with ACE2 in females (r > 0.5) but a negative expression correlation with ACE2 in males (r < 0) (Supplementary Table S4). Fig. 2 shows 25 genes having the strongest positive and negative expression correlation with ACE2 in pan-tissue, female pan-tissue, and male pan-tissue (17). Notably, the gene encoding SLC6A19 (solute carrier family 6 member 19), which interacts with ACE2 (5), had a signi cant positive expression correlation with ACE2 in pan-tissue, female pan-tissue, and male pan-tissue (|r| > 0.3).

Associations of ACE2 expression and immune signatures with the expression of sex hormone receptor genes
Females have a lower disease severity and mortality risk than males infected with SARS-CoV-2 (18,19). A potential explanation is a different host immune response to SARS-CoV-2 infection between females and males (8,14,20,21). We analyzed the correlations between the expression levels of estrogen and androgen receptor genes (ESR1, ESR2, and AR) and ACE2 expression levels in pan-tissue and found that the correlations were consistently positive (Pearson's correlation test, FDR < 1.0 × 10 -60 ) (Fig. 3A).
Interestingly, we found that the expression levels of ESR1 and ESR2 were positively associated with the enrichment levels of immune cells (B cells, CD8+ T cells, and NK cells) (Pearson's correlation test, FDR < 1.0 × 10 -10 ) (Fig. 3B). In contrast, the expression levels of AR inversely correlated with the enrichment levels of these immune cells (FDR < 1.0 × 10 -40 ) (Fig. 3B). We obtained similar results in many individual tissues, including the blood vessel, breast, cervix uteri, colon, esophagus, pituitary, skin, small intestine, stomach, testis, and thyroid (Fig. 3C). These results suggest that females are likely to have a more robust immune defense system against SARS-COV-2 than males.

Discussion
By analyzing the gene expression pro les in human tissues, we identi ed pathways and GO signi cantly associated with ACE2 expression. These pathways and GO were mainly involved in immune response, stromal signature, metabolism, cell growth and proliferation, and cancer and other diseases. In particular, the immune response was positively associated with ACE2 expression, suggesting that the elevated ACE2 expression may boost the immune response. Indeed, the host adaptive and innate immune responses are crucial in ghting off invading SARS-CoV-2, which uses ACE2 as a host cell receptor (3). We found that the enrichment levels of immune cells were positively associated with the expression levels of estrogen receptor genes (ESR1 and ESR2) while they were negatively associated with the expression levels of androgen receptor gene (AR) in human tissues. Meanwhile, these sex hormone receptor genes displayed consistent positive expression correlations with ACE2. Our recent study showed that ACE2 expression levels have no signi cant difference between females and males (8). Collectively, these results suggest that females have a more robust immune defense system against SARS-CoV-2 than males, partially explaining why females have better clinical outcomes of SARS-CoV-2 infections than males (Fig. 4). It also indicates that the supplement with estrogen is a potentially viable approach for treating the patients infected with SARS-CoV-2 (22).
We have also identi ed gene co-expression networks of ACE2. We found that the number of genes having a signi cant positive expression correlation with ACE2 in females far exceeded that in males. Again, this difference may provide cues explaining why females and males have markedly distinct severity and mortality risk of SARS-CoV-2 infection. Thus, the genes with a high expression correlation with ACE2 identi ed in this study warrant further investigation to understand the SARS-CoV-2 infection mechanism.

Declarations
Ethics approval and consent to participate Ethical approval was waived since we used only publicly available data and materials in this study.

Consent for publication
Not applicable.

Availability of data and materials
The GTEx gene expression pro ling datasets for human tissues were downloaded from the GTEx project (https://www.gtexportal.org/home/datasets/).

Competing interests
The authors declare that they have no competing interests.

Funding
This work was supported by the China Pharmaceutical University (grant number 3150120001 to XW).

Authors' Contributions
QF performed data analyses and helped prepare for the manuscript. LL performed data analyses and helped prepare for the manuscript. XW conceived of the research, designed the methods, and wrote the manuscript. All authors read and approved the nal manuscript.    Table S4. 77 genes showing a signi cant positive expression correlation with ACE2 in females (r > 0.5) but a negative expression correlation with ACE2 in males (r < 0). Table S5. The marker genes of immune cells. Figure 1 Pathways and gene ontology associated with ACE2 expression. (A) The KEGG pathways signi cantly associated with ACE2 upregulation identi ed by GSEA (16). (B) Gene ontology signi cantly associated with ACE2 upregulation and downregulation identi ed by WGCNA (25). tissues. The correlation coe cient (r) and the P-value of Pearson's correlation test were shown. * P < 0.05, ** P < 0.01, and *** P < 0.001.

Figure 4
A potential mechanism underlying the signi cantly different clinical outcomes of SARS-CoV-2 infections between females and males.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.