Single-cell transcriptomics analysis of mild cognitive impairment in World Trade Center disaster responders.

INTRODUCTION: Recent research has found that World Trade Center (WTC) responders in their mid-50s have an elevated prevalence of mild cognitive impairment (MCI) that is associated with neural degeneration and subcortical thinning. This article extends our understanding of the molecular complexity of MCI through gene expression profiling of blood.
METHODS: The transcriptomics of 40 male WTC responders were profiled across two cohorts (discovery: nine MCI and nine controls; replication: 11 MCI and 11 controls) using CITE-Seq at single-cell resolution in blood.
RESULTS: Comparing the transcriptomic signatures across seven major cell subpopulations, the largest differences were observed in monocytes in which 226 genes were differentially expressed. Pathway analysis on the genes unique to monocytes identified processes associated with cerebral immune response. DISCUSSION: Our findings suggested monocytes may constitute a key cell type to target in blood-based biomarker studies for early detection of risk of MCI and development of new interventions.

INTRODUCTION

The September 11, 2001 World Trade Center (WTC) terrorist attack was a massive disaster, resulting in long‐term trauma to the responders and survivors. Now, nearly two decades after the attacks, WTC responders are aging, and there has been an increasing awareness that those who worked at the disaster site have problems consistent with mild cognitive impairment (MCI) at midlife. In older individuals, MCI is often an early sign of Alzheimer's disease (AD) or a related dementia (ADRD); however, to date little is known about the etiology of MCI in WTC responders. Consistent with AD, the central deficit in MCI among WTC responders was found to be memory impairment, and neuroimaging studies further found signs of cortical atrophy in responders with possible dementia consistent with AD. Because the potential for neurodegenerative diseases as indicated by MCI are concerning, and because MCI itself is an impairing condition that often is associated with impairment to the capacity to manage daily living, there is a critical need to understand its pathogenesis.

A key benefit of identifying biomarkers of MCI is the potential to facilitate early interventions. Currently, the diagnosis of MCI relies on a combination of clinical assessments, and physical, cognitive, neuropsychological, and neuroimaging exams. These procedures are expensive and time consuming, and do not reliably identify individuals for whom symptoms may be atypical. Thus, there is active research focused on improving the discovery of biomarkers for pre‐dementia screening and diagnosis. Growing evidence on the neuroimmune involvement in dementia has prompted numerous studies to identify putative biomarkers for ADRD using blood, plasma, or serum. , Blood is an attractive source for biomarker profiling because it is inexpensive, minimally invasive, and informative for disease identification.

Gene expression analyses can identify critical downstream biological process (BP) associated with genetic and epigenetic variations and thus can potentially inform efforts to identify biomarkers for ADRD, , of which MCI is a key component. Collection of peripheral blood is far more feasible than collection of brain tissue, and importantly gene expression patterns in blood are consistent with patterns observed in the brain. When considered together, these results suggest that due to strong communication between the central and peripheral immune systems there may be improved chance of seeing molecular signatures consistent with immune reactions to AD in the periphery. Most gene expression studies to date have been performed in AD on whole or peripheral blood samples, including in the large multi‐center Alzheimer's Disease Neuroimaging Initiative (ADNI) study. However, blood is a complex tissue that consists of several populations of cells with distinct gene expression profiles. Alterations in immune regulatory networks are expected to have functional consequences only in some subsets of immune cells. Thus, analyses of whole blood are likely to weaken the signal. One gene expression study on monocytes using the Nanostring panel targeting 255 inflammation‐related genes found differential expression (DE) in several cytokines at different stages of AD. However, this study used a targeted panel of genes, and it remains unclear whether other classes of genes or immune cell subpopulations provide distinct gene expression patterns.

High throughput single‐cell RNA‐Sequation (sncRNA‐Seq) has emerged as the most powerful technique for describing the transcriptomic landscape at single‐cell resolution. Although sncRNA‐Seq is an attractive platform to characterize the heterogeneity in cell populations, it is currently cost‐prohibitive for large sample experiments. In 2018, the New York Genome Center pioneered a cell hashing technique coupled with cellular indexing of transcriptomes and epitopes by sequencing (CITE‐seq), which reduces the experimental costs of sncRNA‐Seq substantially. , CITE‐Seq is a method that combines highly multiplexed protein marker detection with unbiased transcriptome profiling of single cells and has been shown to achieve a more detailed characterization of cellular phenotypes.

Our goal in this study was to characterize the landscape of immune cell subpopulations in peripheral blood mononuclear cells (PBMC) and identify gene expression biomarkers associated with MCI in a population of WTC responders. To the best of our knowledge, this is the first study using the CITE‐Seq technology to understand cellular heterogeneity in PBMC in MCI at midlife in WTC responders to understand cell‐specific gene expression, to identify common and distinct signatures and pathways associated with MCI within each cell subpopulation.

RESEARCH IN CONTEXT

Systematic review: The authors performed a literature review encompassing published articles and abstracts investigating the challenges in biomarker studies for early detection of Alzheimer's disease and related dementia (ADRD). Previous studies suggested that blood‐based gene expression profiling is a promising approach for biomarker discovery in ADRD.

Interpretation: Using single‐cell transcriptomics in peripheral blood mononuclear cells, results revealed that monocytes showed the largest differential gene expression (DGE) in mild cognitive impairment (MCI), consistent with the role of monocytes in brain–immune communication. The authors also found pathways involved in inflammation and viral infections among MCI DGE in monocytes.

Future directions: Further studies include a follow‐up clinical diagnosis of MCI to determine whether the DGE in monocytes can predict disease onset and progression.

METHODS

Participants and clinical assessment

Participants were recruited from the Stony Brook WTC‐Health Program, which conducts the only cognitive monitoring study of WTC responders. We included only males because females show notably different gene expression patterns than males and <10% of responders in the Stony Brook cohort were female. Participants completed a battery of medical and psychological tests annually as part of the Centers for Disease Control (CDC)‐funded monitoring program established in 2002. The initial exams included a detailed WTC exposure history questionnaire. The current study was approved by Stony Brook University Institutional Review Board. Written informed consent was obtained.

MCI was diagnosed following the National Institute on Aging‐Alzheimer's Association criteria as has been previously described in this population. Exclusion criteria for this analysis were cancer of the brain, neurologic diagnoses including AD and other related dementias made prior to initial cognitive assessments, WTC‐related head injuries or head injuries endured while in military service, as well as those with psychotic conditions or substance use disorders. A total of 40 WTC participants were included (20 MCI and 20 controls). All of the participants were non‐smokers and White. The mean age was 55.83 (standard deviation [SD] = 6.07; Table 1).

TABLE 1

Clinical characteristics of samples in discovery and replication sets. The P‐values for comparing age of mild cognitive impairment (MCI) to control were computed from t‐test

Group (cohort)	N	Age mean	Age SD	P
MCI (all)	20	55.15	4.93	.490
Control (all)	20	56.5	7.10
MCI (discovery)	9	53.33	5.29	.866
Control (discovery)	9	53.78	5.72
MCI (replication)	11	56.64	4.30	.438
Control (replication)	11	58.73	7.58

The 40 sets of blood samples were divided into two batches. The first batch (discovery) consisted of nine MCI and nine controls, whereas the second batch (replication) consisted of 11 MCI and 11 controls.

CITE‐Seq profiling

Cryopreserved PBMCs were thawed and washed twice. Cells were stained with barcoded antibodies as previously described for CITE‐seq and Cell Hashing. The cDNA, ADT, and HTO libraries were quantified using Picogreen and Fragment Analyzer then pooled with appropriate percentages. The final library pools were then sequenced on the Illumina NovaSeq. Additional details are provided in supporting information.

CITE‐Seq data preprocessing

Analysis for counting UMIs per cell per gene was done following the Drop‐seq tools protocol, alignment was performed using STAR aligner v2.5.2a, HTO and ADT nUMI per cell barcode per antibody were quantified using CITE‐Seq‐Count. Seurat V3.0.1 was used for all downstream analysis. Cells with >500 genes detected and mitochondrial rate <15% were retained, resulting in a total of 61,046 cells. The HTO and RNA raw counts were normalized. Mixtures of two Gaussian distributions were fitted to each ADT marker to determine the cutoff for identifying cells expressing the marker. The cell identity was determined using combinations of ADT markers for the major cell types BCell, CD4T, CD8T, monocytes, NK, DC, and DNT (Table S1 in supporting information). A subset of monocytes was further divided into classical and non‐classical monocytes. Among the cells which were untagged by the ADT markers, we adapted the method described in Diaz‐Mejia et al. to tag these cells. The UMAP method was used to visualize cell clusters. For each cell type, genes with >90% zero were filtered, resulting in a total of 13,692 genes analyzed. For each cohort, the cells across the samples were pooled for subsequent analyses. Two‐sample t‐test was used to determine whether the proportions differ between MCI and healthy exposed controls, and P < .05 was considered statistically significant. Additional details are provided in supporting information.

Differential expression analysis

DE analysis was carried out using the Wilcoxon rank sum test comparing MCI to control for each cell type within discovery and replication sets, respectively. Genes with Bonferroni‐adjusted P < .05 in discovery were identified. Among these genes, genes with P < .05 in replication and consistent effect size directions were considered replicated. DE genes unique to a specific cell type were defined as those with P > .2 in other cell types in both discovery and replication cohorts, to ensure that the unique DE genes associated with MCI for a specific cell type were not marginally significant in other cell types. The Pearson correlation coefficients were computed on the inverse Gaussian‐transformed P‐values to compare the global DE patterns associated with MCI across cell types.

Sensitivity analyses were conducted via two additional statistical models to assess the effect of age and individual of origin for cell. The consistency across the different statistical methods was compared. Additional details are provided in supporting information.

Pathway and gene ontology analyses

Pathway and gene ontology analyses were carried out using the over‐representation via the Bioconductor package clusterProfiler on the DE genes within each cell type using the functions enrichGO and enrichKEGG. In total, 6402 gene ontologies (GO) including biological processes, molecular functions, and cellular components and 316 KEGG pathways (the range of genes per gene set was 15 to 500) were tested. Statistically significant gene sets corresponded to those with false discovery rate < 0.05 from over‐representation analyses. The significant GO terms were clustered using REVIGO to reduce functional redundancies.

Data availability

The CITE‐Seq data will be available at Synapse (https://www.synapse.org/#!Synapse:syn22855256, https://doi.org/10.7303/syn22855256) upon publication.

RESULTS

Comparison of discovery and replication cohorts

MCI and control patients had comparable age within the discovery and replication cohorts. However, participants in the replication cohort were on average 4.13 years older than the discovery cohort (P = .029; Table 1).

CITE‐Seq cell identity

The average number of cells retained after filtering were 1434 (SD = 495) and 1602 (SD = 278) per sample in discovery and replication set, respectively. Figure 1A,B shows the UMAP plots of the cells by cell type, which indicates that the different cell types were reasonably well segregated in both cohorts. The proportions of cell types aggregated over all samples are provided in Table S1. The proportions of NK, DC, total monocytes, and classical monocytes were higher in the replication cohort, whereas the proportion of CD4T was lower at P < .05. Figure S1 in supporting information compares the proportions of cell types by case (MCI)/control status. The proportion of DNT was lower in MCI in the discovery cohort (P = .029); however, it was not statistically different in the replication cohort. On the other hand, the proportions of total monocytes and classical monocytes were lower in MCI in the replication cohort at P = .026 and .032, respectively, but not significant in the discovery cohort.

FIGURE 1

UMAP depicting the clusters of cells by cell types, (A) discovery set, (B) replication set

DE analysis associated with MCI

DE analysis identified 34, 234, 126, 81, 14, 32, 444, 416, and 55 genes to be associated with MCI in the discovery set in Bcell, CD4T, CD8T, NK, DC, DNT, total monocytes, classical monocytes, and non‐classical monocytes, respectively. Among these, 15/34 (Bcell), 70/234 (CD4T), 27/126 (CD8T), 31/81 (NK), 7/14 (DC), 6/32 (DNT), 226/444 (total monocytes), 215/415 (classical monocytes), and 10/55 (non‐classical monocytes) were replicated as shown in Table 2. All the replicated DE genes within DC and DNT were upregulated in MCI, whereas 76%, 67%, 46%, 48%, and 41% of the replicated DE genes were upregulated in MCI in total monocytes, Bcell, CD4T, CD8T, and NK cells. The complete list of statistically significant DE genes is provided in Table S2 in supporting information.

TABLE 2

Number of genes at Bonferroni‐adjusted P < .05 in discovery cohort, replicated at P < .05 in replication cohort and proportion of upregulated genes

	Significant genes in discovery cohort	Replicated genes	Proportion of upregulated genes
Bcell	34	15	0.667
CD4T	234	70	0.457
CD8T	126	27	0.481
NK	81	31	0.419
DC	14	7	1
DNT	32	6	1
Monocytes	444	226	0.761
Classical monocytes	416	215	0.744
Non‐classical monocytes	55	10	0.7

Among the list of genes associated with MCI, RP11‐742N3.1, SNHG5, HLA‐B, and MT‐ATP6 were identified in monocytes, Bcell, CD4T, CD8T, and NK. SNHG5 was also identified in DNT, whereas MT‐ATP6 was identified in DC. On the other hand, all but one (TRIB1) gene among the 10 DE genes in non‐classical monocytes were also detected in classical monocytes.

Among the 226 DE genes associated with MCI in total monocytes, 85 genes were unique to total monocytes, that is, not DE in other cell types, of which 80 were upregulated (Table S3 in supporting information). The top 10 genes were HLA‐DRB6, CCL3, MARCKS, NR4A1, GBP2, BCL2A1, PLEK, SRGN, ICAM1, and GUSBP3; all were upregulated in MCI.

Results from sensitivity analyses assessing the effect of age and individual of origin for cell are provided in supporting information, and show that these factors had negligible effects on the DE analysis in our study.

Pathways and GO associated with MCI genes

The number of significant KEGG pathways and GO gene sets involved in BP, cellular component (CC) and molecular function (MF) associated with replicated DE genes in MCI for each cell type are provided in Table 3. The complete lists describing the KEGG and GO terms are given, respectively, in Tables S4 and S5 in supporting information.

TABLE 3

Number of KEGG and gene ontologies (GO) biological process (BP), cellular component (CC), and molecular function (MF) gene sets at false discovery rate < 0.05 identified from the replicated differential expression (DE) genes associated with mild cognitive impairment (MCI) within each cell type

	KEGG pathways	GO (BP)	GO (CC)	GO (MF)
Bcell	1	24	13	1
CD4T	1	56	19	1
CD8T	1	34	6	10
NK	22	40	36	10
DC	0	10	1	0
DNT	6	180	0	18
Monocytes	31	438	35	24

Three KEGG pathways were significantly associated with the 85 genes unique to total monocytes, namely Epstein‐Barr virus (EBV) infection, toll‐like receptor (TLR), and NOD‐like receptor (NLR) signaling pathways. On the other hand, 97 GO BP and four GO CC were significantly associated with these 85 genes. Cluster analysis of the significant GO terms using REVIGO to reduce functional redundancies retained 36 GO BP terms summarized in seven clusters as shown in Figure 2, namely cellular response to biotic stimulus, regulation of protein secretion, podosome assembly, viral transcription, cytokine metabolism, antigen processing, and presentation of peptide antigen via MHC class I and ovarian follicle development.

FIGURE 2

Treemap summarizing the clustered gene ontologies (GO) biological process (BP) terms among the 85 genes unique to total monocytes

DISCUSSION

There is a need to develop better blood‐based indicators of cellular functioning in ADRD and also to improve our understanding of the etiology of subtypes of ADRD including WTC‐related cognitive impairment. The present study helped to fill these gaps by profiling transcriptomics in 40 WTC responders with MCI using the CITE‐Seq platform in PBMCs. Using this protocol, we identified DE patterns associated with MCI in seven immune cells subsets, namely, Bcell, CD4T, CD8T, NK, DC, DNT, and total monocytes. Four genes, namely RP11‐742N3.1, SNHG5, HLA‐B, and MT‐ATP6 were common across monocytes; Bcell, CD4T, CD8T, and NK, and all were upregulated in MCI in this study. These results suggest that biomarkers may be developed using the CITE‐Seq platform, and also suggest that changes in the monocyte population warrant more research in this population to determine which indicators of monocytic dysregulation identified here are prognostic of increased risk of incident dementia in WTC responders and in the general population.

Genes differentially expressed in these analyses may play a role in a number of ADRD subtypes. For example, MT‐ATP6 is a mitochondrial gene whose mutations have been previously associated with ataxia, cognitive dysfunction, neuropathy, and seizures, , whereas overexpression of this gene has been observed in blood of early AD and MCI. SNHG5 is a long non‐coding small nucleolar RNA highly expressed in cancer and promotes cell proliferation and apoptosis. Although the link between SNHG5 and AD has yet to be established, this gene was previously found to be implicated in multiple sclerosis. On the other hand, HLA‐B is a histocompatibility complex, class I B gene and plays an important role in the immune system by differentiating the host's proteins from foreign viral or bacterial proteins. Genetic variants in HLA‐B were also found to be associated with AD and were hypothesized to play a protective role by eliminating pathogens that contribute to dementia.

The DE analysis showed that the largest difference occurred in the monocyte cell subsets comparing MCI to control. DE genes in total monocytes were mainly contributed by the dominant classical monocytes subset. A total of 226 genes were differentially expressed in monocytes. These included several cytokines, chemokines, and inflammatory response genes such as IL1B, IL8, CCL3, CCL5, CEBPB, NR4A1, and CD40. All, except CCL5 were upregulated in MCI. Several of these genes have been found to be involved in the inflammatory responses in AD. , For example, upregulation of IL1B and IL8 was associated with amyloid plaque progression and was involved in amyloid beta (Aβ)‐mediated inflammatory processes. The top five genes unique to total monocytes were HLA‐DRB6, CCL3, MARCKS, NR4A1, and GBP2. Genetic variants of HLA‐DRB6, a human leukocyte antigen, have been found to be a risk factor for AD. The cytokine CCL3 was previously found to be expressed by neurons and microglia in AD post mortem brains, and upregulated in both amyloidosis and tauopathy deposits in experimental models (reviewed in Marciniak et al. ). MARCKS is involved in cell morphology, motility, and neural development, and may be implicated in Aβ generation. Specifically, MARCKS is associated with AD pathology via the phosphorylation process, and is a marker of protein kinase C, which was previously found to be activated in microglia and dystrophic neurites by a Aβ in AD brains. NR4A1, a nerve growth factor IB gene, plays a key role in mediating inflammatory responses in macrophages and is involved in Aβ precursor protein metabolism and tau phosphorylation. GBP2 is a gene encoding interferon‐induced guanylate‐binding protein 2. Interferons are cytokines that have antiviral effects. Risk variants of GBP2 were previously found to be potentially associated with AD. Two other GBP genes, namely GBP1 and GBP5, were also unique to monocytes in this study.

EBV infection was one of the three KEGG pathways associated with the genes unique to total monocytes. The other two were TLR and NLR signaling pathways. Elevated levels of EBV‐specific antibodies have been shown to be associated with an increased relative risk for developing AD. On the other hand, the formation of Aβ deposits in the brain, a key characteristics of AD, is linked to a microglial‐mediated inflammatory response. TLRs are involved in this inflammatory response and contribute to AD pathogenesis. Similarly, NLRs are involved in central nervous system inflammation, and are associated with the pathophysiology of AD. Both TLRs and NLRs are important receptors that mediate immune recognition. On the other hand, the GO terms clusters associated with the genes unique to monocytes include podosome assembly, viral transcription, cytokine metabolism, and antigen processing/presentation of peptide antigen via MHC class I. Podosomes consist of a core rich in actin surrounded by adhesion and scaffolding proteins and have been found to be implicated in microglia motility and migration. On the other hand, viral transcription was consistent with the KEGG pathway analysis result, as well as evidence that pathogenic viruses may contribute to the onset and progression of AD. Cytokine metabolism and antigen processing/presentation of peptide antigen are involved in immune and inflammatory responses, all of which are important in the pathogenesis of AD.

Our findings were in line with existing evidence that monocytes play an important role in AD. Monocytes are thought to play a pivotal role in explaining the bidirectional brain–immune relationship and signaling triggered by psychological stress. Emerging evidence also implicates a novel neuroimmune circuit involving microglia activation and sympathetic outflow to the peripheral immune system that reinforces stress‐related behaviors by facilitating the recruitment of inflammatory monocytes to the brain. Additionally, microglial depletion can trigger peripheral macrophage engraftment into the CNS, and these cells retain unique functional identity. This suggests that monocytes could influence disease progression in the CNS, and thus future studies focusing on monocytes might inform potential treatment strategies to slow disease progression.

Strengths and limitations

The current study had several strengths, including the first to profile transcriptomics at single cell resolution in blood of MCI using CITE‐Seq and replication of results in an independent subsample. Nonetheless, our findings must be considered in the context of several limitations. First, our study is cross‐sectional, which can establish concurrent associations between gene expression and MCI, but the direction of the associations cannot be determined, such as whether the identified signatures can predict onset of AD. Longitudinal studies are needed to determine the direction of the effects we observed, as well as follow‐up clinical diagnosis of AD to determine the predictive power of MCI‐transcriptomics in disease onset and progression. Second, our current sample was all male, thus it is unclear to what extent results generalize to females. Third, some cell subsets exist in very small fractions, which may not be captured by our current CITE‐Seq resolution. Finally, another shortcoming is the absence of characterization of the neuropathology of WTC‐related neurocognitive disorders.

Implications

In conclusion, the current study identified common and distinct transcriptomic signatures associated with MCI in seven immune cell subsets at single‐cell resolution, indicating that the cell subpopulations may provide a valuable and to some degree independent source of information to identify the biomarker signature for early detection of AD compared to whole blood. In particular, monocytes showed the most differential gene expression in MCI, in line with evidence that monocytes play a pivotal role in mediating the crosstalk between central and peripheral systems via transduction through the blood‐brain barrier. Together with the findings from pathways analysis, the transcriptomic profiles point to inflammation and viral infections. These results add to growing evidence suggesting that intervention strategies that target inflammatory responses and/or infectious agents early in the disease may help to prevent or slow the progression of WTC‐related neurocognitive disorders.

CONFLICTS OF INTEREST

The authors reported no biomedical financial interests or other potential conflicts of interest.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.