Posterior cingulate cortex reveals an expression profile of resilience in cognitively intact elders

Abstract The posterior cingulate cortex, a key hub of the default mode network, underlies autobiographical memory retrieval and displays hypometabolic changes early in Alzheimer disease. To obtain an unbiased understanding of the molecular pathobiology of the aged posterior cingulate cortex, we performed RNA sequencing (RNA-seq) on tissue obtained from 26 participants of the Rush Religious Orders Study (11 males/15 females; aged 76–96 years) with a pre-mortem clinical diagnosis of no cognitive impairment and post-mortem neurofibrillary tangle Braak Stages I/II, III, and IV. Transcriptomic data were gathered using next-generation sequencing of RNA extracted from posterior cingulate cortex generating an average of 60 million paired reads per subject. Normalized expression of RNA-seq data was calculated using a global gene annotation and a microRNA profile. Differential expression (DESeq2, edgeR) using Braak staging as the comparison structure isolated genes for dimensional scaling, associative network building and functional clustering. Curated genes were correlated with the Mini-Mental State Examination and semantic, working and episodic memory, visuospatial ability, and a composite Global Cognitive Score. Regulatory mechanisms were determined by co-expression networks with microRNAs and an overlap of transcription factor binding sites. Analysis revealed 750 genes and 12 microRNAs significantly differentially expressed between Braak Stages I/II and III/IV and an associated six groups of transcription factor binding sites. Inputting significantly different gene/network data into a functional annotation clustering model revealed elevated presynaptic, postsynaptic and ATP-related expression in Braak Stages III and IV compared with Stages I/II, suggesting these pathways are integral for cognitive resilience seen in unimpaired elderly subjects. Principal component analysis and Kruskal–Wallis testing did not associate Braak stage with cognitive function. However, Spearman correlations between genes and cognitive test scores followed by network analysis revealed upregulation of classes of synaptic genes positively associated with performance on the visuospatial perceptual orientation domain. Upregulation of key synaptic genes suggests a role for these transcripts and associated synaptic pathways in cognitive resilience seen in elders despite Alzheimer disease pathology and dementia.


Introduction
Alzheimer disease is a major public health issue resulting in significant societal and economic burden. 1 Alzheimer disease is considered a spectrum disorder, 2-4 characterized clinically with declining memory, executive function and an inability to perform activities of daily living. 5,6 Neuropathologically, it is characterized by neurofibrillary tangles (NFTs) containing hyperphosphorylated tau, insoluble amyloid plaques, increased production of amyloid-beta peptide (Aβ) species, neuroinflammation and synaptic loss. [7][8][9] Although NFTs are associated with both Alzheimer disease progression and cognitive decline, 10-12 they are not absolute predictors of dementia. At least 15% of adults display NFTs in the medial temporal lobe (MTL) memory circuit characterized as Braak Stage I-II showing an age-associated elevation of tau pathology in cross-sectional health populations. 13,14 Interestingly, elders with a pre-mortem clinical diagnosis of no cognitive impairment (NCI) met criteria for Braak NFT stages ranging from I-VI [15][16][17] suggesting NFT pathology is not necessary for cognitive impairment. Identifying the molecular pathogenesis underlying brain resilience to cognitive decline despite varying stages of NFT pathology will provide new avenues for intervention to delay the onset of Alzheimer disease, an unmet need and a priority for the National Institute on Aging (NIA). 18 Although the MTL is an early site for NFTs, [19][20][21] the posterior cingulate cortex (PCC), a hub of the cortical default mode network 22 (DMN, Fig 1A), that plays a role in autobiographical memory retrieval, attention, salience and emotional context, 23,24 displays metabolic dysregulation during the onset of Alzheimer disease. [25][26][27][28] Neuroimaging studies indicate the DMN monitors the external and/or internal environment. [29][30][31] The PCC is dysregulated at resting state and during attention-demanding tasks in individuals with mild cognitive impairment (MCI) and Alzheimer disease. 32,33 Unlike other DMN hubs (e.g. precuneus, prefrontal cortex), 34 there are no standalone clinical molecular transcriptomic studies of the PCC from elderly people with a pre-mortem clinical diagnosis of NCI and a post-mortem Braak stage of I-IV, which may include a population resilient to the pathogenesis of Alzheimer disease. The lack of PCC transcriptomic information in elders with NCI, but with NFT pathology, impedes discovery science for therapeutics and understanding mechanisms underlying cognitive reserve/resilience that is not possible to model in preclinical animal and cellular preparations.
We performed high-throughput RNA sequencing (RNA-Seq), with subsequent specialized bioinformatic inquiry to assess genes and microRNAs (miRNAs) in association with clinical pathological variables using post-mortem PCC tissue obtained from elderly subjects that came to autopsy with a pre-mortem clinical diagnosis of NCI and received a post-mortem neuropathological Braak score of I-IV from the Rush Religious Order Study (RROS). 35,36 The goal was to identify a transcriptomic baseline profile of the PCC in healthy aged individuals without cognitive impairment but with varying stages of NFT pathology to generate a putative molecular fingerprint of resilience within this key hub of the DMN.

Materials and methods
The study cohort (n = 26) consisted of retired clergy with no signs of dementia at enrolment in the RROS, a longitudinal clinical pathological study. 35,37 Cognitive testing was performed annually during life. Post-mortem brains were examined for neuropathologic features of Alzheimer disease and related disorders. 38 Exclusion criteria included Lewy body dementia, Parkinson disease, hippocampal sclerosis, vascular disease and large strokes. 35,37,39,40 Apolipoprotein E (APOE) genotyping was performed as previously reported 35,37,39,40 and confirmed by RNA-Seq to identify nonsynonymous polymorphisms encoding base substitutions at amino acid positions 112 and 158. 41

Clinical and neuropathological evaluations
Briefly, RROS testing included the Mini-Mental State Examination (MMSE) 42 and a global cognitive score (GCS) compiled from a battery of 19 cognitive tests, which contribute to a cognitive domain score. 35,37 Neuropathological diagnosis was based on Braak NFT staging, NIA-Reagan criteria and the Consortium to Establish a Registry for Alzheimer's disease (CERAD). [43][44][45] In addition, brain slabs containing the PCC were immersion fixed in 4% paraformaldehyde, cryoprotected, cut into 40 µm thick sections and two sections from each case were immunostained with an antibody against the amyloid precursor protein (APP) and Aβ (6E10, 1:400 dilution) and tau (AT8, 1:250 dilution) as previously reported. 46,47 PCC 6E10 and AT8 loads were determined using a semi-quantitative score ranging from no 6E10-positive amyloid plaques and no AT8-positive NFTs, neurites or neuropil threads (0) to mild-to-moderate (2)(3) to moderate-to-severe (4)(5).
Preparation of tissue and RNA-Seq PCC was excised using fiduciary landmarks 48,49 and stored at -80°C until processing at the Collaborative Sequencing Center (Translational Genomics Research Institute, Phoenix, AZ). Total RNA from frozen slabs was extracted (mirVana; Ambion, TX) with enrichment for small RNAs, enabling assessment of mRNAs and non-coding RNAs (ncRNAs) including miRNAs. 50,51 Tapestation (DV200; Agilent, Santa Clara, CA) values ranged from 67.12% to 91.58%. RNA-Seq libraries were prepared using 500 ng of total RNA (TruSeq Stranded RNA Kit; Illumina, CA), ligated with xGen Dual-index UMI adapters (Integrated DNA Technologies, Coralville, IA) and enriched using eight PCR cycles. Libraries were paired-end sequenced (HiSeq4000, Illumina) for 80 base-pair (bp) reads.

Read processing
FastQ files were merged for paired ends before quality filtering and trimming using Fast Read Adjustment of Short reads (FLASH-1.2.11, minimum overlap 10 bases, maximum overlap 80 bases, mismatch allowed 1 in 4). 52 Reads were trimmed (sliding window of 3 bases with an average quality ≥ 32), quality filtered (average trimmed read quality ≥ 30) and size-selected (≥ 50 bases) using Trimmomatic (0.32) 53 resulting in three files per subject converted to fasta: single reads consisting of merged paired-end and R1 or R2 reads without a pair, R1 reads with a pair and R2 reads with a pair (see Supplementary Methods for details). The latter were collapsed into one paired reads file. The resulting two files (paired and unpaired) were mapped to Homo sapiens genome Genome Reference Consortium Human Build 38 patch release 13 (GRCh38.p13, hg38, assembly GCF_000001405.39), retrieved March 2020 (chr1-24, M), in Geneious using a custom annotation-span preference algorithm (v.9.0.1; Biomatters, Inc., CA). This involved a 13-mer index length (reads) and 18-mer word length (genome) and allowed for paired overlaps and gaps in reads as well as intron spanning. The hg38 genome was annotated using feature files for NCBI RefSeq, miRBase, LINC, and SNORD/ miRNA. After mapping to somatic chromosomes 1-22 and X (NC_000001-NC_000023) and mitochondrion (NC_00 0025), unused reads were mapped to the Y chromosome (NC_000024) with no chromosome masking. 54,55 Alignment files were exported and raw counts calculated using StringTie (2.1.1) 56 and the hg38RefSeq gtf attained from UCSC genome table browser August 2020 using default inclusion for All Tracks. Since exon information was not used in generating counts for differential expression, pre-mRNA was not differentiated from spliced mRNA.

Differential expression analysis
StringTie counts were imputed into EdgeR and DESeq2 using three comparison structures: Braak Stage I versus II versus III versus IV (six comparisons); Braak Stage I/II versus III versus IV (three comparisons); and Braak stage I/II versus III/IV (single comparison). Since not all entities were protein encoded, we use 'gene' to refer to both coding and noncoding annotations. A separate miRNA differential expression analysis used a custom reference gtf that included entries from mirBase and RefSeq. 57

Functional annotation clustering and gene enrichment
Each gene list was converted to Gene IDs inputted into Database for Annotation, Visualization and Integrated Discovery (DAVID, version 6.8, release October 2016) 58,59 and processed for annotation clustering (conducted January 2021) using multiple RNA and protein databases with a targeted focus on structure, function and gene ontology (Supplementary Table 1). This software generates an EASE score (one-tail Fisher's exact probability value), P-value 60 and FDR-corrected value 61 for each gene and database link within a group and an overall enrichment score for each grouping based on EASE scores. 62 We used enrichment scores above 1.00 based on the volume of output. Resources used to define gene product interactions and cellular compartment localization included Protein ANalysis THrough Evolutionary Relationships classification system (PANTHER) 63 and SynGO. 64 Protein names are derived from UniProtKB retrieved March 2021.

miRNA and transcription factor binding site databases
A combination of TarBase v7.0, miRBase and TargetScan databases (retrieved August 2020) generated 2,319 miRNA gene features. 57,65,66 Annotations from RefSeq, miRbase, and TarBase were crossed and used for downstream analysis. Genes regulated by miRNAs were determined using curated chromatin immunoprecipitation (ChIP)-Seq and experimental data for nucleic acid interactions, 65,66 and miRNA pathway analysis using the union of genes was performed using DIANA-mirPath. 67 Significant miRNAs as determined by differential expression analysis for protein-coding genes Bar graphs depicting diffuse Aβ (C) and neuritic plaque (D) regional distribution varied across Braak stages. At least 50% of Stage I/II and III cases displayed no or low plaque load, whereas Stage IV varied from absent to severe diffuse and neuritic plaques across the brain regions examined.
were excluded from downstream miRNA analyses to avoid cherry-picking data. Transcription factor binding sites (TFBSs) were examined using a combination of CpG islands <300 nucleotides (nt), transcription factor binding clusters observed with ChIP-Seq from the Encyclopaedia of DNA Elements (ENCODE; source data version 2018), and locations with histone 3 acetylated at lysine 27 (H3K27Ac), a marker of active regulatory sites. [68][69][70] Statistical analysis and data visualization Statistical tests were performed in R (version 4.0.4) and in Excel using custom-designed spreadsheets and scripts. Statistical significance was set at P < 0.05 and a false discovery rate (FDR) correction was applied where indicated. 60 For miRNA differential expression analysis, significance was set at FDR P < 0.10. Kruskal-Wallis and χ 2 tests were used for descriptive statistics and analysis of subject information. Gene counts are presented as counts per million (CPM) based on raw reads normalized within the DESeq2 analysis 71 adjusted for total reads, and as transcripts per million (TPM) based on a reference-guided assembly step in StringTie. 56 Also reported is the percentage of subjects in which a gene was expressed (PE) as found through assembly. Gravity network plots were made in Gephi using a two-step gravity loop that applies a separate algorithm at each step with iterations until a limit cycle or steady state is reached. Separate analysis confirmed uniformity with agglomerative hierarchical clustering using traditional distance metrics (hclust, R).

Data availability
Data are available upon request to the corresponding author and stored in GEO Database. Files document multiple steps in the process to aid further research and cross-validation efforts.

Differential expression data
Differential expression analysis with edgeR software revealed 644 differentially expressed (DE) genes that were significantly altered in either Braak Stage III or IV or when Stages III and IV were combined (Stage III/IV) and compared with Stage I/II. This contrasts with 750 DE genes found using DESeq2. Combining lists, a total of 917 genes were significantly altered in more advanced Braak stages, meaning 477 genes were shared between the two analyses (Fig. 2). None of the shared genes differed in direction of change. DESeq2 was selected for downstream analysis based on normalized counts, FDR correction and total output.
Comparison of DE genes that met inclusion criteria for expression level showed that 750 were statistically significant including 489 (65%) downregulated in advanced Braak stages (Fig. 3). In Stage III, 140, Stage IV, 195, and Stage III/IV, 383 genes were expressed at decreased levels compared with Stage I/II. Of these, 26 genes were shared across all three comparisons (Fig. 2B). Finally, 180 genes were decreased only when Stages III and IV were combined, and no differences were observed between Braak Stage III and IV, supporting the variance across these stages (Fig. 3).
probable n = 2 probable n = 6 possible n = 2 possible n = 1 possible n = 0 No Alzheimer disease n = 4 No Alzheimer disease n = 5 No Alzheimer disease n = 2 NIA-Reagan e Intermediate n = 1 Intermediate n = 2 Intermediate n = 8 P < 0.01 (χ) low n = 7 low n = 6 low n = 2 PCC 6E10 load f 2.6 (n = 8) 2.8 (n = 6) 4.6 (n = 10) Braak staging was deteremined using Bielchowsky silver stain and AT8 immunostaining to identify neurofibrillary tangle (NFT) severity and distribution across the brain. Braak Stages I and II display mild-to-moderate NFTs primarily in the entorhinal cortex; Stages III and IV display a larger involvement into limbic regions including the hippocampus; and stages V and VI revealed moderate-to-severe NFTs across brain regions.

Functional annotation clustering and pathway analysis
Downregulation of structure-related transcripts The majority of DE genes in Braak Stages III and IV were downregulated (Fig. 3). Functional annotation clustering based on protein structure, pathway interactions, shared function and gene ontology revealed a predominance of extracellular matrix (ECM) associated genes. Gene subcategories included coding products involved in basement membrane structure and ECM-cell interactions (e.g. adhesion, signalling and cell-cell junctions). Downregulation was seen for genes encoding classes: membrane proteins, vasculature, and protein metabolism. The latter category involved transcripts for peptidases, collagen digestion/absorption and protease inhibitors (Figs. 4 and 5).

Synapse gene upregulation
An elevated metabolic profile was found in Braak Stage III and IV compared with Stage I/II. Pre-and postsynaptic channel proteins involved in signal propagation, neurotransmitter release and signal summation were upregulated in advanced Braak stages. The 261 upregulated genes, of which 85 were associated with presynaptic signalling and chemical neurotransmission, were assigned to four functional categories: voltage-gated potassium channel and ion transport at the presynaptic compartment including a link with epilepsy ( (Table 2). Genes associated with the post-synapse were also affected, with 50 genes upregulated in Braak Stages III and IV compared with I/II divided into two functional categories: postsynaptic membrane (44 genes, enrichment score 3.23) and voltage-gated sodium channel activity/postsynaptic membrane depolarization (12 genes, enrichment score 2.09). Five genes were shared across groups: calcium voltage-gated   (7 genes, enrichment score 2.52); domain WW (6 genes, enrichment score 2.52); RNA polymerase II TF activity and sequence-specific DNA binding transcription factor forkhead box (FOX) (15 genes, enrichment score 1.90); and positive regulation of transcription from RNA polymerase II promoter (97 genes, enrichment score 1.64; 307 genes, enrichment score 1.54; 134 genes, enrichment score 1.22); and regulation of transcription from RNA polymerase II promoter and negative regulation of protein metabolic process (seven genes, enrichment score 1.40). A total of 17 genes were seen in five of the seven clusters, including two involved in chromatin modelling, 13 involved in gene-specific transcription regulation and 2 coded for cell structure products.

Dimensionality reduction highlighted gene upregulation
PCA of the 750 DE genes explored covariance within individuals. Dimension 1 accounted for 46.2% of the variance and Dimension 2 10.4%. After running the regression calculation using DESeq2 normalized gene expression values, we examined subject factors not used in deriving the PCA results. This process collapses subjects within categorical groupings (e.g., Braak stage) to derive a theoretical variable location and confidence interval, presented as coordinates and ellipsis on a Dimension 1 × Dimension 2 biplots. Although Braak Stage I/II segregated from Stages III and IV, there was no difference between Braak Stages III and IV ( Supplementary Fig. 3). Overlay of male/female (categorical), APOE ϵ status (categorical) and age at death (vector) did not show differences across categories or influence of age on a biplot, or very near the origin, indicating that age lies in a different dimension (Supplementary Fig. 3).
Although 65% of the DE genes were downregulated in Stages III and IV, PCA using the same gene list highlighted upregulated genes as the largest contributor to variance across subject gene expression profiles. Taking the top 10% of contributors to Dimension 1 (75 genes accounting for 22.7% of Dimension 1, and 10.5% of total variance), 69 genes (92%) were upregulated and six genes were downregulated in Braak Stages III and IV compared with Stage I/II, a stark difference from the 35% percent of total DE genes that were upregulated in the more advanced Braak stages. Functional clustering highlighted neuronal cation channel activity (48 genes, enrichment score 2.88), synaptic signalling (10 genes, enrichment score 2.06), and postsynaptic No significant correlations were found for the subject information: years of education, mini-mental state examination, global cognitive score, episodic memory and semantic memory. † FDR P < 0.10; ‡ FDR P < 0.05; * P < 0.05; ** P < 0.01; ns, not significant.
membrane (15 genes, enrichment score 2.00) as pathways and physiological mechanisms enriched in Dimension 1 upregulated genes, whereas phosphoprotein binding (six genes, enrichment score 1.08) and transcription regulation (four genes, enrichment score 1.02) were enriched in the Dimension 1 downregulated genes. PCA revealed neither a contribution by cognitive domain and performance scores nor highlight a difference between Braak stage I/II and Braak stage III or IV. Dimension 1 on the cognitive PCA contributed 27.7%, and Dimension 2 13.4%, which is closer to a random distribution (based on a run of 10 PCA with values from cognitive data replaced with random numbers, average Dimension 1 = 12.5%, Dimension 2 = 11.2%, regression slope 0.64) than to the PCA with genes. There was a difference between male and female theoretical variable overlays on the cognitive PCA ( Supplementary Fig. 3C, D). As with the gene expression PCA, age and APOE ϵ status did not show an influence on Dimension 1 × Dimension 2 biplots.

Synapse-related functional pathways associate with cognitive performance
Significant DE genes correlated with cognitive test scores (Table 3) following FDR correction for multiple comparisons. A cut-off of rho ≥ |0.55| with an uncorrected P-value <0.005 was used to determine associations. Gene expression was not associated with composite cognitive scores for episodic, working, or semantic memory (Table 3). Less than 1% of DE genes (< 8 genes) correlated with each episodic memory test: delayed logical memory II, word list and word list recall; working memory test: alpha span; and semantic memory tests: category fluency and reading test. Ten genes positively and one negatively associated with performance on the Boston naming test of semantic memory. Of note, the 10-item reading test was associated with 5 genes, including neurotrophic receptor tyrosine kinase 1 (NTRK1; CPM <2, TPM <2, PE = 100%), which encodes TrkA, the cognate receptor for nerve growth factor (NGF). 72-74

Figure 6 Association networks showing relationship of regulatory mechanisms and DE genes in Braak Stages I/II compared
with Stages III/IV within the PCC of elderly adults with NCI. Functional annotation clustering of 750 DE genes was performed using a TFBS annotation file that combines information on chromatin structure and ChIP (see Methods) to derive a list of genes associated with a given TFBS. Clustering first matches Braak-stage DE genes to respective associated TFBSs and groupings of TFBS are selected using a calculated enrichment score (based on number of DE genes) to determine significance. Following discovery of significant TFBS clusters, we used a network-based map to illustrate associations. (A) Three TFBS clusters were associated each with DE genes with expression significantly upregulated and downregulated in Braak Stage III/IV compared with I/II. Of note, the direction of expression change refers exclusively to DE genes, and not TFBS factors. This plot is spatially agnostic and no information can be derived from axes; the layout is the consequence of a force-directed algorithm and conveys information only in distance (farther = looser association), not in position relative to any constant (like an axis or grid). Each green dot (node) represents a TFBS cluster (e.g., LHX3/CART/FREAC7); blue, a specific expressed gene (mRNA) with higher (bottom right) or lower (top left) levels in Braak Stages III/IV compared with Braak Stages I/II; and, pink, a specific microRNA. Every line (edge) represents an association as determined from the TFBS database outlined previously, coloured according to component nodes with no information delivered via edge thickness. Node size is based on number of associations but should be considered minimally informative at this resolution. All associated mRNA can be found in Supplementary Fig. 2A. (B) DE microRNA (miRNA), as detected in a separate differential expression analysis, were analyzed for associated genes through a literature and multiple database search (see Methods). This compiled list of associated genes was then crossed with the DE genes. Interestingly, many mRNA nodes have multiple associations with regulatory TFBS and miRNA. The direction of change for miRNA had no consistent association with direction of change in mRNA ( Supplementary Fig. 2B, C) and many miRNA were associated with up-and downregulated mRNA seemingly indiscriminately. Names of genes associated with miRNA can be found in the Supplementary Material.   compared with I/II. Functional annotation clustering based on gene structure, function and gene ontological category using genes positively correlated with the composite perceptual domain score showed enrichment in transcript classes encoding proteins associated with axon activity and postsynaptic membrane potential (Supplementary Tables 5 and 6).

Discussion
We found 489 downregulated and 261 upregulated genes in PCC obtained from elderly subjects that died with a premortem clinical diagnosis of NCI and post-mortem pathological evaluation of Braak Stage I, II, III and IV. Despite predominantly downregulation across Braak stages, upregulation of individual expression profiles was most prevalent in Stage III and IV compared with I/II. Dimension reduction analysis found that upregulated genes primarily contributed to Dimension 1, which accounted for nearly half of the covariance across individuals. Of the top 10% Dimension 1 genes, enrichment was primarily related to excitatory synaptic transmission, which correlated strongly with cognitive performance. Dimension 2, the next highest orthogonal contributor to individual covariance, revealed a decrease of neuromodulatory genes in later Braak stages with differences between Braak Stages III and IV. These novel findings emphasize the profound changes in synaptic and neuromodulatory genes that may underlie a mechanism of resiliency and cognitive reserve in the face of mounting Alzheimer disease pathology with NCI. Commensurate with our post-mortem human brain findings, animal models of aging have been integral in the development of a compendium of possible candidates for cognitive reserve genes (CRGs). [76][77][78] Further, independent studies that support the current results found gene expression alterations between Braak stage I/II compared to III that were related to synaptic plasticity, mitochondrial function, GPCR signalling, electron transport and calcium ion binding, among others in the prefrontal cortex (PFC) DMN hub. 79,80 Upregulation of genes encoding synaptic transmission and cellular energy metabolism observed in the more advanced Braak cases is analogous to increased frontal lobe neuroactivity reported in older adults without cognitive impairment measured by PET imaging. 81,82 These findings suggest that these alterations are involved in the compensatory preservation of cognition despite the increase in neuropathology. Over time, these initial cognitive resilience mechanisms to maintain function may ultimately fail to preserve cognition with advancing age or are overwhelmed by the onslaught of disease pathology. 83 For example, the present findings suggest resiliency at the metabolic level may fail in those with cognitive decline similar to that seen in model organisms. 83 In addition to aging and pathology, resilience likely is influenced by sex, life experiences, education, connectional plasticity, and epigenetics. [83][84][85][86] Whether upregulation in cellular activity genes underlying metabolic dysregulation and altered connectivity patterns found in the PCC across Braak stages is similar or different between hubs of the DMN remains to be determined. 32,33,87,88 Therefore, uncovering the mechanism(s) for increased cortical synaptic activity will have clinical and quality-of-life implications for the elderly and enhance putative therapeutic implications using previously reported novel CRGs. 78,89 There are no PCC transcriptomic datasets in elders across the Alzheimer disease spectrum that offer a tool for comparison, highlighting the importance and novelty of the present findings. Analogous PFC gene expression in pre-middle-aged (≤ 40 years) compared to aged non-demented adults (≥ 70 years) found decreases in genes associated with inhibitory neurotransmission and neuropeptide systems. 90 Interestingly, while GABA marker gamma-aminobutyric acid type A receptor subunit gamma2 (GABRG2) and glutamate marker G protein-coupled receptor 158 (GPR158) expression were increased in the PCC of Braak Stage III/IV compared with I/II, a significant downregulation occurred in the PFC of aged compared with pre-middle-aged adults. 90 Although this may represent regional DMN profile differences, a pathology × age interaction could relate to cohort composition or size. GABRG2 and GPR158 expression levels in the brain 91,92 is linked to aging 93,94 and adult neuropsychiatric conditions, [95][96][97] and GPR158 expression is associated with Alzheimer disease pathology as well as frontotemporal dementia. 98,99 Moreover, GPR158 downregulation is related to hippocampal-mediated cognitive deficits. 93,100,101 Interestingly, glutamatergic presynaptic markers increase in MCI cortex, suggesting a paradoxical inhibitory response to dementia onset. 102 Of the pre-and postsynaptic protein-encoding genes upregulated in Braak Stages III and IV compared with I/II, GRIN2A mRNA is also elevated in the hippocampus in MCI compared with NCI 103 suggesting a target for intervention. 104 Microarray studies also reveal VAMP1 mRNA elevation in the superior frontal gyrus and increased hippocampal STXBP5L mRNA in MCI compared to NCI, while both are decreased in entorhinal cortex, 103 suggesting differential brain regional vulnerability between aging and the onset and progression of Alzheimer disease. A negative association of VAMP1 expression and Braak stage was observed when analysis included Braak Stage V and VI. 105 Upregulation of these genes was found in PCC of Braak Stage III/IV compared with I/II, with no advanced stages for comparison. While altered exocytotic vesicle transcripts along with VAMP1 and STCBP5L occur in the hippocampus and PFC in MCI compared with NCI, 103 similar findings were not seen in our study. Notably, studies using lower organisms report opposing directional changes in transcripts and proteins in response to pathological mutations. 106 We found similar decreases to those reported in MCI compared with NCI including decrements in neocortical expression for ITGB1 and ITGB8. 103 Therefore, ITGB1 may play a role in the progression of Alzheimer disease through alterations in oxidative stress.
Increased expression of postsynaptic genes reveals elevated synapse activity and a decrease in neuromodulatory genes in more advanced Braak stages. Genes encoding vesicular transporters for dopamine, DAT (SLC6A3), and norepinephrine, NET (SLC6A2), involved in the synaptic reuptake of catecholamine neurotransmitters, as well as choline acetyltransferase (ChAT), the synthetic enzyme for acetylcholine were significantly decreased in PCC in Braak Stages III/IV compared with I/II. Although decrements in ChAT activity have been reported in the PCC in Alzheimer disease, ChAT expression remains stable in MCI. 107 Interestingly, we found upregulated excitatory gene profiles even within functional clusters defined by neuromodulatory circuits. For example, we found an increase in VGLUT1 (SLC17A6), a presynaptic transcript that encodes a protein involved in primary excitatory transmitter release, and a decrease in the transcript that encodes a transporter involved in Glu synthesis xCT (SLC7A11) in functional clusters associated with the neurotransmitters, dopamine, noradrenaline, and serotonin. These findings suggest a molecular signature of decreased neuromodulatory activity and elevated excitatory neurotransmission. Examining these changes in light of neuropathology, changes in genes encoding DAT, NET, and ChAT occurred in Braak stage IV, whereas excitatory transmitter changes were seen in Stages III or III/IV. This provides a possible timeline for resilience through molecular mechanisms whereby neuromodulation is altered in response to elevated excitatory neurotransmission. Since the PCC receives neuromodulatory innervation from spatially distinct cell populations, this raises the possibility of a diffuse connectome reorganization in NCI elders with increased NFT pathology. These observations may demonstrate neuroplasticity associated with resilience that may play a role in the ability to perform age-related task completion strategies. 81,82 We found significant differences in miRNAs only in later Braak stages (e.g. IV or III/IV compared with Stage I/II) and no detectable differences at Stage III compared with Stage I/ II. These findings support previous studies suggesting miRNA alterations occur later than alterations in genes they regulate in individuals with MCI compared to NCI. 108 Possible factors contributing to these temporal differences include a secondary regulatory response to disease onset or an inability to regulate homeostasis by post-translational modifications, 109 an idea supported by our finding Braak Stage III changes in gene pathways involved in transcription regulation. Further, Braak Stages III and IV show a marked upregulation in transcripts encoding kinases, downregulation in phosphatases, and an increase in ubiquitin protein-encoding, pathways similarly implicated in Alzheimer disease. 110,111 Alterations of protein metabolic factors also occur in MCI compared with aged controls 103 and are associated with NFTs, 112 differentiating these changes from normal aging. 113 The changes in miRNAs associated with Braak Stage IV indicate an expression imbalance in response to pathogenesis and may provide a viable target for identifying resilience or lack thereof across the Alzheimer disease spectrum. 114 The specific miRNA alterations reported here have not been previously identified across the Alzheimer disease spectrum; however, there is poor consensus and systematization for the evaluation at non-coding regulatory elements, making comparison tenuous. 115,116 When we crossed the miRNA list with known gene interactions, we found DE genes more highly expressed in Braak Stages III and IV were those involved in synaptic activity, especially with regards to regulatory elements hsa-miR-8061 and hsa-miR-548a-3p, miRNA decreased in Braak Stage IV, and hsa-miR-5692b and hsa-miR-134-5p, miRNA increased in Braak Stage IV compared with Stage I/ II. Moreover, comparison of miRNA expression with cognitive function revealed high association with the same visuospatial domains associated with DE genes found in a separate analysis. The precise molecular pathogenic role that miRNAs play during the progression of Alzheimer dementia remains to be defined. As exploration into CRGs continues, these regulatory mechanisms may prove insightful for defining a timeline and therapeutic targets for the treatment of cognitive decline in the elderly and those with dementia. 78 PCC expression profiling revealed a significant association with NFTs but not amyloid or neuritic plaque pathology. Differential expression analysis using CERAD or NIA-Reagan neuropathological scores as grouping factors was indistinguishable from random grouping. This corresponds with prior investigation of individuals with MCI or Alzheimer disease that demonstrated minimal correlation between parenchymal plaque pathology and cognitive impairment. 17,117 Similarly, ApoE allele as a grouping factor was indistinguishable from a random grouping factor on differential expression analysis and PCAs with either gene expression or cognitive performance. However, a study of ApoE status and brain glucose metabolism in non-demented adults aged 30-95-years-old found an age-related significant decline with greater uptake in ϵ4 noncarriers compared with carriers in DMN hubs including the PCC. 118 Moreover, in participants older than 70 years, there was no interaction between Pittsburgh Compound B amyloid binding status and APOE ϵ4 genotype with respect to glucose metabolism. 118 These findings indicate the PCC has a unique vulnerability to reductions in glucose metabolic rate as a function both of age and APOE allele status, perhaps due to its role as a hub of the DMN that deactivates when mental effort is required but is less efficient in deactivation during the progression of Alzheimer disease. 118 Since ApoE genotype represents a life-long state, persons with a higher level of education or a lifestyle that involves frequent cognitive engagement may be less likely to have detectable differences on cognitive tests that correlate with ApoE allele status. A potential limitation in the present study is that the small number of ϵ4 carriers may mask PCC genotype changes associated with ApoE ϵ status. It is possible that educational level affects the expression of various classes of transcripts including the upregulation of synaptic genes. Level of education has been suggested to play an important role in preventing the onset of dementia through brain reserve. 119 Interestingly, the higher Braak stage group had an average education level 4 years greater than the lower Braak cases, suggesting the intriguing concept that educational level plays an active role in the upregulation of synaptic transcripts found in the current high Braak cases. 119 More detailed investigations of the interaction between education, Braak stage, brain resilience, and gene expression are warranted.
Studies indicate TOMM40 variants are associated with estimating onset of Alzheimer disease and interaction with ApoE status can increase disease onset, which may be geographically dependent. 120 We found that TOMM40 expression was significantly increased in Braak stage III compared to I/II with no difference in transcript variants based on reference-guided assembly. Previously, blood analysis revealed a significant association between longer TOMM40 poly-T lengths and neuroimaging higher medial temporal cortex plaque and NFT burden in non-demented older adults. 121 TOMM40 '523 polymorphism affects expression levels of APOE, and TOMM40 mRNAs in the temporal and occipital cortices of late-onset Alzheimer disease and non-demented controls. 122 The molecular and biochemical mechanism(s) underlying the effect of increased TOMM40 expression upon Alzheimer disease pathophysiology remains to be investigated. However, structural DNA variations, especially those in intronic or intergenic regions such as TOMM40 '523, may alter gene transcription efficiency, timing of transcription, transcript stability, transcript splicing and/or epigenomic modifications. 123 While we have studied transcript variants from reference-based assembly, we have not yet investigated polymorphisms. This is in progress for all DE genes and will help to clarify the possible role of TOMM40 in CRG-related processes. However, it is possible that TOMM40 is part of a resilience mechanism that is specific to a select group of variants and not the generalized elderly population.
We provide evidence for putative brain cognitive reserve as a mechanism for resiliency based upon differential molecular expression profiling of the PCC genes derived from elders with NCI but with different Braak scores. 114, 124 Although the present definition is similar to that established by the Collaboratory on Research Definition for Reserve and Resilience, it also is reminiscent of 'potential cognitive reserve genes', in which genes are selected depending upon whether they display differential expression 78 based upon Braak stage. In the present report, brain resilience and cognitive reserve suggest that a population older individuals have functional and structural physiological changes, such as increased synapse number or size, or adjusted cognitive strategies which allow the brain to tolerate a greater degree of pathology without suffering decline on cognitive tasks. 114,125 Along this line, resiliency/reserve may also involve recruitment of other brain areas resulting in increased cortical innervation from regions not severely affected to aid in task performance. Our findings suggest that cognitive reserve and resilience likely involves synaptic and metabolic pathway expression that increases across Braak Stages III and IV as a potential compensatory response to age-related cortical denervation. 126 In this regard, it has been proposed that reserve can be measured or inferred either through increased brain structural and/or physiological premorbid capacity. 127 Interestingly, a disconnect between the Alzheimer disease proteome and transcriptome in the PFC was reported, 128 suggesting the importance of investigating proteins in addition to their coding transcripts that likely play a role in brain resilience, especially within hubs of the DMN including the PFC and PCC. Interestingly, a mathematical assessment of the transcriptome from different aging studies found in relevant animal models one in ∼six age-related genes were considered poor behavioural predictors, highlighting expression variability and biological variance 78 that may be applicable to defining CRGs and exploiting them for therapeutic interventions.
Finally, it is important to consider study limitations. Tissue was obtained from a subpopulation of the RROS cohort with lifestyle elements that differ from a secular community-based cohort, 129,130 which likely affect the bidirectional relationship between cognitive stimulation and cognitive status. 125,[131][132][133] Since we examined individuals who aged into their 9th decade without cognitive impairment, natural limitations affect cohort size and applications of computational detection allowing for clustering into expression between successful agers versus those progressing to MCI. 134 However, a strength of this population is homogeneity and low rate of subject attrition over time. Importantly, regional brain dissections consist of an admixture of different cell types resulting in an expression profile that masks changes in specific cells at the sequencing and computational stage. Notwithstanding these caveats, we uncovered mRNAs in human PCC that were differentially expressed between Braak Stages I/II and III/IV in addition to associated miRNAs and TFBSs. Inputting significantly different gene/ network data into a functional annotation clustering model revealed elevated presynaptic, postsynaptic and ATP-related expression in Braak Stages III and IV compared with Stages I/II, suggesting these pathways are integral for cognitive resilience seen in elderly non-demented cases. Braak stage was not associated with cognitive function but upregulation of synaptic genes positively correlated with visuospatial perceptual orientation tasks. These findings suggest increased synaptic expression, in part, underlies cognitive resilience in elders despite Alzheimer disease pathology.

Funding
Supported by National Institute of Aging grants PO1AG014449, PO1AG017617, P30AG072931 and P30AG10161 and Arizona Alzheimer's Consortium at Barrow Neurological Institute.