Front. Genet. Frontiers in Genetics Front. Genet. 1664-8021 Frontiers Media S.A. 875762 10.3389/fgene.2022.875762 Genetics Original Research Single-Cell Transcriptomics-Based Study of Transcriptional Regulatory Features in the Non-Obstructive Azoospermia Testis Tang et al. Specific TFs in iNOA Tang Xiao-juan 1 2 Xiao Qiao-hong 2 Wang Xue-lin 2 He Yan 2 3 Tian Ya-nan 2 3 Xia Bin-tong 4 Guo Yang 2 Huang Jiao-long 2 Duan Peng 2 * Tan Yan 1 5 6 * 1 Department of Andrology, Renmin Hospital, Hubei University of Medicine, Shiyan, China 2 Key Laboratory of Zebrafish Modeling and Drug Screening for Human Diseases of Xiangyang City, Department of Obstetrics and Gynaecology, Xiangyang No. 1 People’s Hospital, Hubei University of Medicine, Xiangyang, China 3 Postgraduate Training Basement of Jinzhou Medicical University, Xiangyang No.1 People’s Hospital, Hubei University of Medicine, Xiangyang, China 4 Department of Urology Surgery, Xiangyang No.1 People’s Hospital, Hubei University of Medicine, Xiangyang, China 5 Hubei Key Laboratory of Embryonic Stem Cell Research, Hubei University of Medicine, Shiyan, China 6 Biomedical Engineering College, Hubei University of Medicine, Shiyan, China

Edited by: Xiao Chang, Children’s Hospital of Philadelphia, United States

Reviewed by: Qi-En Yang, Northwest Institute of Plateau Biology (CAS), China

Zhao Liangyu, The Fifth Affiliated Hospital of Sun Yat-sen University, China

 *Correspondence: Peng Duan, meduanpeng@163.com; Yan Tan, tanyan-1@163.com

This article was submitted to Computational Genomics, a section of the journal Frontiers in Genetics

 20 05 2022 2022 13 875762 14 02 2022 04 05 2022 Copyright © 2022 Tang, Xiao, Wang, He, Tian, Xia, Guo, Huang, Duan and Tan. 2022 Tang, Xiao, Wang, He, Tian, Xia, Guo, Huang, Duan and Tan

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

Non-obstructive azoospermia (NOA) is one of the most important causes of male infertility. Although many congenital factors have been identified, the aetiology in the majority of idiopathic NOA (iNOA) cases remains unknown. Herein, using single-cell RNA-Seq data sets (GSE149512) from the Gene Expression Omnibus (GEO) database, we constructed transcriptional regulatory networks (TRNs) to explain the mutual regulatory relationship and the causal relationship between transcription factors (TFs). We defined 10 testicular cell types by their marker genes and found that the proportion of Leydig cells (LCs) and macrophages (tMΦ) was significantly increased in iNOA testis. We identified specific TFs including LHX9, KLF8, KLF4, ARID5B and RXRG in iNOA LCs. In addition, we found specific TFs in iNOA tMΦ such as POU2F2, SPIB IRF5, CEBPA, ELK4 and KLF6. All these identified TFs are strongly engaged in cellular fate, function and homeostasis of the microenvironment. Changes in the activity of the above-mentioned TFs might affect the function of LCs and tMΦ and ultimately cause spermatogenesis failure. This study illustrate that these TFs play important regulatory roles in the occurrence and development of NOA.

 non-obstructive azoospermia spermatogenesis leydig cells testicular macrophages transcription factors National Natural Science Foundation of China10.13039/501100001809 China Scholarship Council10.13039/501100004543 Natural Science Foundation of Hubei Province10.13039/501100003819 Introduction

Male infertility accounts for almost half of all infertility cases and is considered a common, multifactorial pathological condition resulting from a combination of genetic, environmental and lifestyle factors (Kuroda et al., 2020; Sharma et al., 2021). The genetic landscape of male infertility has not been well defined likely due to technical challenges. A substantial proportion of male infertility is accompanied by azoospermia, usually manifested as non-obstructive azoospermia (NOA), which affects about 1% of men in the general population (Tournaye et al., 2017; Zhao et al., 2020; Salas-Huetos et al., 2021). NOA is defined as the complete absence of spermatozoa during ejaculation as a result of failed spermatogenesis. Currently, NOA remains the most challenging and clinically severe form of male infertility and it is primarily associated with genetic abnormalities (Kuroda et al., 2020; Han et al., 2021). A relatively small proportion of NOA is caused by congenital factors such as Klinefelter syndrome (KS) and microdeletions in the azoospermia factor (AZF) region of the Y chromosome (Kuroda et al., 2020; Sharma et al., 2021). However, a majority of the remaining NOA cases have idiopathic (unknown) causes and are diagnosed as idiopathic NOA (iNOA), which accounts for more than 70% of cases (Dabaja and Schlegel, 2013; Zhao et al., 2020; Salas-Huetos et al., 2021). Although there have been some advances in understanding the aetiology and pathogenesis of NOA (including inherent testicular injury or gonadotropin deficiency), the transcriptional regulatory network (TTRN) in spermatogenesis as well as NOA remains not fully clear.

Spermatogenesis is regulated by cross-talk between somatic cells and germ cells in the testis (Meinhardt et al., 2018; Winge et al., 2018; Zhou et al., 2019; Zhao et al., 2020; Hashimoto et al., 2021). Testicular somatic cells, including Leydig cells (LCs) and macrophages, interact with each other to create a supportive microenvironment for germ cell development and self-renewal of spermatogonial stem cells (SSCs), both of which are indispensable for spermatogenesis and male fertility (DeFalco et al., 2015; Teerds and Huhtaniemi, 2015; Mossadegh-Keller et al., 2017; Jauregui et al., 2018; Meinhardt et al., 2018; Winge et al., 2018; Zhou et al., 2019; Lokka et al., 2020; Figueiredo et al., 2021; Hashimoto et al., 2021). In particular, LCs and testicular macrophages (tMΦ) are both located in the testicular interstitial compartment and are functionally related (DeFalco et al., 2015; Mossadegh-Keller et al., 2017; Lokka et al., 2020; Hashimoto et al., 2021). LCs are the primary cells responsible for synthesising and releasing androgens, and these hormones regulate both spermatogenesis and the development of male-specific secondary sex characteristics. tMΦ not only sustain an immune-privileged microenvironment but also engage in collaborative interactions with LCs (Meinhardt et al., 2018; Mossadegh-Keller and Sieweke, 2018; Figueiredo et al., 2021). Increasing evidence indicates that alteration in somatic cell function or the somatic microenvironment could hinder spermatogenesis and lead to NOA (Zhao et al., 2020; Zheng et al., 2021b; Hauptman et al., 2021; Zhu et al., 2021). Nevertheless, the mechanism by which somatic cells contribute to spermatogenesis, particularly the somatic cells in the testicular iNOA microenvironment, is still poorly understood.

More recently, it has emerged that, even after cell development and differentiation, the maintenance of adult somatic cell identity and function relies on the continuous activity of transcription factors (TFs) (Shan et al., 2017). Some TFs could affect not only cell function through expression, but could also regulate cell differentiation (Ieda et al., 2010; Riddell et al., 2014; Han et al., 2018). An increasing number of TFs in testicular somatic cells have been found to regulate a variety of fundamental cell functions during spermatogenesis and testis development (Meroni et al., 2019; Zhou et al., 2019; Jia et al., 2020; Uchida et al., 2020; Wang et al., 2020; Sarkar et al., 2021). However, it remains largely unclear how TFs in LCs and tMΦ play an active role in the development of NOA. In this study, we used the complete atlas of human testicular single cell data to construct transcriptional regulatory networks (TRNs) in the iNOA testis (Zhao et al., 2020). For this endeavour, we used single-cell transcriptome data in conjunction with the gene regulatory network approach. We initially defined TFs influencing the testicular somatic TRNs, revealing that LCs and tMΦ have specific TRNs in the NOA testis.

 Materials and Methods

As our study was based on a conjoint analysis of existing data and no additional patients were included. Hence, ethical approval was not required.

Datasets

We downloaded a human testicular single cell RNA-Seq data set (GSE149512) from the Gene Expression Omnibus (GEO). We selected three normal adults (GSM4504189, GSM4504187 and GSM4504184) and three patients with iNOA (GSM4504195, GSM4504196 and GSM4504197). The quality control (QC) of single-cell RNA-Seq data was performed by using the scater package in R (McCarthy et al., 2017). Genes expressed in at least 2 cells were retained. Mitochondrial (MT) genes were set as the internal reference. For cells with total counts <25,000 or total genes >6,000, the percentage of MT genes >40 were removed. The scImpute package in R was used for imputation, and normalisation was conducted by using the scran package in R (Vieth et al., 2019). RNA-Seq data were normalised by using the transcripts per kilobase million (TPM) method for further analysis.

 Dimensional Reduction and Clustering

We performed principal component analysis (PCA) together with JackStraw and PCEIbow-Plot functions by using the Seurat package (version 3.2.2) in R (version 4.0.2), to select important principal components (PCs) (Lin et al., 2017; Butler et al., 2018). Seurat’s Find All-Markers function was used to identify specific genes for each cell subpopulation. The Run TSNE function was then used for cell clustering and visual analysis of t-distributed stochastic neighbour embedding (t-SNE). The marker genes were thereafter annotated with the singleR package and corrected with CellMarker according to their characteristics (Aran et al., 2019; Zhang et al., 2019). Then, heatmaps were made up of the first 50 marker genes in each cell population, and Gene Ontology (GO) terms were selected to represent the function of each cell type with p < 0.05 among top 30 terms.

 Gene Set Variation Analysis

We performed GSVA to reveal the underlying changes in signalling mechanisms using R (Hanzelmann et al., 2013). The differentially expressed Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were identified from the data of the normal and iNOA groups. The gene set c2. cp.kegg.v7.2. symbols.gmt was downloaded from the Molecular Signatures Database (MSigDB) and set as the reference gene list (Liberzon et al., 2015). After inputting the gene expression profile matrix, the GSVA algorithm transformed the genes of the matrix into scores that represented the activity of each KEGG pathway based on the reference gene. Then, the differentially activated pathways between the normal and iNOA groups were determined by using the limma package in R with |log2fold change| ≥ 2 and p < 0.05 (Ritchie et al., 2015).

 Inference of Regulons and Activity

A number of methods have been developed to predict genetic regulatory networks (GRNs) from single-cell gene expression data. We adopted the SCENIC method as previously described with slight modification (Fiers et al., 2018). The SCENIC analytic process, comprised three steps. First, a gene co-expression network was constructed through gene co-expression analysis. Second, possible TF-target regulatory relationships were established based on the gene co-expression network. In this step, the direct regulatory relationship was established by using motif analysis. Any direct downstream genes occurring for each TF were profiled as regulons. In particular, SCENIC could only support transcriptional positive regulation analysis. Third, based on the results of step 2, a regulon activity score (RAS) was calculated for each cell. As described in previous studies, the Avg20 method was repeated three times to assess the variability of random sampling. Thereafter, a t-test was used to assess whether the Avg20 method performed better than using all individual cells (Cao et al., 2017).

 Functional Validation

As in a previous study, we used Search-Based Exploration of Expression (SEEK) analysis to determine whether the predicted regulons correlated with their cell type (Zhu et al., 2015). In brief, we used the human version of SEEK to assess whether genes in the regulons were co-expressed. Significantly co-expressed genes in multiple data sets associated with a particular cell type scored positive for high relevance of the function of the regulon to that cell type.

 Regulon Module Analysis and Quantifying Cell Type Relationship

To identify regulon modules, we employed two main steps (Aibar et al., 2017). First, each pair of regulatory relationships was analysed by calculating the Pearson correlation coefficient. To systematically describe the regulatory relationships of TFs, we compared the regulatory activity scores of each regulatory pair based on the connection-specificity index (J. I. F. Bass et al., 2013). For each pair of regulatory relationships, we defined a regulatory specificity score (RSS) based on the Jensen–Shannon scatterplot (Cabili et al., 2011). Next, we selected the specific regulator with the highest RSS value and further examined its functional characteristics. The activity score of each regulon module in relation to a cell type was then defined as the average of the activity scores of its regulon members in all cells of that cell type. The highest ranked units were then filtered for each module. We quantified the relationship between different cell types based on the similarity of overall regulon activity. A pair of cell types was linked if their Spearman correlation coefficient was >0.8. Finally, we used the Markov Clustering Algorithm (MCL) to identify related cell types (Van Dongen and Abreu-Goodger, 2012).

 GO, KEGG Enrichment Analysis and Protein-Protein Interaction Networks

We conducted GO enrichment analysis and KEGG analysis on the TF of each regulatory module example. GO analysis depicts the unique biological significance based on differentially expressed genes (DEGs) between groups. We used the KEGG database to determine important pathways. The ‘p < 0.05’ and the ‘|log2 fold change| ≥ 2’ conditions were used as the cut-off criteria for GO and KEGG enrichment analyses. The genes of each module are then incorporated into a search tool (STRING) (Jensen et al., 2009) that retrieves interacting genes/proteins. In the multivariate analysis, the confidence levels were 0.4, 0.7 and 0.9. Then, we input the gene network file into Cytoscape to draw PPIN diagrams.

 Results Single-Cell Maps Define the Heterogeneity of Normal and iNOA Testicular Cells

To investigate the heterogeneity of patients with iNOA, we used the Seurat package to perform quality control and t-SNE analysis on single-cell data from the GEO (GSE149512) data set. 32,048 cells, each with 500–8,000 genes, were reserved for subsequent analysis. We divided the normal and iNOA groups into 10 cell populations based on marker genes of each cell population (Figures 1A,B, Supplementary Figure S1). The cellular composition of the iNOA group was different compared qwith the normal group. Almost no sperm cells were observed in the iNOA group, while the proportion of LCs and tMΦ were increased significantly in the iNOA group compared with the normal group (Figure 1C).

Single-cell RNA sequencing (scRNA-Seq) analysis shows the testicular cell lineages in normal individuals and patients with idiopathic non-obstructive azoospermia (iNOA). (A) t-Distributed random neighbour embedding (t-SNE) isolation (left) and combination (middle and right) of single testicular cells in the normal and iNOA groups. Ten main cell types were defined, including endotheliocytes, T cells, mast cells, Leydig cells, Sertoli cells, peritubular myoid (PTM) cells, macrophages, spermatogonia, spermatocytes and spermatids. (B) The t-SNE map shows the expression level distribution of marker genes in cell types, including FGFR3 (spermatogonia), SYCP3 (spermatocytes), TNP1 (spermatids), VWF (endotheliocytes), CD163 (macrophages), DLK (Leydig cells), MYH11 (PTM cells), SOX9 (Sertoli cells), CCL5 (T cells) and TPSB2 (mast cells). (C) The mean cell number and relative proportion of testicular subsets from various sample sources.

 Marker Gene GO Analysis

We identified unique characteristics of the iNOA group. For example, normal LC GO terms include ‘Regulation of multicellular organismal development’ and ‘Negative Regulation of multicellular cells’ organismal process (Figure 2A), while iNOA group LCs GO terms include ‘Extracellular matrix structural constituent’ and ‘Glycosaminoglycan binding’ (Figure 2B). The normal tMΦ GO terms include ‘Myeloid leukocyte activation’, while iNOA tMΦ GO terms include ‘Immune effector process’ (Figures 2A,B). These results indicate that the functions of LCs and tMΦ are different between the normal and iNOA groups.

Gene Ontology (GO) enrichment analysis of marker genes in (A) the normal group and (B) the idiopathic non-obstructive azoospermia (iNOA) group. Left: heatmap showing expression signatures of the top 50 specifically expressed genes in each cell type; the value for each gene is the row-scaled Z score. Right: representative GO terms.

 GSVA Analysis of Testicular Cells

A direct comparison of the iNOA and normal groups revealed that ‘G2M checkpoint’, ‘MTORC1 signaling’ and ‘TGF beta signaling’ as the top enriched signatures in iNOA LCs, are associated with proliferation (Figure 3A). Moreover, the GSVA scores of the ‘mitotic spindle’, ‘MTORC1 signaling’ and ‘TGF beta signaling’ were obviously increased in iNOA tMΦ, and these pathways are also associated with proliferation (Figure 3B). We found that the spermatogenesis pathways were obviously downregulated in the iNOA group including LCs, tMΦ, peritubular myoid (PTM) cells and Sertoli cells compared with the normal group (Figures 3A,B, Supplementary Figure S2).

Gene set variation analysis (GSVA) of testicular somatic cells in the normal group and the idiopathic non-obstructive azoospermia (iNOA) group. (A) Differences in pathway activities scored per cell by GSVA between normal and iNOA Leydig cells (n = 968 and 9,412 cells, respectively; six patients per group). (B) Differences in pathway activities scored per cell by GSVA between normal and iNOA macrophages (n = 144 and 763 cells, respectively; six patients per group).

 Analysis of Specific Regulation of Somatic Cell Types

Our network analysis found that NR2F1, CREB3l1, HC1, GLI2 and KLF4 were specific TFs associated with LCs in normal human (Supplementary Figure S3A). The t-SNE diagrams further demonstrated the highly specific activity of NR2F1 in LCs (Supplementary Figure S3B, C). To test the validity of the above analysis, SEEK analysis was used to identify GEO data sets with significant co-expression of NR2F1. SEEK analysis did not show significant co-expression of the regulatory gene NR2F1 (Supplementary Figure S3D, Fisher’s exact test, p = 0.0858). CEBPA, MEF2A, IRF5, HIVEP2 and EOMES were found to be specific TFs in normal tM φ (Supplementary Figure S3E). CEBPA activity was highly specific in normal tMΦ (Supplementary Figure S3F, G). SEEK analysis did not show significant co-expression of the regulatory gene CEBPA (Supplementary Figure S3H, Fisher’s exact test, p = 1). Other specific regulation of cell types is shown in Supplementary Figure S3I–P; Furthermore, our network analysis identified LHX9, ARID5B, KLF8, RXRG and KLF4 as specific TFs associated with iNOA LCs, and found that the activity of these TFs was elevated in iNOA patients (Figure 4A). The t-SNE diagram further demonstrated the highly specific activity of LHX9 in iNOA LCs (Figures 4B,C). To test the validity of the above analysis, we performed SEEK analysis to identify GEO data sets with significant co-expression of LHX9 (Figure 4D, Fisher’s exact test, p = 0.0345). We identified POU2F2, SPIB, IRF5, CEBPA and CREM as specific TFs in iNOA tMΦ, and the activity of these TFs was also elevated in iNOA patients (Figure 4E). The activity of POU2F2 in iNOA tMΦ was highly specific (Figures 4F,G). SEEK analysis did not show significant co-expression of the regulatory gene POU2F2 (Figure 4H, Fisher’s exact test, p = 0.0608). For other specific regulation of cell types, see Figure 4I–P; for specific regulation of cell types in the normal group, see Supplementary Figure S3.

Analysis of cell type-specific regulation in testis of idiopathic non-obstructive azoospermia (iNOA) patients. (A–D) Leydig cells (LCs). (A) The regulation of testicular cells was ranked according to the regulation specificity score. (B) LCs are marked by red dots on the t-distributed random neighbour embedding (t-SNE) diagram. (C) The expression values of the genes with the highest regulatory activity score are displayed in the t-SNE diagram. (D) Search-Based Exploration of Expression (SEEK) analysis was used to determine the co-expression of the highest regulatory genes in different Gene Expression Omnibus (GEO) data sets. The x-axis represents the different data sets, and the y-axis represents the co-expression significance of the target gene in each data set. Data sets with a significant correlation (p < 0.05) are highlighted with yellow dots. (E–H) The same as (A–D), but for testicular macrophages. (I–L) The same as (A–D), but for peritubular myoid (PTM) cells. (M–P) The same as (A–D), but for Sertoli cells.

 Organizing Regulons Into Combinatorial Modules

Based on the identified connection specificity index (CSI) matrix model (M1–M7), we mapped the average activity of each module to t-SNE (Figure 5A). The iNOA testicular cells were then ranked according to the regulatory specificity score (Figure 5B). Each module occupies a different region, and all the highlighted regions indicate the high transcriptional activity of different modules (Figure 5A). As shown in Figure 6A, the M1 and M5 modules showed high transcriptional activity mainly in iNOA tMΦ, while the M2 and M4 modules showed high transcriptional activity mainly in iNOA LCs. Figure 6A showed the determination of the regulation module based on the regulation CSI matrix, along with associated cell types, corresponding binding motifs, and representative transcription factors. The PPINs of TFs in each module were shown in Figure 6B. The M1 module contains MYB, POU2F2 and PBX4 in iNOA tMΦ. In addition, in iNOA LCs the M2 and M4 modules contain LHX9, KLF8, KLF4 and FOXC2. We performed GO functional enrichment analysis on genes of the M1, M2, M4 and M5 modules, as shown in Supplementary Figure S4.

The activity of regulatory modules in testicular somatic cells of different idiopathic non-obstructive azoospermia (iNOA) types. (A) Identification of the regulatory modules (M1–M7) according to the regulatory connection specificity index (CSI) matrix and a map of the average activity of each module based on t-distributed random neighbour embedding (t-SNE). (B) Ranking regulation in iNOA testicular cells based on the regulation specificity score. The y-axis represents the regulation activity score. The x-axis represents the cell type.

Identification of combinatorial regulon modules. (A) Determination of the regulon modules based on the regulation connection specificity index (CSI) matrix, along with associated cell types, corresponding binding motifs, and representative transcription factors (TFs). (B) Protein-protein interaction networks (PPINs) of regulatory factors in each module.

 Discussion

Spermatogenesis is a highly sophisticated process with well-organised cellular and molecular events involving gene networks regulated by TFs in testicular somatic cells (Bettegowda and Wilkinson, 2010). Our results also support the view that the somatic cellular microenvironment provides favourable biochemical and biophysical components for spermatogenesis. We focussed on LCs and tMΦ in the testicular microenvironment, which are critical for sperm production, and identified some important TFs in LCs and tMΦ that could play an important role in NOA genesis.

Unlike obstructive azoospermia, the pathogenesis of iNOA remains complex and multifactorial. To study the molecular mechanism of iNOA pathogenesis, we analysed single-cell data sets from the GEO database, and clarified the heterogeneity of different cell types of the human testis in more detail by using single-cell transcriptome sequencing analysis. We found that the proportion of LCs and tMΦ in iNOA testis is higher than the normal group. Furthermore, the function of LCs and tMΦ in the iNOA testis is different from the normal group. We suggest that changes in the proportions and functions of LCs and tMΦ could alter the testicular microenvironment, contributing to spermatogenesis failure and thereby leading to male infertility.

Testicular LCs are the major producers of circulating testosterone, which is essential for testis development and spermatogenesis. LC dysfunction can lead to testosterone deficiency and impair male fertility (Zhou et al., 2019). The function of LCs is correlated with the developmental stages of LC lineage specification and differentiation, both as Sertoli cells and tMΦ. LC development involves at least three steps: the proliferation of LC precursors (also known as stem LCs), their differentiation into immature LCs and their final maturation into adult LCs (Peak et al., 2016). We suppose that LCs of patients with iNOA remain at the stage of proliferation without transitioning to the differentiation and maturation stages. These immature and non-functional LCs are unable to (fully) fulfil their steroidogenic function to maintain spermatogenesis (O'Hara et al., 2015; Guan et al., 2019). The successive stages of LC development are regulated by an array of lineage-specific TFs (Winge et al., 2018). We identified specific TFs including LHX9, KLF8, KLF4, ARID5B and RXRG in LCs of patients with iNOA. Among these TFs, LHX9 is an important steroidogenesis-related TF and indispensable for testis development (Hu et al., 2018). KLF4 and KLF8 are Krüppel-like factors known to regulate several biological processes, such as cell proliferation, differentiation and metabolism (Chu et al., 2016; Kult et al., 2021; Kumar et al., 2021). Similarly to KLF4, ARID5B also plays a pivotal role in adipogenesis and lipid metabolism, which might be closely related to cell differentiation and development (Hata et al., 2013; Claussnitzer et al., 2015). RXRG is a rexinoid receptor that participates in the regulation of cell differentiation (Gely-Pernot et al., 2015; Cheng et al., 2018). All these identified TFs play direct/indirect roles in cell development. However, the occurrence of NOA caused by TFs has not been reported. This is our new finding. We conducted GO analysis on LCs of normal persons and iNOA patients, and found that LCs of iNOA patients were inhibited in the ‘Regulation of Multicellular Development’ pathway. GSVA enrichment analysis showed that LCs of iNOA patients were significantly enriched in the proliferation pathway, indicating that the development of LCs in iNOA patients was inhibited in the proliferation stage, but did not enter the stages of differentiation and maturity. We infer that changes in the activity of these TFs potentially affect the maturation and function of LCs and impair the microenvironment of spermatogenesis, dysfunctions that may eventually cause azoospermia in humans.

TMΦ are the principal immune cell population of the mammalian testis, and together with LCs and Sertoli cells, they maintain testicular immune privilege (Heinrich and DeFalco, 2020; Rehman et al., 2021). tMΦ produce several growth and differentiation factors for LC development (Heinrich and DeFalco, 2020). In addition to testicular immunosuppression, tMΦ also locally regulate LC steroidogenesis (Fang et al., 2021). Previous studies have delineated the phenotype of tMΦ in the normal testis, and tMΦ have been implicated in the development of azoospermia (Hussein et al., 2005; Duan et al., 2011). We identified specific TFs including POU2F2, SPIB, IRF5, CEBPA, ELK4 and KLF6 in iNOA tMΦ. Among these identified TFs, both POU2F2 and SPIB are essential for cell proliferation, differentiation and functional maturation of immune cells (Hodson et al., 2016; Klisuric et al., 2019). IRF5 has been shown to act as a master switch that promotes proinflammatory cytokine production from macrophages and thus contributes to the plasticity of macrophage polarisation (Banga et al., 2020). CEBPA is required for the regulation of cell proliferation and terminal differentiation and participates in the control of immune and inflammatory processes (Bristol et al., 2009; Minner et al., 2019). ELK4 has been implicated in maintaining cellular homeostasis, but also in macrophage M2 polarisation (Zheng et al., 2021a). Similarly to KLF4 and KLF8, KLF6 signalling engages in various cellular processes, including cell differentiation and development. Moreover, KLF6 impedes macrophage polarisation to the M2 phenotype (Zhao et al., 2021). These TFs in macrophages are newly discovered and may affect NOA. GSVA enrichment analysis of tMΦ in normal persons and iNOA patients showed that tMΦ of iNOA patients was significantly enriched in the proliferation pathway, which also indicated that thetMΦ of iNOA patients was overproliferated, which was consistent with the previous research results of Wenzhong Zheng (Zheng et al., 2021). The changes of TFs activity found in iNOA tMΦ may regulate the development and function of tMΦ and the differentiation of tM into M2 macrophages. When macrophages are activated and elaborate inflammatory mediators, the function of LCs secretion of testosterone may be impaired. (Hales, 2002), and thus fail to create an optimal immune microenvironment for spermatogenesis.

In this study, the transcriptional regulatory network of iNOA testicular cells was established to provide a new idea for understanding the regulatory mechanism and functional relationship of iNOA testicular cells. However, the roles and mechanisms of these TFs in iNOA pathogenesis in different testicular cell types need to be further investigated experimentally.

 Conclusion

With the recent application of single-cell sequencing technology in the human testis, the understanding of testicular cell heterogeneity has greatly improved. It is clear that we need more in-depth investigation of the mechanism of cellular heterogeneity during iNOA development. This study has provided a new approach to dissect the regulatory mechanisms and functional relationships by establishing the TRNs and PPINs of TFs in iNOA LCs and tMΦ. Several of the identified TFs, such as Krüppel-like factors, are predicted to regulate the differentiation and function of both LCs and tMΦ. We demonstrated that aberrant regulation of TFs identified in iNOA LCs and tMΦ potentially affects the testicular microenvironment and germ cell development. This study should improve knowledge regarding TFs involved in the regulatory landscape of LC and tMΦ development and the crosstalk among cell types. However, additional experiments are needed to investigate the function and mechanism of TFs in different testicular cell types that are involved the pathogenesis of NOA.

 Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.

 Author Contributions

PD, YT and X-JT launched the study. X-JT performed the data analysis. Q-HX, X-LW, YH, Y-NT, B-TX, YG, and J-LH participated in reference collecting. X-JT, Q-HX, and PD completed the manuscript. All authors contributed to the article and approved the submitted version.

 Funding

This study was supported by the research grants from the National Natural Science Foundation of China (81901567), the China Scholarship Council (201808420351), the Natural Science Foundation of Hubei Provincial Department of Education (Q20202105), the Scientific and Technological Project of Xiangyang City of Hubei Province (2021YL29 and 2021YL30), and the Innovative Research Program for Graduates of Hubei University of Medicine (YC2020001).

 Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

 Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

 Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2022.875762/full#supplementary-material

References Aibar S. González-Blas C. B. Moerman T. Huynh-Thu V. A. Imrichova H. Hulselmans G. (2017). SCENIC: Single-Cell Regulatory Network Inference and Clustering. Nat. Methods 14 (11), 1083–1086. 10.1038/nmeth.4463 Aran D. Looney A. P. Liu L. Wu E. Fong V. Hsu A. (2019). Reference-based Analysis of Lung Single-Cell Sequencing Reveals a Transitional Profibrotic Macrophage. Nat. Immunol. 20 (2), 163–172. 10.1038/s41590-018-0276-y Banga J. Srinivasan D. Sun C.-C. Thompson C. D. Milletti F. Huang K.-S. (2020). Inhibition of IRF5 Cellular Activity with Cell-Penetrating Peptides that Target Homodimerization. Sci. Adv. 6 (20), eaay1057. 10.1126/sciadv.aay1057 Bass J. I. F. Diallo A. Nelson J. Soto J. M. Myers C. L. Walhout A. J. M. (2013). Using Networks to Measure Similarity between Genes: Association Index Selection. Nat. Methods 10 (12), 1169–1176. 10.1038/nmeth.2728 Bettegowda A. Wilkinson M. F. (2010). Transcription and Post-transcriptional Regulation of Spermatogenesis. Phil. Trans. R. Soc. B 365 (1546), 1637–1651. 10.1098/rstb.2009.0196 Bristol J. A. Morrison T. E. Kenney S. C. (2009). CCAAT/enhancer Binding Proteins α and β Regulate the Tumor Necrosis Factor Receptor 1 Gene Promoter. Mol. Immunol. 46 (13), 2706–2713. 10.1016/j.molimm.2009.05.024 Butler A. Hoffman P. Smibert P. Papalexi E. Satija R. (2018). Integrating Single-Cell Transcriptomic Data across Different Conditions, Technologies, and Species. Nat. Biotechnol. 36 (5), 411–420. 10.1038/nbt.4096 Cabili M. N. Trapnell C. Goff L. Koziol M. Tazon-Vega B. Regev A. (2011). Integrative Annotation of Human Large Intergenic Noncoding RNAs Reveals Global Properties and Specific Subclasses. Genes Dev. 25 (18), 1915–1927. 10.1101/gad.17446611 Cao J. Packer J. S. Ramani V. Cusanovich D. A. Huynh C. Daza R. (2017). Comprehensive Single-Cell Transcriptional Profiling of a Multicellular Organism. Science 357 (6352), 661–667. 10.1126/science.aam8940 Cheng S. Wang M. Wang Y. Zhang C. Wang Y. Song J. (2018). RXRG Associated in PPAR Signal Regulated the Differentiation of Primordial Germ Cell. J. Cell. Biochem. 119 (8), 6926–6934. 10.1002/jcb.26891 Chu L.-F. Leng N. Zhang J. Hou Z. Mamott D. Vereide D. T. (2016). Single-cell RNA-Seq Reveals Novel Regulators of Human Embryonic Stem Cell Differentiation to Definitive Endoderm. Genome Biol. 17 (1), 173. 10.1186/s13059-016-1033-x Claussnitzer M. Dankel S. N. Kim K.-H. Quon G. Meuleman W. Haugen C. (2015). FTO Obesity Variant Circuitry and Adipocyte Browning in Humans. N. Engl. J. Med. 373 (10), 895–907. 10.1056/NEJMoa1502214 Dabaja A. A. Schlegel P. N. (2013). Microdissection Testicular Sperm Extraction: an Update. Asian J. Androl. 15 (1), 35–39. 10.1038/aja.2012.141 DeFalco T. Potter S. J. Williams A. V. Waller B. Kan M. J. Capel B. (2015). Macrophages Contribute to the Spermatogonial Niche in the Adult Testis. Cell Rep. 12 (7), 1107–1119. 10.1016/j.celrep.2015.07.015 Duan Y.-G. Yu C.-F. Novak N. Bieber T. Zhu C.-H. Schuppe H.-C. (2011). Immunodeviation towards a Th17 Immune Response Associated with Testicular Damage in Azoospermic Men. Int. J. Androl. 34 (6 Pt 2), e536–e545. 10.1111/j.1365-2605.2010.01137.x Fang Y. Su Y. Xu J. Hu Z. Zhao K. Liu C. (2021). Varicocele-mediated Male Infertility: from the Perspective of Testicular Immunity and Inflammation. Front. Immunol. 12, 729539. 10.3389/fimmu.2021.729539 Fiers M. W. E. J. Minnoye L. Aibar S. Bravo González-Blas C. Kalender Atak Z. Aerts S. (2018). Mapping Gene Regulatory Networks from Single-Cell Omics Data. Brief. Funct. Genomics 17 (4), 246–254. 10.1093/bfgp/elx046 Figueiredo A. F. A. Wnuk N. T. Vieira C. P. Gonçalves M. F. F. Brener M. R. G. Diniz A. B. (2021). Activation of C-C Motif Chemokine Receptor 2 Modulates Testicular Macrophages Number, Steroidogenesis, and Spermatogenesis Progression. Cell Tissue Res. 386 (1), 173–190. 10.1007/s00441-021-03504-w Gely-Pernot A. Raverdeau M. Teletin M. Vernet N. Féret B. Klopfenstein M. (2015). Retinoic Acid Receptors Control Spermatogonia Cell-Fate and Induce Expression of the SALL4A Transcription Factor. PLoS Genet. 11 (10), e1005501. 10.1371/journal.pgen.1005501 Guan X. Chen F. Chen P. Zhao X. Mei H. Liu J. (2019). Effects of Spermatogenic Cycle on Stem Leydig Cell Proliferation and Differentiation. Mol. Cell. Endocrinol. 481, 35–43. 10.1016/j.mce.2018.11.007 Hales D. B. (2002). Testicular Macrophage Modulation of Leydig Cell Steroidogenesis. J. Reproductive Immunol. 57, 3–18. 10.1016/s0165-0378(02)00020-7 Han B. Yan Z. Yu S. Ge W. Li Y. Wang Y. (2021). Infertility Network and Hub Genes for Nonobstructive Azoospermia Utilizing Integrative Analysis. Aging 13 (5), 7052–7066. 10.18632/aging.202559 Han X. Wang R. Zhou Y. Fei L. Sun H. Lai S. (2018). Mapping the Mouse Cell Atlas by Microwell-Seq. Cell 173 (5), 1307. 10.1016/j.cell.2018.05.012 Hänzelmann S. Castelo R. Guinney J. (2013). GSVA: Gene Set Variation Analysis for Microarray and RNA-Seq Data. BMC Bioinforma. 14, 7. 10.1186/1471-2105-14-7 Hashimoto M. Kimura S. Kanno C. Yanagawa Y. Watanabe T. Okabe J. (2021). Macrophage Ubiquitin-specific Protease 2 Contributes to Motility, Hyperactivation, Capacitation, and In Vitro Fertilization Activity of Mouse Sperm. Cell. Mol. Life Sci. 78 (6), 2929–2948. 10.1007/s00018-020-03683-9 Hata K. Takashima R. Amano K. Ono K. Nakanishi M. Yoshida M. (2013). Arid5b Facilitates Chondrogenesis by Recruiting the Histone Demethylase Phf2 to Sox9-Regulated Genes. Nat. Commun. 4, 2850. 10.1038/ncomms3850 Hauptman D. Perić M. H. Marić T. Bojanac A. K. Sinčić N. Zimak Z. (2021). Leydig Cells in Patients with Non-obstructive Azoospermia: Do They Really Proliferate? Life 11 (11), 1266. 10.3390/life11111266 Heinrich A. DeFalco T. (2020). Essential Roles of Interstitial Cells in Testicular Development and Function. Andrology 8 (4), 903–914. 10.1111/andr.12703 Hodson D. J. Shaffer A. L. Xiao W. Wright G. W. Schmitz R. Phelan J. D. (2016). Regulation of Normal B-Cell Differentiation and Malignant B-Cell Survival by OCT2. Proc. Natl. Acad. Sci. U.S.A. 113 (14), E2039–E2046. 10.1073/pnas.1600557113 Hu F. Zhu Q. Sun B. Cui C. Li C. Zhang L. (2018). Smad Ubiquitylation Regulatory Factor 1 Promotes LIM‐homeobox Gene 9 Degradation and Represses Testosterone Production in Leydig Cells. FASEB J. 32 (9), 4627–4640. 10.1096/fj.201701480R Hussein M. Aboudeif E. Bedaiwy M. Said T. Mustafa M. Nada E. (2005). Phenotypic Characterization of the Immune and Mast Cell Infiltrates in the Human Testis Shows Normal and Abnormal Spermatogenesis. Fertil. Steril. 83 (5), 1447–1453. 10.1016/j.fertnstert.2004.11.062 Ieda M. Fu J.-D. Delgado-Olguin P. Vedantham V. Hayashi Y. Bruneau B. G. (2010). Direct Reprogramming of Fibroblasts into Functional Cardiomyocytes by Defined Factors. Cell 142 (3), 375–386. 10.1016/j.cell.2010.07.002 Jauregui E. J. Mitchell D. Topping T. Hogarth C. A. Griswold M. D. (2018). Retinoic Acid Receptor Signaling Is Necessary in Steroidogenic Cells for Normal Spermatogenesis and Epididymal Function. Development 145 (13), dev160465. 10.1242/dev.160465 Jensen L. J. Kuhn M. Stark M. Chaffron S. Creevey C. Muller J. (2009). STRING 8--a Global View on Proteins and Their Functional Interactions in 630 Organisms. Nucleic Acids Res. 37 (17), D412–D416. 10.1093/nar/gkn760 Jia G.-X. Lin Z. Yan R.-G. Wang G.-W. Zhang X.-N. Li C. (2020). WTAP Function in Sertoli Cells Is Essential for Sustaining the Spermatogonial Stem Cell Niche. Stem Cell Rep. 15 (4), 968–982. 10.1016/j.stemcr.2020.09.001 Klisuric A. Thierry B. Delon L. Prestidge C. A. Gibson R. J. (2019). Identifying Human and Murine M Cells In Vitro . Exp. Biol. Med. (Maywood) 244 (7), 554–564. 10.1177/1535370219838674 Kult S. Olender T. Osterwalder M. Markman S. Leshkowitz D. Krief S. (2021). Bi-fated Tendon-To-Bone Attachment Cells Are Regulated by Shared Enhancers and KLF Transcription Factors. Elife 10, e55361. 10.7554/eLife.55361 Kumar S. Behera A. Saha P. Kumar Srivastava A. (2021). The Role of Krüppel-like Factor 8 in Cancer Biology: Current Research and its Clinical Relevance. Biochem. Pharmacol. 183, 114351. 10.1016/j.bcp.2020.114351 Kuroda S. Usui K. Sanjo H. Takeshima T. Kawahara T. Uemura H. (2020). Genetic Disorders and Male Infertility. Reprod. Med. Biol. 19 (4), 314–322. 10.1002/rmb2.12336 Liberzon A. Birger C. Thorvaldsdóttir H. Ghandi M. Mesirov J. P. Tamayo P. (2015). The Molecular Signatures Database Hallmark Gene Set Collection. Cell Syst. 1 (6), 417–425. 10.1016/j.cels.2015.12.004 Lin P. Troup M. Ho J. W. K. (2017). CIDR: Ultrafast and Accurate Clustering through Imputation for Single-Cell RNA-Seq Data. Genome Biol. 18 (1), 59. 10.1186/s13059-017-1188-0 Lokka E. Lintukorpi L. Cisneros-Montalvo S. Mäkelä J.-A. Tyystjärvi S. Ojasalo V. (2020). Generation, Localization and Functions of Macrophages during the Development of Testis. Nat. Commun. 11 (1), 4375. 10.1038/s41467-020-18206-0 McCarthy D. J. Campbell K. R. Lun A. T. L. Wills Q. F. (2017). Scater: Pre-processing, Quality Control, Normalization and Visualization of Single-Cell RNA-Seq Data in R. Bioinformatics 33 (8), btw777–1186. 10.1093/bioinformatics/btw777 Meinhardt A. Wang M. Schulz C. Bhushan S. (2018). Microenvironmental Signals Govern the Cellular Identity of Testicular Macrophages. J. Leukoc. Biol. 104 (4), 757–766. 10.1002/JLB.3MR0318-086RR Meroni S. B. Galardo M. N. Rindone G. Gorga A. Riera M. F. Cigorraga S. B. (2019). Molecular Mechanisms and Signaling Pathways Involved in Sertoli Cell Proliferation. Front. Endocrinol. 10, 224. 10.3389/fendo.2019.00224 Minner S. Lutz J. Hube‐Magg C. Kluth M. Simon R. Höflmayer D. (2019). Loss of CCAAT‐enhancer‐binding Protein Alpha (CEBPA) Is Linked to Poor Prognosis in PTEN Deleted and TMPRSS2:ERG Fusion Type Prostate Cancers. Prostate 79 (3), 302–311. 10.1002/pros.23736 Mossadegh-Keller N. Gentek R. Gimenez G. Bigot S. Mailfert S. Sieweke M. H. (2017). Developmental Origin and Maintenance of Distinct Testicular Macrophage Populations. J. Exp. Med. 214 (10), 2829–2841. 10.1084/jem.20170829 Mossadegh-Keller N. Sieweke M. H. (2018). Testicular Macrophages: Guardians of Fertility. Cell. Immunol. 330, 120–125. 10.1016/j.cellimm.2018.03.009 O'Hara L. McInnes K. Simitsidellis I. Morgan S. Atanassova N. Slowikowska‐Hilczer J. (2015). Autocrine Androgen Action Is Essential for Leydig Cell Maturation and Function, and Protects against Late‐onset Leydig Cell Apoptosis in Both Mice and Men. FASEB J. 29 (3), 894–910. 10.1096/fj.14-255729 Peak T. C. Haney N. M. Wang W. DeLay K. J. Hellstrom W. J. (2016). Stem Cell Therapy for the Treatment of Leydig Cell Dysfunction in Primary Hypogonadism. Wjsc 8 (10), 306–315. 10.4252/wjsc.v8.i10.306 Rehman A. Pacher P. Haskó G. (2021). Role of Macrophages in the Endocrine System. Trends Endocrinol. Metabolism 32 (4), 238–256. 10.1016/j.tem.2020.12.001 Riddell J. Gazit R. Garrison B. S. Guo G. Saadatpour A. Mandal P. K. (2014). Reprogramming Committed Murine Blood Cells to Induced Hematopoietic Stem Cells with Defined Factors. Cell 157 (3), 549–564. 10.1016/j.cell.2014.04.006 Ritchie M. E. Phipson B. Wu D. Hu Y. Law C. W. Shi W. (2015). Limma Powers Differential Expression Analyses for RNA-Sequencing and Microarray Studies. Nucleic Acids Res. 43 (7), e47. 10.1093/nar/gkv007 Salas-Huetos A. Tuttelmann F. Tüttelmann F. Wyrwoll M. J. Kliesch S. Lopes A. M. (2021). Disruption of Human Meiotic Telomere Complex Genes TERB1, TERB2 and MAJIN in Men with Non-obstructive Azoospermia. Hum. Genet. 140 (1), 217–227. 10.1007/s00439-020-02236-1 Sarkar R. K. Sen Sharma S. Mandal K. Wadhwa N. Kunj N. Gupta A. (2021). Homeobox Transcription Factor Meis1 Is Crucial to Sertoli Cell Mediated Regulation of Male Fertility. Andrology 9 (2), 689–699. 10.1111/andr.12941 Shan Q. Zeng Z. Xing S. Li F. Hartwig S. M. Gullicksrud J. A. (2017). The Transcription Factor Runx3 Guards Cytotoxic CD8+ Effector T Cells against Deviation towards Follicular Helper T Cell Lineage. Nat. Immunol. 18 (8), 931–939. 10.1038/ni.3773 Sharma A. Minhas S. Dhillo W. S. Jayasena C. N. (2021). Male Infertility Due to Testicular Disorders. J. Clin. Endocrinol. Metab. 106 (2), e442–e459. 10.1210/clinem/dgaa781 Teerds K. J. Huhtaniemi I. T. (2015). Morphological and Functional Maturation of Leydig Cells: from Rodent Models to Primates. Hum. Reprod. Update 21 (3), 310–328. 10.1093/humupd/dmv008 Tournaye H. Krausz C. Oates R. D. (2017). Novel Concepts in the Aetiology of Male Reproductive Impairment. Lancet Diabetes & Endocrinol. 5 (7), 544–553. 10.1016/s2213-8587(16)30040-7 Uchida A. Sakib S. Labit E. Abbasi S. Scott W. Underhill M. (2020). Development and Function of Smooth Muscle Cells Is Modulated by Hic1 in Mouse Testis. Development 147 (13), dev185884. 10.1242/dev.185884 Van Dongen S. Abreu-Goodger C. (2012). Using MCL to Extract Clusters from Networks. Methods Mol. Biol. 804, 281–295. 10.1007/978-1-61779-361-5_15 Vieth B. Parekh S. Ziegenhain C. Enard W. Hellmann I. (2019). A Systematic Evaluation of Single Cell RNA-Seq Analysis Pipelines. Nat. Commun. 10 (1), 4667. 10.1038/s41467-019-12266-7 Wang S. Wang P. Liang D. Wang Y. (2020). BRG1 Is Dispensable for Sertoli Cell Development and Functions in Mice. Ijms 21 (12), 4358. 10.3390/ijms21124358 Winge S. B. Dalgaard M. D. Belling K. G. Jensen J. M. Nielsen J. E. Aksglaede L. (2018). Transcriptome Analysis of the Adult Human Klinefelter Testis and Cellularity-Matched Controls Reveals Disturbed Differentiation of Sertoli- and Leydig Cells. Cell Death Dis. 9 (6), 586. 10.1038/s41419-018-0671-1 Zhang X. Lan Y. Xu J. Quan F. Zhao E. Deng C. (2019). CellMarker: a Manually Curated Resource of Cell Markers in Human and Mouse. Nucleic Acids Res. 47 (D1), D721–D728. 10.1093/nar/gky900 Zhao L. Yao C. Xing X. Jing T. Li P. Zhu Z. (2020). Single-cell Analysis of Developing and Azoospermia Human Testicles Reveals Central Role of Sertoli Cells. Nat. Commun. 11 (1), 5683. 10.1038/s41467-020-19414-4 Zhao Y. Yu Z. Ma R. Zhang Y. Zhao L. Yan Y. (2021). lncRNA-Xist/miR-101-3p/KLF6/C/EBPα axis Promotes TAM Polarization to Regulate Cancer Cell Proliferation and Migration. Mol. Ther. - Nucleic Acids 23, 536–551. 10.1016/j.omtn.2020.12.005 Zheng L. Xu H. Di Y. Chen L. Liu J. Kang L. (2021). ELK4 Promotes the Development of Gastric Cancer by Inducing M2 Polarization of Macrophages through Regulation of the KDM5A-PJA2-KSR1 axis. J. Transl. Med. 19 (1), 342. 10.1186/s12967-021-02915-1 Zheng W. Zhang S. Jiang S. Huang Z. Chen X. Guo H. (2021). Evaluation of Immune Status in Testis and Macrophage Polarization Associated with Testicular Damage in Patients with Nonobstructive Azoospermia. Am. J Rep Immunol 86 (5), e13481. 10.1111/aji.13481 Zhou R. Wu J. Liu B. Jiang Y. Chen W. Li J. (2019). The Roles and Mechanisms of Leydig Cells and Myoid Cells in Regulating Spermatogenesis. Cell. Mol. Life Sci. 76 (14), 2681–2695. 10.1007/s00018-019-03101-9 Zhu F. Luo Y. Bo H. Gong G. Tang R. Fan J. (2021). Trace the Profile and Function of Circular RNAs in Sertoli Cell Only Syndrome. Genomics 113 (4), 1845–1854. 10.1016/j.ygeno.2021.04.022 Zhu Q. Wong A. K. Krishnan A. Aure M. R. Tadych A. Zhang R. (2015). Targeted Exploration and Analysis of Large Cross-Platform Human Transcriptomic Compendia. Nat. Methods 12 (3), 211–214. 10.1038/nmeth.3249