Identification of hub genes associated with cognition in the hippocampus of Alzheimer's Disease.

Alzheimer's Disease (AD) is a neurodegenerative disease featured by cognitive impairment. This bioinformatic analysis was used to identify hub genes related to cognitive dysfunction in AD. The gene expression profile GSE48350 in the hippocampus of AD patients aged >70 years was obtained from the Gene Expression Omnibus (GEO) database. A total of 96 differentially expressed genes (DEGs) were identified, and subjected to Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses; a protein-protein interaction (PPI) network was constructed. The DEGs were enriched in synapse-related changes. A protein cluster was teased out of PPI. Furthermore, the cognition ranked the first among all the terms of biological process (BP). Next, 4 of 10 hub genes enriched in cognition were identified. The function of these genes was validated using APP/PS1 mice. Cognitive performance was validated by Morris Water Maze (MWM), and gene expression by RT-qPCR, Cholecystokinin (CCK), Tachykinin precursor 1 (TAC1), Calbindin 1 (CALB1) were downregulated in the hippocampus. These genes can provide new directions in the research of the molecular mechanism of AD.

Introduction

Alzheimer’s Disease (AD) is a neurodegenerative disease featured by cognitive decline, disorientation, abnormal behavior, impaired mobility, hallucinations, seizures, and emotional instability [1] [2]. The prevalence of AD is positively related to age, with an incidence of 10% to 30% in patients aged over 65 and doubling every 10 years after 60 [3]. However, the mechanism underpinning AD remains controversial. Given its close relevance with age, amyloid-beta (Aβ) and tau protein accumulated in the aging brain may serve as critical players in AD [4].

AD patients show short-term memory deficits. In the central nervous system (CNS), the hippocampus plays a vital role in the formation and consolidation of memory by long-term potentiation (LTP) and synaptic remodeling [5]. Therefore, it has been researched extensively in the field of cognition. In the pathology of AD, the hippocampus suffers damage arising from abnormal amyloid beta-peptide and p-tau.

Microarray analysis or high-throughput sequencing have been employed to mine the potential biomarkers or therapeutic targets of AD. Berchtold NC et al. used a microarray method to obtain RNA expression profiles in different postmortem brain regions of AD patients, including the postcentral gyrus, superior frontal gyrus, entorhinal cortex and hippocampus [6]. GSE48350 has also been used in some other studies [7-12]. In these studies, however, other factors (e.g., immunity, aging, toxicity) and samples (e.g., cortex, blood) are analyzed.

Here we aim to identify genes responsible for the cognitive decline in the hippocampus.

Method

Data source

The original microarray dataset (.CEL file) was downloaded from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) in May 2020. GSE48350 was based on GPL570 platform (HG-U133_Plus_2), composed of postmortem brain samples from the AD patients and normal people. The samples meeting the following criteria were included: (1) from the hippocampus; (2) from aged patients (over 70-year-old). As a result, a total of 40 samples from the hippocampus were restricted into further analysis, including 19 samples from AD patients and 21 samples from normal people as a control group. (AD group: GSM1176221, GSM1176222, GSM1176223, GSM1176224, GSMx1176225, GSM1176226, GSM1176227, GSM1176228, GSM1176229, GSM1176211, GSM1176212, GSM1176213, GSM1176214, GSM1176215, GSM1176216, GSM1176217, GSM1176218, GSM1176219, GSM1176220; Control group: GSM300169, GSM300170, GSM300178, GSM300182, GSM300185, GSM300205, GSM300215, GSM300235, GSM300268, GSM300333, GSM300168, GSM300171, GSM300172, GSM300190, GSM300193, GSM300197, GSM300223, GSM300227, GSM300239, GSM300329, GSM300339).

Data preprocessing

Data processing was conducted using R (version 4.0.2) [13]. The R package Affy was utilized to read raw data [14], then the Robust Multichip Averaging (RMA) algorithm was used to normalize the data about the samples [15], insuring that each sample was comparable. Probe ID was transferred to gene name according to the platform’s annotation. When a gene corresponded to multiple probes, the median expression value in various probe IDs was taken to represent the duplicated gene’s expression.

Identification of differentially expressed genes (DEGs)

The LIMMA package is a R’s extension for the analysis of gene expression data via linear models [16]. It was used to compare the gene expression between AD group and control group in our study. Differentially expressed genes (DEGs) were identified according to P value < 0.05 and the absolute value of Fold Change > 2 (log2 |Fold Change| > 1).

Enrichment analysis of DEGs

The clusterProfiler is a package for functional and comparative study in R [17]. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed via clusterProfiler package in our research. GO analysis covers three aspects about cellular component (CC), biological function (BP) and molecular function (MF). KEGG enrichment analysis covers pathways responsible for molecular interactions and biochemical reactions. A threshold of P < 0.05 and enrichment counts > 2 were set as thresholds. Replicated terms with different descriptions were removed.

PPI network of DEGs

Search Tool for the Retrieval of Interacting Genes (STRING) database(https://string-db.org/) is an online tool to calculate PPI connections and identify the potential relationship among DEGs [18]. The Cytoscape software was used to visualize the results [19]. The clusters of PPI were identified by MCODE tool, a plugin that defines highly interconnected regions in a network. The hub genes were screened out by the plugin cytoHubba, a tool that identifies hub genes through multiple algorithms [20]. The intersections of four ranking algorithms, including Maximal Clique Centrality (MCC), Density of Maximum Neighborhood Component (DMNC), Maximum Neighborhood Component (MNC), Degree, were considered as hub genes.

Animal experiments

All experiments were conducted observing related laws and guidelines. All procedures were approved by the Nanjing University Institutional Animal Care and Use Committee. All processing methods were in line with Directive 2010/63/EU. Male APP/PS1 transgenic mice (aged 18 months) and wild-type littermates were purchased from Ziyuan Laboratory Animal Technology Co. Ltd. (Hangzhou, China) in June 2020. All mice were housed in separated cages with free access to food and water in a room under a 12 h light/dark cycle, with a temperature of 22 ± 2°C and a humidity of 50–70%. The behavioral test initiated after all the mice adjusted to the environment for 2 weeks.

Morris water maze

The MWM test contained a training phase and a probe test phase [21]. The training phase was performed for five consecutive days. Each mouse completed four trials in different quadrants. In each trial, a mouse was allowed to search the platform within 60s and stay on the platform for 15s afterward. The mouse that failed to find the platform within 60s was manually guided to the platform and stayed at the platform for 30s. 60s was recorded as the escape latency. The probe test was performed on the 6th day, the platform was removed, and the mice were placed into 2 selected quadrants that did not contain the platform and were allowed to swim for 60s. The percentage of time spent in the target quadrant and platform crossings were recorded and analyzed.

RNA extraction and RT-qPCR

Mice were sacrificed after 2 hours of MWM. Total RNA was extracted from the hippocampus of the mice using TRIzol reagent (Vazyme, China), then cDNA was synthesized using HiScript Q RT SuperMix (Vazyme, China) according to manufacturer’s protocol. The SYBR qPCR master mix (Vazyme, China) was used to detect the relative mRNA expression in the QuantStudio 5 instrument (Applied Biosystems Inc.). The relative expression level was determined by the 2−ΔΔCt method and normalized to loading control β-actin. Primers of four hub genes were shown as follow:

β-actin forward primer 5ʹ-CTGTCCCTGTATGCCTCTG-3ʹ, reverser primer 5ʹ-ATGTCACGCACGATTTCC-3ʹ;

CCK forward primer 5ʹ-ATACATCCAGCAGGTCCGCAA-3ʹ, reverse primer 5ʹ- CAGACATTAGAGGCGAGGGGT-3ʹ;

TAC1 forward primer 5ʹ-CTAAATTATTGGTCCGACTG-3ʹ, reverse primer 5ʹ-TTCTGCATTGCGCTTCTTTC-3ʹ;

CALB1 forward primer 5ʹ-TCTGGCTTCATTTCGACGCTG-3ʹ, reverse primer 5ʹ-ACAAAGGATTTCATTTCCGGTGA-3ʹ;

CNR1 forward primer 5ʹ-AAGTCGATCTTAGACGGCCTT-3ʹ, reverse primer 5ʹ-TCCTAATTTGGATGCCATGTCTC-3ʹ.

Statistical analysis

All data were expressed as the mean ± SD. Data collected from the MWM training phase were analyzed using a repeated-measures ANOVA, followed by Student’s t-test for between-group comparison. The RNA expression of hub genes was analyzed using Student’s t-test. All data were analyzed with SPSS 26.0 software, and P < 0.5 was considered a significant difference.

Results

In this study, we identified hub genes associated with cognitive decline in AD patients by bioinformatic tools, including GO and KEGG enrichment analyses, PPI network analysis. Then we validated our findings by RT-qPCR on APP/PS1 mice. The aim was to delve into the molecular mechanism of cognitive dysfunction in AD and provide new options of diagnosis and treatment.

Identification of DEGs

The dataset GSE48350 contained the expression profiles of genes in the hippocampus of 21 AD patients and 19 normal people. In total, 96 DEGs were identified between the two groups, with 2 upregulated and 94 downregulated. The expression levels in two upregulated genes increased by over twofold (ABCA6, CP). Of the 94 downregulated genes, the expression levels decreased by over fourfold in 2 (CALB1, TAC) and over onefold in 92 (PDYN, BDNF, etc.). All DEGs were visualized in volcano plot and heatmap (Figure 1).

Figure 1.

The visualization of differential gene expression profile (a) A volcano plot of the expression of DEGs in the whole gene profile, significantly upregulated genes are indicated in red, downregulated in green (b) Cluster pattern of Heat map concerning DEGs, high level expressions are indicated in red, low level in blue

GO and KEGG enrichment analyses

The 96 DEGs were then submitted to GO and KEGG pathway enrichment analyses. In the aspect of cellular component (CC), the DEGs were enriched in presynapse, axon part, distal axon (Figure 2(a)). In the aspect of molecular function (MF), the DEGs were enriched in receptor ligand activity, structural constituent of synapse, neurotransmitter receptor activity (Figure 2(b)). In the aspect of biological process (BP), the DEGs were significantly enriched in the regulation of membrane potential, regulation of trans-synaptic signaling, and cognition (Figure 2(c)). In KEGG pathway enrichment analysis, the DEGs were enriched in neuroactive ligand–receptor interaction (Figure 2(d)).

Figure 2.

GO and KEGG pathway enrichment results of DEGs: (a) cellular component, (b) molecular function, (c) biological process, (d) KEGG pathways

PPI

After removing isolated proteins from the PPI network, 69 proteins were identified from the STRING database. Using MCODE, a cluster composed of PENK, CNR1, CALB1, TAC1, SST, GNG, CCK, BDNF, GAD1 and PDYN was identified (Figure 3(a)). The cluster was then subjected to GO and KEGG analyses. In the aspect of CC, the DEGs in the cluster were enriched in presynapse, distal axon, axon part (Figure 3(b)), neuronal cell body. In the aspect of MF, these DEGs were enriched in receptor ligand activity, hormone activity, G protein-coupled receptor binding (Figure 3(c)). In the aspect of BP, these DEGs were enriched in cognition, sensory perception of pain (Figure 3(d)). In KEGG pathway enrichment analysis, the DEGs were enriched in neuroactive ligand–receptor interaction (Figure 3(e)).

PPI network and GO/KEGG pathway enrichment results of genes in the cluster (a) DEGs in Protein-protein interaction network. Nodes in red represent upregulated genes and nodes in blue represent downregulated genes, the internal circle indicates the cluster (b) cellular component, (c) molecular function (d) biological process, (e) KEGG pathway

Continued

Hub genes

The tool cytoHubba in Cytoscape was used to screened out the hub genes. The top 20 genes in four ranked algorithms (MCC, DMNC, MNC, Degree) were figured out, and the intersections between the four algorithms were considered as the hub genes (Table 1, Figure 4(a)). Finally, 10 genes were screened out. Of them, 4 genes, including Cholecystokinin (CCK), Cannabinoid receptor type 1 (CNR1), Tachykinin precursor 1 (TAC1), Calbindin 1 (CALB1), were identified related to cognitive dysfunction in BPs of Gene Ontology (Figure 4(b)).

Table 1.

Top 20 genes and the intersection in the cytoHubba. MCC, Maximal Clique Centrality. DMNC, Density of Maximum Neighborhood Component. MNC, Maximum Neighborhood Component

		Hub genes
Algorithms	MCC	CCK, SST, TAC1, PDYN, BDNF, GAD1, PENK, GNG2, CNR1, CALB1, GRP, SYN2, CHRM1, SCG2, HTR2C, GABRG2, GAP43, SYNPR, PTGER3, CHGB
	DMNC	CCK, TAC1, PDYN, SST, PENK, GNG2, CNR1, GRP, GAD1, CALB1, SCG2, PTGER3, CHGB, CHRM1, NEFH, INA, NEFM, CALY, SVOP, OLFM3
	MNC	CCK, SST, TAC1, BDNF, PDYN, GAD1, PENK, GAP43, SYN2, SYNPR, CALB1, CNR1, GABRG2, GNG2, GRP, NEFL, STMN2, SCG2, SNAP91, RAB3C
	Degree	CCK, SST, BDNF, TAC1, PDYN, GAD1, PENK, GAP43, SYN2, SYNPR, CNR1, GABRG2, CALB1, STMN2, GNG2, NEFL, CHGB, SCG2, SNAP91, RAB3C
Intersection		CNR1, GAD1, PENK, CALB1, GNG2, CCK, PDYN, TAC1, SST, SCG2

Figure 4.

Venn diagrams for identification of hub genes (a) 10 hub genes were identified via the intersection of four algorithms (b) 4 genes related to the cognition according to biological process of Gene Ontology database. CALB1, Calbindin 1. TAC1, Tachykinin precursor 1. CCK, cholecystokinin. CNR1, cannabinoid receptor type 1

Impaired cognition of APP/PS1 mice in the Morris water maze

In MWM test, repeated measures ANOVA indicated a significant difference between the two groups (F[1,8] = 44.134, P < 0.0001). Wild-type mice showed a more pronounced downward trend than APP/PS1 in five consecutive days. In each day, APP/PS1 mice spent more time in finding the platform than wild-type mice: D1(60 ± 0 vs 47.90 ± 10.86, P = 0.038), D2(57.90 ± 6.50 vs 35.56 ± 12.48, P = 0.009), D3(56.48 ± 6.30 vs 28.30 ± 9.20, P < 0.0001), D4(49.54 ± 6.92 vs 22.38 ± 5.35, P < 0.0001), D5(49.17 ± 10.06 vs 16.46 ± 8.94, P = 0.001). In the probe test, there was a decline in the number of crossing removed platform (P = 0.018) and the time stayed in the platform quadrant (P = 0.011) (Figure 5). So, APP/PS1 mice showed a cognitive impairment compared with wild-type littermates.

Figure 5.

Validation of APP/PS1 mice. A. APP/PS1 mice spent more time finding the platform than wildtype mice in training phase, Day1(60 ± 0 vs 47.90 ± 10.86, P = 0.038), Day2(57.90 ± 6.50 vs 35.56 ± 12.48, P = 0.009), Day3(56.48 ± 6.30 vs 28.30 ± 9.20, P < 0.0001), Day4(49.54 ± 6.92 vs 22.38 ± 5.35, P < 0.0001), Day5(49.17 ± 10.06 vs 16.46 ± 8.94, P = 0.001). B. In probe test, the number of crossing removed platform location declined in APP/PS1 mice (P = 0.018) C. In probe test, the time stayed in the platform quadrant declined in APP/PS1 mice (P = 0.011)

Validation of hub genes

To validate the reliability of bioinformatic analysis, RT-qPCR was performed for hippocampal samples harvested from mice (five mice and control littermate). The expression of four hub genes was measured (Figure 6). Compared with those in the normal tissue, the expression levels of four genes were decreased in the tissue samples of AD mice, with significance in three genes (P < 0.05).

Figure 6.

mRNA expression of four hub genes in two groups. CALB1, CCK, TAC1 showed a significant difference between the two groups. CALB1(P = 0.01), TAC1(P = 0.039), CCK (P = 0.02), CNR1 (P = 0.29)

Discussion

Therapeutic targets for human AD remain to be discovered. Berchtold NC et al. made use of a microarray method to analyze RNA expression profiles in different brain regions of postmortem AD patients and normal people, including the postcentral gyrus, superior frontal gyrus and entorhinal cortex [6]. Compared with other AD microarray analyses, theirs contained the data of RNA expression in brain samples from AD patients. Previous studies have shown an increasing incidence of cognitive disorders in patients aged over 70 years [22-24], so we chose the data of these patients in our bioinformatic analysis, aiming to identify AD-related genes in the hippocampus.

After data processing, we found 96 DEGs involved in AD, including 94 downregulated and 2 upregulated (Figure 1). Then, GO and KEGG analyses were performed for these DEGs (Figure 2). Three aspects of results (CC, MF, BP) were all related to the synapse, a communicational unit that participates in neuronal interactions. Numerous basic and clinical studies have demonstrated the direct role of synapse in cognitive impairment. At the biological level, the synaptic function is distorted in the settings of neurogenerative diseases, especially the ‘synaptic plasticity’ in AD [25]. Synaptic plasticity is closely associated with the formation of learning and memory [26]. Several autopsy studies have reported that AD patients have a lower level of synapses than healthy people [27-29]. Neurogenesis declines in the hippocampus of AD patients [30]. Synaptic communications involve a series of neuropeptides that act in receptor–ligand interaction to regulate the postsynaptic membrane potential [31]. These neuropeptides also show neuroprotective functions against AD [32]. A message is transmitted through neurotransmitters, such as Acetylcholine (ACh), GABA, NMDA, which is stored in presynaptic vesicles. Any abnormal transmission results in disordered trans-synaptic signaling and membrane potential. The membrane potential in the CNS may contribute to memory formation [33]. Then the protein–protein interaction (PPI) analysis was conducted, A protein cluster including PENK, CNR1, CALB1, TAC1, SST, GNG, CCK, BDNF, GAD1 and was identified and was transferred to the GO and KEGG analyses (Figure 3(a-e)). These 10 genes were enriched in CC and MF terms similar to those of the 96. Interestingly, in BP terms, ‘cognition’ ranked the first, indicating that the cluster is most closely related to cognitive function. As a consequence, we intended to identify hub genes related to cognition in the whole PPI. Four genes were enriched in the term of cognition, including CCK, CALB1, TAC1 and CNR1 (Figure 4). Then we confirmed the reliability of bioinformatic analysis using APP/PS1 mice with a cognitive impairment phenotype (Figure 5), and the expression of genes in the hippocampus were confirmed via RT-qPCR method (Figure 6). Taken together, these hub genes play important roles in AD-related cognitive impairment.

CCK encodes the cholecystokinin and works as a pancreatic enzyme in the digestive system as known. However, it is also a kind of neuropeptide that’s abundant in the CNS [34]. In AD, CCK may play a protective role. A higher CCK level predicts a better cognitive performance in a linear regression study [35]. Anti-LINGO-1 reverses AD-induced cognitive decline via enriching CCK-GABA interneuron. On the contrary, CCK-KO mice showed worse cognition than wild-type littermates [36]. CNR1 encodes cannabinoid receptor type 1 (CB1) and is distributed in the brain. Lower CB1 availability is correlated with cognitive decline in Parkinson’s Disease [37]. A CB1 knockout model manifests an altered synaptic ingredient and memory impairment [38], and cannabinoid receptor type1 agonist alleviates the cognitive dysfunction [39]. CALB1, distributed in excitatory neurons and inhibitory interneurons in the hippocampus, encodes the calcium-binding protein Calbindin 1. In a CALB1-KD model, reduced level of hippocampal CALB1 contributed to spatial memory loss [40]. In AD, the level of calbindin is related to the severity of cognitive decline, due to the epigenetic modification at the promotor of CALB1 [41]. Furthermore, calbindin is upregulated in the early stage and downregulated in the late stage, making it an ideal biomarker [42] . Tachykinin precursor 1 (TAC1) encodes a precursor protein to produced substance P, neurokinin A, etc. Among them, substance P has been studied most widely. It acts as a neurotransmitter to regulate the activities of the CNS, and also modulates learning and memory consolidation in the hippocampus. Substance P also decreases in the hippocampus of postmortem AD patients [43-45]. It has been reported that the substance P can be modulated to improve cognitive performance in an Alzheimer’s model [46]. Substance P also induces neuroinflammation. An antagonist of Substance P relives hippocampal inflammation and spatial memory impairment [47]. Besides, we found that both CCK and CNR1 were expressed at CCK-expressing basket cells (CCKBC) located in the hippocampal CA1 region. CCKBC targets postsynaptic pyramidal cells with the different neural circuits in different layers. CCK receptors are expressed in postsynaptic pyramid cells. When CCK binds to the receptors, a cascade molecular pathway is activated, ending up with the synthesis of cannabinoids. Cannabinoids interact the receptor located at the presynaptic CCK’s terminal CB1 receptor [48], a process known as the retrograde signaling pathway. Delta9-THC, an exogenous cannabinoid, can disrupt hippocampal network oscillations by binding to CB1R, resulting in a cognitive decline [49] [50]. The electrophysiological activation of CCKBC enhances multiple cognitive behaviors [51]. On the contrary, an abnormal CCKBC wiring pattern disrupts spatial memory [52]. The reduction of CCKBC causes alternations in spatial learning and memory [53].

CCK, CNR1, CALB1, TAC1 were neuropeptides that play roles in the neurophysiology. Interestingly, all these four genes were packed into a PPI cluster, indicating that they interact with each other closely. It is worth noting that these four genes, expressed at the GABAergic neuron, can regulate the primary neuron circuit. Kenneth et al. have classified GABAergic neurons into 10 major clusters by single-cell transcriptomics, in which these genes were co-expressed in one of these clusters. CCK, CALB1 and CNR1 were mainly expressed in CCK interneurons while TAC1 was in PV interneurons. Both CCK and PV interneurons shared the same fast-spiking pattern [54]. AD animal models showed decreased activity of GABAergic interneurons because Aβ decrease the activity of GABAergic via cholinergic system [55] [56].

Based on these reports and our findings, we argue that the downregulation of CCK, TAC1, CALB1 is associated with cognitive impairment in AD patients.

Conclusion

CCK, TAC1 and CALB1 is associated with cognition in the hippocampus of AD patients. These genes may serve as potential therapeutic targets for cognitive dysfunction in AD.