Abstract
Vascular endothelial cells play an important role in maintaining brain health, but their contribution to Alzheimer's disease (AD) is obscured by limited understanding of the cellular heterogeneity in normal aged brain and in disease. To address this, we performed single nucleus RNAseq on tissue from 32 human AD and non-AD donors (19 female, 13 male) each with five cortical regions: entorhinal cortex, inferior temporal gyrus, prefrontal cortex, visual association cortex, and primary visual cortex. Analysis of 51,586 endothelial cells revealed unique gene expression patterns across the five regions in non-AD donors. Alzheimer's brain endothelial cells were characterized by upregulated protein folding genes and distinct transcriptomic differences in response to amyloid β plaques and cerebral amyloid angiopathy. This dataset demonstrates previously unrecognized regional heterogeneity in the endothelial cell transcriptome in both aged non-AD and AD brain.
SIGNIFICANCE STATEMENT In this work, we show that vascular endothelial cells collected from five different brain regions display surprising variability in gene expression. In the presence of Alzheimer's disease pathology, endothelial cell gene expression is dramatically altered with clear differences in regional and temporal changes. These findings help explain why certain brain regions appear to differ in susceptibility to disease-related vascular remodeling events that may impact blood flow.
Introduction
Cerebral vasculature plays a central role in the maintenance of brain health, and in Alzheimer's disease (AD), changes that impact blood flow and blood–brain barrier function may exacerbate neurodegeneration and cognitive decline (Iturria-Medina et al., 2016; Montagne et al., 2016; Benedictus et al., 2017; Leijenaar et al., 2017; Nation et al., 2019). A growing body of literature suggests that diverse blood vessel changes, including aberrant angiogenesis, vascular pruning, inflammation, senescence, and other remodeling events, occur alongside AD-related pathologic aggregation of amyloid β and tau (Brown and Thore, 2011; Bennett et al., 2018; Bryant et al., 2020; Lau et al., 2020). However, our current understanding is limited by the lack of characterization of the regional heterogeneity observed in the human brain, particularly in non-neuronal cell types, such as endothelial cells. Indeed, few human studies have deeply characterized vasculature from multiple brain regions, and comparisons with mouse models may be inadequate (Hodge et al., 2019). Thus, we have, at best, an incomplete picture of how these vascular mechanisms underlie functional alterations that affect disease progression.
Previous work has highlighted the heterogeneity of the endothelial transcriptome and translatome between the brain and other organs, such as lung, heart, liver, and kidney (Nolan et al., 2013; Jambusaria et al., 2020; Kalucka et al., 2020; Paik et al., 2020). Moreover, advances in single-cell sequencing technologies have enabled the resolution of brain endothelial cell subtypes (Vanlandewijck et al., 2018; Kalucka et al., 2020; Yang et al., 2022). However, such studies have only examined either 1 or 2 brain regions or a pooled cortical homogenate, which may obscure cerebrovascular heterogeneity at the region-specific level (Bernier et al., 2021). Other groups have recently demonstrated region-specific cellular subtypes in murine astrocytes (Batiuk et al., 2020; Hasel et al., 2021) and microglia (Tan et al., 2020). We reasoned that the endothelial cell transcriptome may also vary by location in the brain, given both developmental and morphologic differences in the microvasculature across brain regions.
Understanding regional heterogeneity in the human brain is important in AD because amyloid β, neurofibrillary tangles, and neurodegeneration do not affect all brain regions to the same extent. Numerous neuropathological studies have outlined the spatial distributions of these changes. Amyloid deposits in region-specific patterns are distinct for plaques and for cerebral amyloid angiography (CAA). The temporal and frontal association cortices accumulate more plaques than primary sensory areas, such as visual cortex (Haroutunian et al., 1998; Serrano-Pozo et al., 2011), whereas CAA has much higher prevalence in parietal and occipital cortices (Vinters and Gilbert, 1983). Neurofibrillary tangles show yet another distinct hierarchical regional pattern. Braak staging describes the progressive severity of tau neurofibrillary tangle pathology, with early stages affecting primarily the transentorhinal and entorhinal cortex (EC, stages I-II), then structures including the medial and inferior temporal cortex (ITG, stages III-IV), followed by spread to isocortical areas including the pre frontal cortex and association visual cortex regions (PFC, V2, stage V) before appearing in primary sensory areas (stage VI) including visual cortex (Braak and Braak, 1991; Braak et al., 2006). Notably, neuron loss correlates better with neurofibrillary tangle burden than with amyloid β plaque burden (Arriagada et al., 1992; Gómez-Isla et al., 1997). For this reason, we focused these studies on cortical regions included in the Braak staging scheme that represent the hierarchical accumulation of tau pathology in AD, to examine the effect of disease progression on the endothelial cell population.
Considering that cortical brain areas are differently affected by amyloid β and tau pathologies, and that vascular changes are intimately related to AD, the goals of this work were twofold: (1) to characterize baseline heterogeneity in the endothelial cell transcriptome in the normal aged brain; and (2) to elucidate how regional differences in disease progression impact endothelial cell gene expression in AD.
Materials and Methods
Donor tissue collection
Thirty-two donors were selected from the Neuropathology Core of the Massachusetts Alzheimer's Disease Research Center, which obtained written informed consent from next-of-kin for brain donation. Tissues were assessed by a neuropathologist for AD-related neuropathology burden according to the National Institute on Aging and Alzheimer's Association guidelines (Hyman et al., 2012) (Extended Data Table 1-1) and included high-pathology Braak V/VI donors (n = 16, 10 F/6 M; B3), intermediate-pathology Braak III/IV donors (n = 5, 3 F/2 M; B2), and low-pathology Braak 0-II donors (n = 11, 6 F/5 M; B0/1). For each donor, frozen postmortem cortical tissue was obtained from the following five brain regions where available: EC, inferior temporal gyrus (ITG), PFC, secondary visual cortex (V2), and primary visual cortex (V1). Formalin-fixed, paraffin embedded tissue blocks from the opposite hemisphere were also obtained for histologic studies.
Table 1-1
Clinical and neuropathological characteristics of the n = 32 donors included in this study. In the five rightmost columns, an X indicates that nuclei were sequenced for the given donor and region, while a blank space means tissue was unavailable for sequencing. Download Table 1-1, XLSX file.
Histology
For validation of targets identified by transcriptomic analysis, paraffin-embedded EC and V1/V2 were sliced on a microtome at 8 μm, and sections were mounted on Superfrost slides. Two sets of sections were prepared: one containing EC and ITG at the level of the lateral geniculate nucleus, and the other containing occipital cortex, including primary and association visual areas. Sections were rehydrated in xylene and a descending concentration of ethanol. Epitope retrieval was conducted in citrate buffer (pH 6.0, containing 0.05% Tween-20) with microwave heating at 95°C for 20 min. Sections were then blocked with 5% BSA in TBS containing 0.25% Triton-X for 1 h at room temperature and then incubated with primary antibodies overnight at 4°C. Primary antibodies used in this study targeted CREM (Invitrogen catalog #PIPA581971), HSPH1 (Proteintech, catalog #13383-1-AP), HSP90AA1 (Invitrogen catalog #PA3-013), IGFBP3 (Invitrogen catalog #PA5-27190), MT2A (Invitrogen catalog #PA5-102549), and GLUT1 (EMD Millipore catalog #07-1401). The following day, sections were rinsed 3 times with TBS and then incubated in Alexa dye-conjugated secondary antibodies in TBS containing 0.25% Triton-X for 1 h at room temperature. Afterward, sections were rinsed 3 times with TBS and Trueblack (Biotum) was applied for 30 s to quench autofluorescence. Sections were rinsed with H2O and then coverslipped with Fluoromount G with DAPI. An Olympus VS120 slide scanner was used to capture digital images using a 20× lens.
For measurement of plaque burden in tissues, cryosections (10 µm thickness) adjacent to those used for single-nucleus (sn) RNAseq were subjected to immunohistochemistry with mouse monoclonal anti N-terminal abeta antibody clone 3D6. The staining was performed on a Leica BOND Rx automated stainer using DAB-based detection (Leica). Sections were scanned in a slide scanner (3DHistech, Pannoramic 250) and area fraction (i.e., % area of tissue section occupied by 3D6-immunoreactive plaques) was measured using HALO software (Indica Labs).
For assessing CAA, tissue sections were processed for amyloid β (clone 6F/3D, Agilent catalog #M0872) DAB histology by the Massachusetts Alzheimer's Disease Research Center. Four tissue blocks were assessed which encompassed all five brain regions; V1 and V2 were located in the same block and received the same score. Slides were digitized and evaluated by a blind researcher (B.W.) using a 0-4 CAA grading scheme (Olichney et al., 1995). An average score for each donor was used for analysis.
ELISA
For quantification of phospho(T231)tau and total tau, 25 mg of tissue was dissected before nuclei isolation. Tissue was Dounce homogenized with 5 v/w of PBS and then centrifuged to pellet insoluble material at 3000 × g for 5 min at 4°C. The resulting supernatant was reserved for analysis using MesoScale Discovery Kit #K15121D according to the manufacturer's instructions.
Isolation of nuclei from postmortem frozen brain tissue
Nuclei isolation and sorting were performed as previously described with minor modifications (Gerrits et al., 2021). Briefly, fresh frozen tissue was cryosectioned (∼40 sections of 40 µm thickness each) and lysed in sucrose lysis buffer (10 mm Tris-HCl, pH 8.0; 320 mm sucrose; 5 mm CaCl2; 3 μm Mg(Ac)2; 0.1 mm EDTA; 1 mm DTT and 0.1% Triton X-100). Lysates were filtered through a 70 µm cell strainer. Nuclei were purified by ultracentrifugation (107,000 × g for 1.5 h at 4°C) through a sucrose cushion (10 mm Tris-HCl, pH 8.0; 1.8 m sucrose; 3 μm Mg(Ac)2; 0.1 mm EDTA and 1 mm DTT). Supernatants were removed and pellets were resuspended in 2% BSA/PBS containing RNase inhibitor (0.2 U/μl) (Roche). Nuclei were incubated with fluorescently conjugated antibodies against the neuronal marker NEUN (RBFOX3/NEUN (1B7) AF647 mouse mAB, Novus Biologicals, NBP1-92693AF647) and the oligodendrocyte transcription factor OLIG2 (Anti-OLIG2 clone 211F1.1 mouse mAb, Merck Millipore, MABN50A4). Samples were kept on ice throughout the isolation and staining procedure. Nuclei were stained with Sytox blue (Thermo Fisher Scientific) and sorted on a BD FACSAria Fusion. For each sample, we collected Sytox+NeuN+Olig2– and Sytox+NeuN–Olig2– nuclei. This second group, depleted for neurons and oligodendrocytes, was examined in this study.
snRNAseq library preparation
Single-nucleus cDNA libraries were constructed using the Chromium Single Cell 3′Reagents Kit version 3. Samples were pooled and sequenced targeting at least 30,000 reads per cell on an Illumina NovaSeq.
snRNAseq data preprocessing and quality control
Sytox+NeuN–Olig2– nuclei raw data were processed with 10× Genomics cellranger software version 4, based on GRCh38 pre-mRNA reference, and quality control was performed to identify potential sample outliers. Subsequently, raw counts were filtered per region to include only cells with >100 exonic counts, and >800 features and counts (both intronic and exonic) as well as a mitochondrial gene percentage <15%. A sample specific subsequent filter removed additional potential low-quality cells, only including cells higher/lower than the median ± 3*MAD (median absolute deviation) of gene and UMI count numbers per sample. This yielded final matrices of 231085 (ITG) + 202514 (EC) + 280665 (PFC) + 253403 (V1) + 234567 (V2) nuclei and 35 423 features.
UMI counts were normalized and integrated across donors for each brain region separately, following the recommended Seurat integrative workflow (Stuart et al., 2019). For each brain region, the filtered count matrices were log-normalized and the top 1000 highly variable genes were identified for each donor using FindVariableFeatures() with selection.method = “vst” Anchoring features for cross-sample integration were identified using FindIntegrationAnchors() with reduction = “rpca” for reciprocal principal component analysis (PCA), as is recommended for larger datasets (see https://satijalab.org/seurat/articles/integration_rpca.html). The total number of reads (nCount_RNA) and percent mitochondrial reads (percent.mito) were regressed out. PCA was applied for dimensionality reduction using RunPCA() with npcs = 50. The PCA-reduced dataset was further passed to uniform manifold approximation and projection (UMAP) for nonlinear dimensionality reduction using the RunUMAP() function, with dims = 1:30. Graph-based shared nearest neighbor clustering was applied using FindClusters() with resolutions 0.1 through 1. Subcluster resolution 0.1 was selected to identify clusters of the following cell types: astrocyte, microglia, neuron, endothelial cell, oligodendrocyte, or other. Cell types were manually identified using DimPlot and heatmap visualization of canonical cell type markers (for markers used, see Extended Data Table 1-2).
Table 1-2
Marker genes used to identify cell types. (Data shown in Figure 4). Download Table 1-2, XLSX file.
The Seurat objects for each brain region were subsetted to only the clusters identified as endothelial cells, and the endothelial subclusters were reintegrated using canonical correlation analysis within each brain region. Specifically, anchor features were identified based on the top 2000 highly variable genes and were passed to FindIntegrationAnchors() with k.filter = 20. The data were canonical correlation analysis-integrated using IntegrateData() with k.weight = 20 and scaled using ScaleData(), with the total number of reads (nCount_RNA) and percent mitochondrial reads (percent.mito) regressed out. PCA was applied for dimensionality reduction using RunPCA() with npcs = 50. The resulting scree plot was visualized to select the number of components to retain for graph-based clustering for each brain region. We chose the following numbers of components: EC, 8; ITG, 9; PFC, 7; V2, 8; V1, 9. UMAP dimensionality reduction and shared nearest neighbor clustering were performed as described above. The final region-wise endothelial cell matrices for each brain region had the following cell counts: EC, 4271; ITG, 10,850; PFC, 12,638; V2, 15,553; V1, 18,235.
These endothelial matrices were then reciprocal PCA-integrated across brain regions using FindIntegrationAnchors() with the top 2000 highly variable genes and IntegrateData() with k.weight = 20. The inter-region integrated dataset was scaled using ScaleData() with nCount_RNA and percent.mito as variables to regress out. Dimensionality reduction was performed using RunPCA() with npcs = 50 and and RunUMAP() with dims = 1:30. Graph-based shared nearest neighbor clustering was applied using FindClusters() with resolutions 0.1 through 1. Residual contamination from nonendothelial cell types was evaluated using the DoHeatmap() function with known cell type markers (Extended Data Table 1-2) using subcluster resolution 0.2. Subclusters 4, 5, 7, and 10 were removed based on upregulation of genes related to one or more nonendothelial cell types and/or general downregulation of endothelial marker genes.
Vessel segment identification
To examine AD-related changes along the cerebral arteriovenous axis (Vanlandewijck et al., 2018; Chen et al., 2020; L. Zhao et al., 2020), we used endothelial zonation-specific gene markers from a recent large-scale brain endothelial sequencing study (Yang et al., 2022) (Extended Data Table 5-1). This set included 62 arterial markers, 129 capillary markers, and 87 venous markers. For each donor and brain region, we tabulated the geometric mean of the log-normalized expression across the arterial, capillary, and venous gene markers, respectively. Each endothelial nucleus was classified according to the segment with the largest geometric mean expression of the three types.
Differential expression (DE) analysis
All DE analysis was performed using either the FindMarkers() or the FindAllMarkers() function in Seurat and the model-based analysis of single-cell transcriptomics (MAST) algorithm (Finak et al., 2015), as has been applied in previous AD single-cell sequencing studies (Morabito et al., 2020; Yang et al., 2022), unless otherwise specified. We included donor ID as a latent variable to be regressed out in all MAST modeling unless otherwise specified.
Endothelial cells versus nonendothelial cells
We examined gene expression across brain regions in low-pathology endothelial cells using FindMarkers(), opting to retain a gene as an endothelial-specific marker if the average log2 fold change (log2FC) was ≥1, if the gene was detected in at least 10% of endothelial nuclei, the p value was <0.05, and the Benjamini–Hochberg (BH)-adjusted p value was <0.1. The resulting set of endothelial marker genes was resupplied to FindMarkers() with no thresholding to calculate the donor-regressed average log2FC across all five cell types for comparison.
Region-specific low-pathology endothelial transcriptome analysis
To examine regional endothelial heterogeneity in normal aged brain, we subsetted the integrated endothelial dataset to only the low-pathology subjects (Braak 0/I/II; Extended Data Table 1-1) and performed DE analysis using FindMarkers(), with the brain region as the identity variable to compare. Genes were retained as significant for each region if the absolute average log2FC was >0.1, if the gene was detected in at least 10% of endothelial nuclei in the given region, the nominal p value was <0.05, and the BH-adjusted p value was <0.1. These thresholds were chosen based on previously published criteria for endothelial snRNAseq in AD (Lau et al., 2020).
Whole-brain high- versus low-pathology endothelial transcriptome analysis
We compared the endothelial transcriptome in the high- versus low-pathology donor cortex, first by combining nuclei from all five brain regions into one analysis. To focus on donors at either end of the neuropathological spectrum, we omitted intermediate-pathology donors (B2, n = 5) and kept low (B0/1; n = 11) and high-pathology donors (B3; n = 16). DE analysis was performed using FindMarkers() with the pathology level as the identity variable to compare (i.e., high vs low). Genes were retained as significant if the absolute average log2FC was >0.1, if the gene was detected in at least 10% of endothelial nuclei in the given region, the nominal p value was <0.05, and the BH-adjusted p value was <0.1.
Region-specific high- versus low-pathology endothelial transcriptome analysis
Upon identifying genes with DE based on pathology level across the whole brain, we further examined genes with region-specific DE. We reran FindMarkers() for each brain region separately, keeping donor ID as the latent variable to regress out. A gene was retained as significant for a given region if the absolute average log2FC was >0.1, if the gene was detected in at least 10% of endothelial nuclei in the region, the nominal p value was <0.05, and the BH-adjusted p value was <0.1.
Identification of gene sets that are temporally related to AD disease stages
To identify genes that were related to the pathologic progression of AD, we performed DE analysis between adjacent pathology groups (i.e., Group 1 vs Group 2; Group 2 vs Group 3; Group 3 vs Group 4; Extended Data Table 1-1) using FindMarkers() with region as the latent variable and logistic regression (instead of MAST). Genes with logFC > 0.25 and BH-adjusted p value < 0.05 and with expression >10% in all groups were included. The standardized gene expression was determined by scaling across all samples for each gene and computing the average score for each pathology group. The standardized gene expression scores were grouped using spectral clustering with SNFTool (version 2.3.1; k = 6) (Wang et al., 2014).
Plaque burden, pTau/Tau, and apolipoprotein E (APOE) allele endothelial transcriptome analysis
To directly investigate the association between the endothelial transcriptome and AD neuropathology burden, gene expression was modeled as a function of (1) amyloid β plaque burden, (2) CAA score, (3) pTau/Tau ratio, or (4) APOE E4 allele status, respectively. For all four analyses, data were filtered to only include genes with reads detected in at least 20% of endothelial cells and standardized variance >0.5. For the CAA analysis (Model 2), data were additionally filtered to remove samples with missing CAA scores (as indicated by “NA” in the CAA column in Extended Data Table 1-1). After filtering the data, for each of the four pathology variables, a zero-inflated MAST generalized linear mixed effects model was fit to regress endothelial cell gene expression on the variable of interest using the zlm() function. The donor ID was included as a random effect and the cellular detection rate (“cngeneson”) and percent mitochondrial DNA (“pc_mito”) variables were included as fixed effects, with the method set to “glmer” for a generalized linear mixed-effects model. The ebayes flag was set to FALSE to turn off variance regularization.
Functional pathway annotation
Pathway overrepresentation analysis
To examine which Gene Ontology (GO) biological processes (Ashburner et al., 2021; Gene Ontology Consortium, 2021) are enriched in endothelial genes of interest, we performed overrepresentation analysis using PANTHER with Fisher's exact test and the Calculate False Discovery Rate option (Mi et al., 2013). The resulting json files were downloaded and imported into R using the jsonlite package (Ooms, 2014). The 36,601 genes assayed in this study were provided as the reference gene set for PANTHER. We defined significance as nominal p value < 0.05 and adjusted for multiple comparisons in R using the p.adjust() function with BH correction.
Gene set enrichment analysis (GSEA)
To complement pathway overrepresentation analysis, we examined the enrichment of gene sets based on the log2FC ranking of all genes using GSEA (Mootha et al., 2003; Subramanian et al., 2005). This method is advantageous in that it does not rely on any specific fold change or significance threshold since it takes the full list of genes (Mootha et al., 2003; Subramanian et al., 2005). A log2FC-ranked list of genes was compiled in the .rnk format and GO biological processes were curated in the .gmt format by msigdb (Liberzon et al., 2011). We performed preranked GSEA using the fgseaMultiLevel() function from the fgsea package in R (Korotkevich et al., 2021), with the following parameters: minSize = 10, maxSize = 500, eps = 0.0, nPermSimple = 1000. The observed enrichment score (ES) is compared with an empirical null distribution of ESs generated by randomly shuffling phenotype labels and gene ordering 1000 times; the nominal p value is then calculated as 1 minus the proportion of permutations for which the observed ES was greater or lesser than the null distribution, for positive or negative ES, respectively. A pathway was considered significantly enriched (positive or negative) if the p value was <0.05 and the BH-adjusted p value was <0.25, as per GSEA recommendations (Subramanian et al., 2005).
Results
Identification of endothelial cells from multiple brain regions using single-nucleus RNA sequencing
To examine the region-specific transcriptome of brain endothelial cells in normal control brain and in AD, we performed snRNAseq on postmortem brain tissue in five regions (EC, ITG, PFC, V2, and V1) available from 16 high-pathology (B3) AD donors, 5 intermediate-pathology (B2) donors, and 11 low-pathology (B0/1) donors (Figs. 1–3; Extended Data Table 1-1). An overview of this experimental and analytical workflow is depicted in Figure 1A, B. We separated out most neurons and oligodendrocytes using fluorescence-based sorting, then clustered nuclei and identified endothelial cell nuclei based on von Willebrand factor (VWF) expression, and finally, removed subclusters likely to reflect contamination (see Materials and Methods; Figs. 3B, 4A,B; Extended Data Table 1-2). After these steps, we retained 19,271 genes measured in 51,586 endothelial nuclei (counts per brain region: EC, 3484; ITG, 9259; PFC, 10,377; V2, 13,110; and V1, 15,356). In addition to snRNAseq analysis, pathology was confirmed in each tissue block by tau ELISA and amyloid β histology (Fig. 2A,B). There was no significant difference in donor postmortem interval between B0/1 versus B3 donors (Fig. 2C), and nuclei from the five regions exhibited comparable read counts per nucleus, RNA integrity number, and total number of genes detected (Fig. 3C). Last, we confirmed that the number of endothelial cell nuclei did not differ by pathology group (Fig. 3D). Of note, nonendothelial cell types collected are concurrently being analyzed in separate studies (Serrano-Pozo et al., 2022; Alzheimer DataLENS, 2023).
Study design overview. A, NeuN-/Olig2-negative nuclei were isolated from five brain regions (EC, ITG, PFC, V2, V1) and sequenced (10× Genomics) to examine endothelial cell-specific transcriptome throughout the AD brain. B, Endothelial cells were identified in each brain region individually and isolated for integration, enabling the examination of region-specific and vessel segment-specific cerebrovascular transcriptomes. Undefined and donor-specific clusters removed from plots. C, DE analysis using Findmarkers() between low-pathology donor endothelial nuclei and low-pathology donor nonendothelial nuclei identified a total of 527 unique endothelial cell markers, 347 of which were shared across all five regions. D, Heatmap visualization of the 347 core endothelial marker genes across cell types depicts the endothelial-specific upregulation. The top eight genes across the regions are highlighted on the right (FLT1, CLDN5, ABCB1, ATP10A, ELOVL7, ANO2, MECOM, and VWF). Data corresponding to this figure are included in Extended Data Tables 1-1, 1-2, 1-3, and 1-4.
Neuropathological features are shown for each donor and brain region. A, The region-wise pTau:Tau ratio is plotted for each donor, labeled by ABC-B. B, The region-wise amyloid plaque area is plotted for each donor, labeled by ABC-C. Donors are ordered from greatest to least severity (left to right) in each panel. C, Postmortem interval (PMI) for donors by ABC-B pathology group. Values were compared using the Kruskal–Wallis nonparametric ANOVA (p = 0.57).
Quality control for endothelial cells sequenced. A, Summary of the number of donors with tissue for each brain region, stratified by ABC-B group. B, VWF expression overlaid onto the UMAP visualization of all nuclei highlights the endothelial nuclei cluster for each brain region. Undefined and donor-specific clusters removed. C, Standard quality control metrics plotted for each brain region confirm comparable quality across regions. D, Comparison of the total number of endothelial nuclei per donor (top row) and the proportion of endothelial cells to total nuclei isolated per donor (bottom row) across the Braak stages per brain regions. Values were compared using the Kruskal–Wallis nonparametric ANOVA; not significant (ns) indicates p > 0.05.
Quality control for Seurat integration. A, Elbow plots for the five brain regions show the SD explained by each of the first 50 principal components (PCs). B, Cell type marker visualization for the subclusters identified in endothelial cells integrated from all five regions identifies subclusters with heightened nonendothelial contamination (4, 5, 7, 10). Ast, Astrocyte; mic, microglia; oligo, oligodendrocyte; end, endothelial cell; opc, oligodendrocyte precursor cell; smc/per, smooth muscle cell/pericyte. Data corresponding to this figure can also be found in Extended Data Table 1-2.
To confirm that the endothelial cell populations were enriched for canonical cell type markers, we compared endothelial nuclei versus nonendothelial nuclei per region (Fig. 1C; Extended Data Table 1-3). Across the five brain regions, a total of 527 unique genes were identified as endothelial cell markers, of which n = 347 (65.8%) were shared across all five brain regions. The set of 347 pan-region endothelial marker genes was visualized in a heatmap, including other cell types to visually confirm enrichment of endothelial cells in each brain region (Fig. 1D; Extended Data Table 1-4). These include established endothelial cell markers, such as FLT1 (VEGF receptor-1), CLDN5 (tight junction protein claudin 5), ABCB1 (ATP-binding cassette subfamily B member 1; P-glycoprotein), and VWF. We then characterized the endothelial nuclei by vessel segment zonation, identifying 5046 arterial, 41,517 capillary, and 5023 venous nuclei, across the five brain regions (Extended Data Table 5-1). The proportion of capillary nuclei identified here (80.5%) is consistent with findings that capillary endothelial cells comprise ∼85% of the blood–brain barrier (Montagne et al., 2017).
Table 1-3
Genes significantly upregulated in endothelial cell nuclei relative to other cell types. (Data shown in Figure 1C). Download Table 1-3, XLSX file.
Table 1-4
Comparison of endothelial marker genes expression in other cell types (Data shown in Figure 1D). Download Table 1-4, XLSX file.
Regional transcriptomic variability in endothelial cells from low-pathology control brains
To investigate regional heterogeneity in the endothelial cell transcriptome, we examined nuclei isolated from aged donors with no to minimal AD-related neuropathology (B0/1 donors in Extended Data Table 1-1). UMAP visualization of the integrated low-pathology endothelial nuclei at 0.3 resolution did not show any clear differences between regions (Fig. 5A). However, in terms of vessel segment composition, there were clear differences in the relative proportions of arterial, capillary, and venous endothelial cells across the five regions (Fig. 5B). For example, the EC had more arterial nuclei than the other regions, while V1 had more capillary nuclei than the other regions. Whether this reflects genuine differences in vascular bed composition across the regions or experimental variability remains to be clarified.
Endothelial cells from low-pathology donor brain regions exhibit differences in transcriptomic pathways and vessel segment composition. A, UMAP visualization of the integrated low-pathology donor nuclei show that the brain regions are evenly distributed in UMAP space. B, UMAP visualization showing endothelial cell arterial-venous zonation of all nuclei. The proportions of nuclei across brain regions differed from expected. Right, χ2 test results. *p < 0.05. **p < 0.005. ***p < 0.001. C, Differential gene expression analysis showing upregulated genes in each of the five brain regions. D, Heatmap showing upregulated genes for each brain region are increased relative to other areas and GO biological process (BP) overrepresentation analysis identified the pathways that are enriched in each region's marker gene sets. Data corresponding to this figure are included in Extended Data Tables 5-1, 5-2, and 5-3.
Table 5-1
Markers used to identify vessel segments. Gene list adapted from Yang et al., 2022(Yang et al., 2022). (Data shown in Figure 5B). Download Table 5-1, XLSX file.
Table 5-2
Differentially upregulated endothelial cell markers between brain regions in low-pathology subjects. (Data shown in Figure 5D). Download Table 5-2, XLSX file.
Table 5-3
Region-specific GO biological processes identified from pathway overrepresentation analysis (PANTHER) on expression levels from Extended Data Table 5-2. (Data shown in Figure 5D). Download Table 5-3, XLSX file.
We then identified whether individual genes varied across brain regions and identified between 192 and 379 genes (or 1%-2% of all genes detected in endothelial cells) that were uniquely upregulated in one of the regions relative to the other four, making them region-specific marker genes (Fig. 5C). Critically, donor ID was included as a latent variable in this, and all subsequent, analyses. Relevant region-specific genes included APOE, C1QB, and CD93 in EC; CTNAA2, LSAMP, and CADM1 in ITG; CALM2, ANXA3, and TM6SF1 in PFC; RNASE1, ICAM2, and SLC2A1 in V2; and HSPH1, HMCN1, and FLT1 in V1. These regional marker genes were then grouped into biological processes (Fig. 5D; Extended Data Tables 5-2 and 5-3). EC marker genes were enriched for cytokine production, regulation of lipid localization, and oxidative stress response. ITG marker genes were enriched for Wnt signaling, endothelial cell migration, and axonogenesis. PFC marker genes were enriched for microtubule organization, epigenetic gene expression, and histone modification. V2 marker genes were enriched for cytoplasmic translation, downregulation of cell migration, and vasculogenesis. V1 marker genes were enriched for virus response, “de novo” protein folding, and angiogenesis regulation. These data together suggest that differences in baseline gene expression exist in each of the five regions examined in the aging brain even before AD pathology develops.
AD is associated with unique sets of dysregulated gene networks in endothelial cells
To examine AD-related changes in each brain region, we first analyzed AD-associated transcriptomic differences in endothelial cells from all five brain regions combined by comparing high-pathology (B3) versus low-pathology (B0/I) donors. We identified 936 genes that were upregulated in high-pathology donor endothelial nuclei and 962 that were downregulated (Fig. 6A; Extended Data Table 6-1). Among the top 10 upregulated genes, eight of them are heat shock proteins (HSPs) with log2 fold change (log2FC) values of 2.18-1.04, representing more than double the expression of these genes in high-pathology AD individuals compared with low-pathology controls. Less extensive fold changes were observed among the top 10 downregulated genes, which had log2FC values of −0.8 to −0.51. To validate the results of this analysis, we performed immunohistochemistry in sections from the temporal and visual cortices of donors used in this study. To validate expression changes, we selected AD-upregulated genes MT2A, CREM, IGFBP3, HSPH1, and HSP90AA1, which represent a range of small to large effect sizes, and confirmed expression at the protein level in capillary endothelial cells (Fig. 6B). No specific cortical layer differences were observed in the expression of these proteins by histology.
Whole-brain DE analysis reveals dysregulated pathways that are specific to endothelial cells in the AD brain. A, In high- versus low-pathology endothelial nuclei, there are 936 upregulated genes and 962 downregulated genes. B, Immunofluorescence images of five representative genes that are upregulated in AD brain endothelial cells (white). Glut1 colabeling is used to identify endothelial cells (magenta). All images are from V1/V2, except for MT2A, which is from EC. C, GSEA reveals that high-pathology endothelial cells are enriched for several pathways, including those related to protein folding, IL-8 production, leukocyte adhesion, and blood–brain barrier maintenance. Only a subset of GO biological processes are shown in this figure. NES, Normalized ES. D, DE analysis in each cell type shows that 522 of 936 of the upregulated genes and 501 of 962 of the downregulated genes are specific to endothelial cells. E, Heatmap visualization of these 1023 genes depicts their specificity to endothelial cells in AD. Data corresponding to this figure are included in Extended Data Tables 6-1, 6-2, and 6-3.
Table 6-1
Differential expression analysis for high- versus low-pathology endothelial cells, all brain regions together. (Data shown in Figure 6A). Download Table 6-1, XLSX file.
Table 6-2
Pathology-associated GO biological processes identified from GSEA on expression levels in Extended Data Table 6-1. (Data shown in Figure 6C). Download Table 6-2, XLSX file.
Table 6-3
Expression levels across cell types for genes uniquely up- or down-regulated in high-pathology endothelial cells. (Data shown in Figure 6D). Download Table 6-3, XLSX file.
Next, we examined whether endothelial cells from high-pathology donors were enriched for different biological processes compared with low-pathology donors using GSEA. For each GO biological process, we calculated its normalized ES, which indicates how much the genes in each gene set are overrepresented at the top or bottom of the entire ranked list of genes, normalized to gene set size. Among all 19,271 genes, a total of 660 biological processes were significantly positively enriched in the high-pathology endothelial cells, including upregulated pathways related to protein folding, cytokine production, blood–brain barrier maintenance, and leukocyte adhesion (Fig. 6C; Extended Data Table 6-2). Of the 660 biological processes, 68 pathways were downregulated, which included those related to lipid and glycoprotein metabolism.
Previous studies that examined the AD cerebrovascular transcriptome have not clarified the endothelial cell specificity of resulting gene expression changes. To address this, we calculated the average high- versus low-pathology log2FC for the n = 1898 genes from Figure 6A for each of the other four cell types (astrocytes, microglia, neurons, and oligodendrocytes). This revealed that 522 of the upregulated genes and 501 of the downregulated genes are unique to endothelial cells (Fig. 6D,E; Extended Data Table 6-3). Of the genes that were not unique to endothelial cells, most of the upregulated genes were shared with microglia (n = 122) and most of the downregulated genes with oligodendrocytes (n = 170). Of note, 32 of the upregulated genes unique to endothelial cells were identified in Uniprot to be secreted proteins and may thus represent potential plasma or CSF biomarkers for brain endothelial cell changes. These include ADM, BMP6, EGFR, IL32, IGFBP3, TGFB1, and WNT2B.
Brain regions exhibit different endothelial transcriptomes amid AD pathology
Since the five brain regions assayed are at different stages along the Braak continuum (Braak and Braak, 1991; Braak et al., 2006), we next asked whether biological processes that were identified from analysis of total endothelial nuclei were more specific to individual regions. We reexamined these 16 pathways for each region (shown in Fig. 6C), comparing their normalized ESs for high-pathology versus low-pathology subjects (Fig. 7A; Extended Data Table 7-1). Of the 16 biological processes, seven were significantly enriched in all five regions: protein folding, oxidative phosphorylation, cellular response to heat, regulation of intrinsic apoptotic signaling pathway, vasculogenesis, aging, and positive regulation of cytokine production. Of note, maintenance of the blood–brain barrier was significantly enriched in PFC, ITG, and to a lesser extent in V1, but not in V2 or EC.
Regional differences in endothelial cell responses to AD pathology. A, The same pathways highlighted in Figure 3C are shown in region-wise GSEA comparing high- versus low-pathology donors, showing that some pathways are enriched in all five regions (e.g., protein folding, vasculogenesis) while others are more region-specific (e.g., blood–brain barrier maintenance, regulation of vascular permeability). NES, Normalized ES. B, The number of differentially expressed genes increases from EC to V1 along the Braak stage continuum. C, Euler diagram visualization of the genes of interest across regions shows that most genes are uniquely dysregulated in one or two regions, while a few are commonly dysregulated in all five regions (38 upregulated and 21 downregulated). Data corresponding to this figure are included in Extended Data Tables 7-1 and 7-2.
Table 7-1
Pathology-associated endothelial-specific GO Biological processes identified from GSEA on expression levels in Extended Data 6-3, analyzed at the regional level. (Data shown in Figure 7A). Download Table 7-1, XLSX file.
Table 7-2
Differential expression analysis for high- versus low-pathology endothelial cells, separated by brain region. (Data shown in Figure 7B). Download Table 7-2, XLSX file.
We further examined the AD-related endothelial transcriptome within each region separately. Region-wise gene expression changes appeared to increase along the Braak stage neuroanatomical continuum, with the lowest number in the EC (228 up, 219 down) and the highest number in the V1 (981 up, 1061 down) (Fig. 7B; Extended Data Table 7-2). Nearly all genes were protein-coding genes; of these, 13%-33% were exclusive to one region, while others were shared between many or all regions (Fig. 7C). These included the upregulation of 38 genes enriched for detoxification and metal ion response (metallothionein genes MT1E, MT1F, MT1M, MT1X, and MT2A) and notably included PICALM, which has been previously associated with AD risk in genome-wide association studies and is involved in vascular clearance of Aβ pathology (Harold et al., 2009; Z. Zhao et al., 2015; Kisler et al., 2023). The 21 common downregulated genes included ones related to inhibition of TGF-β signaling (SMAD6, SMAD7), endocytosis (RFTN1, MYOF), and enzymes that modify basement membrane and glycocalyx (HGSNAT, CHST15, SULF2, GALNT17). These AD-related genes showed limited overlap with the region-specific genes (i.e., a gene identified in low-pathology nuclei as a marker of V1 endothelial cells was not also altered throughout the brain with AD pathology). Together, these findings indicate that endothelial cells within specific brain regions express diverse biological pathways in AD.
Endothelial cell gene expression programs exhibit diverse temporal relationships with disease progression
One of the limitations of human tissue studies is the difficulty in resolving the temporal relationships between disease stage and transcriptomic changes. This study was designed a priori to help address this by including donor samples that fall into four equally sized (n = 8) groups based on Braak stage (0-VI) and CERAD score (neuritic plaques; Extended Data Table 1-1). Pairwise comparisons between adjacent groups yielded genes that were increased or decreased with pathologic burden, and clustering of these genes revealed six common trends in gene expression with disease progression (Fig. 8; Extended Data Table 8-1). Some of these genes showed trajectories that decreased linearly with pathology, including ion transport and intracellular signal transduction (Gene Set 1), while others were enriched for pathways that increased linearly, including genes related to metal ion homeostasis, apoptosis, and cell cycle regulation (Gene Set 2). By comparison, genes related to lipid localization and cytoskeleton organization appeared relatively consistent across disease stages (Gene Set 3). Gene Set 4, which includes protein folding, cellular response to heat, and regulation of proteolysis genes, also showed increases in intermediate stage pathology donors but then returned to baseline in late-stage disease. Other transient changes in intermediate states were observed in Gene Set 5, which included antigen processing, calcium ion response, and vasculature development genes. Lastly, Gene Set 6 included vascular transport and cell signaling genes that were only reduced in an intermediate stage. In all, these data indicate that endothelial cell gene expression changes follow complex, nonlinear trajectories in response to disease progression.
Temporal differences in endothelial cell responses to AD pathology. A, Six gene sets identified by pair-comparisons between tissue from four distinct donor categories: 1 = Braak 0-II, no neuritic plaques; 2 = Braak II/III, sparse-moderate neuritic plaques; 3 = Braak V, moderate neuritic plaques; 4 = Braak VI, severe neuritic plaques. Graphs represent standardized gene expression on the y axis across each of the four disease stage categories. Gray lines indicate individual genes. Colored line indicates the mean trend for each gene set. B, Selected top GO biological process terms plotted on the left with example genes from each set listed on the right. Data corresponding to this figure are included in Extended Data Table 8-1.
Table 8-1
Gene sets associated with AD pathology temporal dynamics. (Data shown in Figure 8). Download Table 8-1, XLSX file.
Distinct endothelial gene changes are related to amyloid plaque burden, CAA, tau and APOE genotype
To examine how these endothelial gene expression changes relate to AD neuropathology, we performed linear regression analysis of individual gene expression on histologic amyloid plaque burden, phosphorylated tau content, and APOE genotype. With regards to amyloid β plaque load (measured by 3D6-positive area fraction via immunohistochemistry on adjacent tissue sections), 864 genes were found to be positively or negatively associated with plaque burden, including 89 that were significant after correcting for multiple comparisons (Fig. 9A). In this subset, top genes that were positively associated with amyloid plaque burden include TMSB10, EEF1A1, TPT1, KLF2, NFKBIA, and ICAM2, and genes that were negatively associated with amyloid plaque burden include CTNNA2, GLIS2, FBXL7, PDZRN3, FHIT, and ABCA1 (Extended Data Table 9-1). Biological processes that were overrepresented in this gene set were related to metabolism, inhibition of apoptosis, zinc homeostasis, translation, and endosomal transport (Extended Data Table 9-2).
Genes associated with amyloid β plaques and CAA. A, Examples of low and high amyloid plaque burden across brain regions and volcano plots showing genes that are positively or negatively associated plaque burden. B, Examples of CAA grading severity (0-3) across 4 donors and volcano plots showing genes that are positively or negatively associated CAA severity. Data corresponding to this figure are included in Extended Data Tables 9-1, 9-2, 9-3, and 9-4.
Table 9-1
Genes associated with amyloid plaque burden. (Data shown in Figure 9A). Download Table 9-1, XLSX file.
Table 9-2
Table containing amyloid plaque-associated GO Biological processes obtained from enrichment analysis. (Data shown in Figure 9A). Download Table 9-2, XLSX file.
Table 9-3
Genes associated with CAA pathological burden. (Data shown in Figure 9B). Download Table 9-3, XLSX file.
Table 9-4
Table containing CAA-associated GO Biological processes obtained from enrichment analysis. (Data shown in Figure 9B). Download Table 9-4, XLSX file.
We also examined whether endothelial gene expression was affected by cerebral amyloid angiopathy (CAA; scored 0-4 based on severity using paraffin-embedded sections from the opposite hemisphere). This yielded 684 genes associated with CAA, including positively associated genes BACE2 and LRP1B and negatively associated genes SLC7A5, SLC16A1, SLC2A1, TFRC, FTH1, and CLU. Overall, 37 remained significant after corrections for multiple comparisons, such as CADPS22, IQCJ-SCHIP1, PTGDS, and MTUS1 (Fig. 9B; Extended Data Table 9-3). Of the CAA-related genes with nominal p value s <0.05, only 11 positively associated genes and 34 negatively associated genes overlapped with those related to amyloid β plaques. CAA-related genes included overrepresented biological processes for regulation of amyloid β formation, transport across the blood–brain barrier, iron homeostasis, and response to VEGF stimulus (Extended Data Table 9-4). Collectively, this suggests that CAA-related changes are distinct from changes related to amyloid β plaques.
Finally, we examined the relationship between gene expression and phosphorylated tau (measured biochemically as a ratio between p231 tau and total tau in the same bulk tissue used for the snRNAseq) or APOE genotype (E3 vs E4). This resulted in 154 and 177 genes associated with these variables, respectively, but none was statistically significant after correction for multiple comparisons (Extended Data Tables 9-5 and 9-6). Further, no biological processes were overrepresented in these gene sets. Together, few of the genes associated with amyloid plaque burden, CAA, phosphorylated tau, and APOE overlapped, and no genes correlated with all four variables. This indicates that AD-related neuropathological factors have distinct effects, and that amyloid β exerts the greatest impact on endothelial cell gene expression independent of brain region.
Table 9-5
Genes associated with ptau/tau levels. Download Table 9-5, XLSX file.
Table 9-6
Genes associated with presence of the APOE4 allele. Download Table 9-6, XLSX file.
High-tau endothelial cells recapitulate published enrichment of senescence, endothelial aging, and CAA pathways
To anchor the findings presented herein with previous studies examining the cerebral vasculature in health and disease, we curated a set of genes identified in stroke, CAA, aging, and small vessel disease and compared them with this dataset (Extended Data Table 10-1). At the whole-brain level, high- versus low-pathology endothelial cells were significantly enriched for genes related to CAA (Moursel et al., 2018), endothelial aging (L. Zhao et al., 2020), senescence response (Dehkordi et al., 2021), and the endothelial senescence-associated secretory phenotype (Kiss et al., 2020) (Fig. 10A,B). These cells also showed a trend (nominal p value < 0.05, adjusted p value < 0.1) for enrichment of genes involved in senescence initiation (Dehkordi et al., 2021). At the regional level, two pathways (CAA and endothelial aging, Fig. 10C) remained significantly enriched across all five regions, whereas senescence-related pathways were more region-specific, particularly in the ITG and PFC. Genes associated with stroke, white matter hyperintensities, small vessel disease, and other disorders (progressive supranuclear palsy and frontotemporal dementia) were not enriched. Within these donors, these data confirm key endothelial changes observed in previous datasets and highlight the importance of considering regional differences when examining specific pathways.
Table 10-1
List of gene sets related to aging and neurodegenerative disease that were compared to this dataset. (Data shown in Figure 10). Download Table 10-1, XLSX file.
Enrichment of aging and disease-associated pathways in endothelial cells. Comparison with previous publications highlights that high-pathology donor endothelial cell nuclei are enriched for gene sets related to CAA, endothelial aging, senescence response, and the senescence-associated secretory phenotype (SASP). A, Normalized ESs (NES) are shown for each of the 12 gene sets derived from recent AD and/or vascular disease-related publications. All brain regions were included in this analysis. B, GSEA plot visualization depicts the genes' rank in the ordered dataset versus running ES in each of the significantly enriched modules, respectively. C, Region-wise GSEA reveals that CAA and endothelial aging gene sets are enriched in endothelial cells from all five regions, while senescence response, SASP, and senescence initiator pathways are region-specific. *p < 0.05. **p < 0.005. ***p < 0.001. #p < 0.05. BH-adjusted p > 0.25. Data corresponding to this figure are included in Extended Data Table 10-1.
Discussion
Using a large dataset of endothelial cell nuclei, the transcriptomic analysis that we report reveals new insights into regional gene expression differences in normal aged brain and AD. These data indicate that endothelial cells from different cortical areas are clearly distinguishable and are characterized by the expression of different marker genes and biological pathways. We also observed that AD pathology affects gene expression in brain endothelial cells, particularly those that make up capillaries. By comparing across donor samples from multiple disease stages, we observed the temporal relationship between pathology and endothelial transcriptome changes. Further, using these data, we compared the extent that AD-related factors, including plaques, CAA, tau, and APOE genotype, impact this cell type and examined how previously reported transcriptomic signatures vary between brain regions.
A fundamental biological question we set out to answer was whether the vascular transcriptome differs across brain regions. Previous endothelial cell transcriptomic datasets have been derived from mouse models (He et al., 2018; Vanlandewijck et al., 2018), where examination of distinct cortical areas is challenging, and existing human datasets include only a single cortical area per donor (Lau et al., 2020; Yang et al., 2022). By incorporating multiple cortical areas per donor, we showed that significant heterogeneity exists in endothelial cells, with 1%-2% of total genes detected in endothelial cells differentially expressed in each region. Previous assays of transcriptomic differences that used mouse endothelial cells from brain, heart, and lung showed that DE of 3%-9% of the total genes was responsible for tissue specificity (Jambusaria et al., 2020). Similarly, mouse and human brain endothelial cells show a 4.4% divergence in total gene expression (Yang et al., 2022). This suggests that the number of genes with expression differences across brain regions in the normal aged brain could have a relatively large impact on endothelial cell phenotypes. This warrants further study to better understand functional implications of this regional heterogeneity in the healthy endothelial transcriptome.
A better understanding of these transcriptomic changes could help explain why relatively similar cortical brain areas vary in susceptibility to CAA or neurodegeneration. Interestingly, visual cortex areas, which are affected late in AD progression and experience less neurodegeneration, expressed more genes related to vasculogenesis and angiogenesis. This transcriptomic difference might explain the relatively greater vascularization and distinct organization of visual cortex versus other brain areas (Duvernoy et al., 1981; Schmid et al., 2019). By comparison, highly vulnerable areas, such as the EC, expressed more oxidative stress-related genes in normal aged brain, suggesting endothelial dysfunction in this region even in the absence of severe AD pathology. A limitation of these studies is that analysis was restricted to cortical areas and vulnerable subcortical structures were not examined. In particular, the hippocampal vasculature is known to be disrupted in aging and Alzheimer's disease. Yang et al. (2022) recently reported differences in the endothelial cell transcriptome of healthy adult hippocampus relative to the superior frontal cortex (see also Montagne et al., 2015). In sum, these data support the idea that important differences in endothelial beds exist among brain regions. Whether these are because of development (i.e., origination of capillaries from the posterior cerebral artery vs middle cerebral artery), local structure, or age-related changes in various brain areas remains to be examined.
In comparing donor tissue from high- and low-pathology brains, the most highly upregulated genes were members of the HSP family. HSPs are important components of translational machinery in endothelial cells, and such upregulation could be a response to stress stimuli, including inflammatory cytokines, oxidized LDL, elevated blood glucose, and fluid sheer stress (Hochleitner et al., 2000). Exposure of cultured endothelial cells to metal ions is also known to result in upregulation of HSPs (Wagner et al., 1999). Increased expression of these proteins could also indicate a switch from quiescence to an activated, proliferative state (Pawlowska et al., 2005). Importantly, these genes appeared to have distinct trajectories and increased in intermediate disease stages only to return to baseline in late-stage disease. Because HSPs are a critical component of proteostasis mechanisms, we speculate that their elevation throughout the AD brain raises the possibility of brain-wide proteostatic stress that impacts cells without visible protein aggregates and that failure of this pathway in late-stage disease may contribute to severe phenotypes. Of note, inhibitors of HSP90 and its co-chaperones have been shown to attenuate both tau and amyloid pathology in animal models (Blair et al., 2014). Considering that these genes were increased in all brain areas examined, and not significantly associated with amyloid plaque burden, CAA, phosphorylated tau levels, or APOE genotype, it is possible that this is a response to neurodegeneration more generally. Future studies should expand the scope to non-AD neurodegenerative diseases, such as frontotemporal dementia or Lewy body disease, to clarify the specificity of AD in dysregulated expression of these pathways.
We also explored regional heterogeneity in response to AD pathology and observed that endothelial responses to pathology are not uniform across brain regions, with a greater number of gene expression differences in the visual cortex compared with other areas. One explanation could be more age- and/or disease-related pathology accumulation in other brain areas, such as the EC, resulting in more similarities and fewer gene expression differences between normal aged brain and AD. It is notable that the Moursel et al. (2018) CAA gene set was prominently enriched in AD endothelial cells, and that CAA is more frequently encountered in occipital regions, such as the visual cortex, compared with anterior areas (see also Vinters and Gilbert, 1983). This is also in line with our finding that both CAA and amyloid plaque burden are more strongly associated with endothelial gene expression changes than tau pathology and/or the presence of the APOE4 allele in this dataset. However, future studies that include greater numbers of APOE4 carriers, including AD APOE4 homozygotes and control subjects without AD, will be required to specifically address the contribution of this risk allele to endothelial cell transcriptomic changes.
Last, it was surprising that genes associated with amyloid β plaques minimally overlapped with genes associated with CAA in this dataset. Both CAA and amyloid β plaques frequently, but not always, co-occur within regions (Vinters et al., 1996; Thal et al., 2003). This indicates that, while both result in local increases in pathologic amyloid β, they appear to exert different effects on endothelial transcription programs. The form of pathologic amyloid β, vessel segment, and specific location within cortex may contribute to this observation. We note that CAA predominantly affects leptomeningeal vessels (excluded from this analysis via gross dissection) particularly near the occipital cortex, while plaques are often found more diffusely in deeper cortical layers (Rogers and Morrison, 1985; Charidimou et al., 2017). As such, future studies that examine both pial vessels and parenchyma using spatial transcriptomics could determine how local microenvironments contribute to the observed gene expression programs in endothelial cells. Additional analysis of other cell types, including smooth muscle cells and pericytes, which were not examined here, and their interactions with endothelial cells could further shed light on these amyloid-related changes.
The current study extends and reinforces our earlier work examining gene expression changes in isolated vasculature, and work from other groups that have performed single-cell analysis using similar preparations (Lau et al., 2020; Yang et al., 2022). Previously, we had observed that, with increasing severity of AD pathology, isolated blood vessels expressed more cellular senescence-related changes as measured by qPCR (Bryant et al., 2020). The present work shows that senescence-related gene signatures are increased across several brain regions, and confirms these changes in endothelial cells in the absence of other vascular cell types. Senescent cells are known to have paracrine effects on surrounding tissue, often inducing sprouting of new dysfunctional blood vessels. Consistent with this, other changes in the proliferative capacity of endothelial cells were also noted, including the enrichment of genes related to vasculogenesis at the whole-brain level. Previous reports have observed altered angiogenesis, blood–brain barrier maintenance, and leukocyte adhesion gene changes in samples from PFC (Lau et al., 2020). Our analyses add to this and show that vasculogenesis- and angiogenesis-related genes have distinct regional and temporal expression profiles.
Together, these data can be used to guide future studies exploring disease processes in endothelial cells and indicate that careful tissue selection is warranted, as cortical regions are not all equivalent and show key differences at baseline in the aged brain and in response to AD pathology. We found that amyloid β plaques and CAA have the greatest impact on gene expression programs in endothelial cells from AD donors. Last, while endothelial cells are not typically associated with protein aggregation, upregulated protein folding pathways suggest that proteostatic stress is a key pathway in this cell type.
Data availability
All raw sequencing data is deposited under the NCBI sequence read archive (SRA), with the accession number PRJNA916657. Data may also be explored via an interactive web browser (http://alzdatalens.partners.org/; Alzheimer DataLENS, 2023). Detailed instructions for accessing this tool are provided in Figure 11.
Interactive data visualization. A, Data pertaining to this publication can be viewed at http://alzdatalens.partners.org. From the landing page, select log in on the right-hand tab. B, Log in using user-specific credentials. C, After log in, select “AD Progression Study” and “Sample-Level Analysis.” D, On the left-hand side, the five ROIs and different cell types can be selected. A gene of interest is input into the search box and selected from the drop-down menu—Occludin (OCLN) is shown here as an example. Data may be visualized as correlation plots with amyloid plaque burden, ptau/tau levels, or as box plots (by Braak stage, CERAD score, pathology groups, or by ABC scores).
Footnotes
This work was supported by National Institutes of Health Grant R01AG071567 to R.E.B., Grant P30AG062421 to B.T.H. and S.D., Grant K08AG064039 to A.S.-P.; and MASSCATS. We thank Patrick Dooley and Teresa Connors for help with tissue selection; Derek Oakley and Matthew Frosch for neuropathological expertise; Simon Dujardin and Sarah Hopp for intellectual support; and the many donors and their families who have contributed to the Massachusetts Alzheimer's Disease Research Center.
M.E.W., A.W., G.L., T.K., R.V.T., K.B., and E.H.K. are employees of AbbVie. The design, study conduct, and financial support for this research were provided by AbbVie. AbbVie participated in the interpretation of data, review, and approval of the publication. B.T.H. has received research funding from AbbVie as part of a collaboration agreement with the General Hospital Corporation, d/b/a Massachusetts General Hospital. B.T.H. has a family member who works at Novartis and owns stock in Novartis; he serves on the SAB of Dewpoint and owns stock. He serves on a scientific advisory board or is a consultant for AbbVie, Avrobio, Axon, Biogen, BMS Cell Signaling, Genentech, Ionis, Novartis, Seer, Takeda, the U.S. Department of Justice, Vigil, and Voyager. His laboratory is supported by Sponsored research agreements with AbbVie, F Prime, and research grants from the National Institutes of Health, Cure Alzheimer's Fund, Tau Consortium, and the JPB Foundation. A.B., Z.L., A.S.-P., S.D., and R.E.B. work on the AbbVie-Hyman Collaboration and have no other funding to disclose. The remaining authors declare no competing financial interests.
- Correspondence should be addressed to Rachel E. Bennett at rebennett{at}mgh.harvard.edu
This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International license, which permits unrestricted use, distribution and reproduction in any medium provided that the original work is properly attributed.
