A widespread decrease of chromatin accessibility in age-related macular degeneration

Age-related macular degeneration (AMD) is a leading cause of blindness in the elderly. The extent to which epigenetic changes regulate AMD progression is unclear. Here we globally profiled chromatin accessibility in the retina and retinal pigmented epithelium (RPE) from AMD patients and controls. Global decreases in chromatin accessibility occurr in RPE in early AMD, and in the retina with advanced disease, suggesting that dysfunction in RPE cells drives disease progression. Footprints of photoreceptor and RPE-specific transcription factors are enriched in differentially accessible regions (DARs). Genes associated with DARs show altered expression in AMD. Cigarette smoke treatment of RPE cells recapitulates epigenomic changes seen in AMD, providing an epigenetic link between the known risk factors for AMD and AMD pathology. Finally, overexpression of HDAC11 is partially responsible for the reduction in chromatin accessibility, identifying potential new targets for treatment of AMD.

Age-related macular degeneration (AMD) is by far the most common cause of irreversible visual impairment in people over 60 1 . The estimated number of people with AMD in 2020 is 196 million, and will increase substantially with the aging of the global population 2 . The disease is characterized by the early appearance of drusen, pigmentary abnormalities of the retinal pigment epithelium (RPE), and progressive photoreceptor dysfunction that is restricted primarily in the macula, a 6 mm diameter region of the fundus 3 . Although treatments aimed at inhibiting blood vessel growth can effectively slow the progression of the "wet" AMD, no useful treatments exist for the atrophic ("dry") form of the disease, which account for 90% of all AMD cases 4 .
Currently, GWAS analysis has identified at least 34 AMD genetic risk loci involved in the regulation of complement pathway and inflammation 5,6 . However, these genetic variants can only explain a subset of AMD cases, suggesting a substantial role of environmental factors in the pathogenesis of AMD. Indeed, studies have linked variables such as cigarette smoking and obesity to AMD susceptibility, both of which are known to induce cellular stress and inflammation in a wide range of tissues 7,8 . However, no comprehensive studies have yet been reported on global epigenetic changes associated with AMD progression. This in part, reflects the lack of widely-accepted animal models for AMD 4 , as well as the difficulty in obtaining sufficient amounts of human pathological tissue for analysis. Here we observed global and progressive decreases in chromatin accessibility associated with AMD progression. Both cigarrette smoking treatment and overexpression of the epigenetic regulator HDAC11 in human iPSC-derived RPE recapitulated the changes in chromatin accessibility. These findings suggest that global decreases in chromatin accessibility may play a critical role in the onset and progression of AMD.

Landscape of chromatin accessibility in the retina and RPE
In this study, we obtained 16 eyes from 5 controls and 5 AMD donors (Extended Data Table 1). We collected retina and RPE tissues from the macular and peripheral regions of each donor eye, which altogether yielded a total of 19 normal, 9 early AMD, and 17 late AMD-derived samples (Table 1). Although the procurement time is slightly longer for normal samples, major characteristics including gender and age are comparable among normal and AMD samples. Disease status was confirmed with visual examination by an expert observer (J.T.H.).
To study the global epigenetic landscape of AMD, we used the assay for transposaseaccessible chromatin using sequencing (ATAC-Seq) to detect genomic chromatin accessibility, which depicts the active (i.e. open) chromatin and inactive (i.e. condensed) chromatin 9 . We obtained an average of 78.5% mappability and 35.8 million qualified fragments per sample (Extended Data Table 2). ATAC-Seq data from two replicate samples, obtained from adjacent regions in peripheral retina of the same eye, showed high correlation (R = 0.98, Extended Data Fig. 1a), indicating that ATAC-Seq can reliably and reproducibly measure chromatin accessibility in these samples. In total, 78,795 high-confidence open chromatin regions (or peaks) were identified across all retinal samples, and 49,217 peaks were identified across all RPE samples, representing a total of 93,863 distinct peaks. Chromatin accessibility in the retina is overall higher than that of the RPE, potentially reflecting the much greater diversity of cell types in the retina relative to the RPE (Fig. 1a). Comparison of samples from the macular vs. peripheral retina, as well as the macular vs. peripheral RPE, showed broadly similar profiles of chromatin accessibility (Fig. 1a).
The data revealed categories of peaks that are either specific to, or shared between, the retina and RPE. For instance, a peak associated with RLBP1 is shared by the retina and RPE, whereas a peak associated with SLC1A2 is specific to the retina, and another peak within the SLC45A2 gene is RPE-specific (Fig. 1b). Furthermore, a small number of region-specific peaks (e.g. peaks in KCNC2 and ZIC1) are selectively detected in the macular and peripheral retina, respectively, while others (e.g. peaks in PKD1L2 and ALDH1A3) are selectively accessible in macular and peripheral RPE (Fig. 1b) (Fig. 1c). We identified 5,855 increased and 2,689 decreased peaks in the macular retina, relative to the paired retinal samples from the peripheral region (Extended Data Fig. 1b). Meanwhile, 432 increased and 959 decreased peaks were detected in macular RPE relative to peripheral RPE (Extended Data Fig. 1c). We observed that a great majority (81.7%) of the peaks that are shared between the retina and RPE are also detected in other tissues ( Fig. 1c). In contrast, only 4,626 (12%) of retina-specific peaks and 7,644 (48.2%) of RPE-specific peaks are detected in other tissues, implying that these peaks potentially represent highly tissue-specific cis-regulatory elements.
We then calculated the overall similarity of the ATAC-Seq profiles among all the samples using multidimensional scaling. As expected, this analysis showed that the samples are clustered into two groups, one from the retina and the other from the RPE (Fig. 1d). Moreover, most AMD samples are clearly separated from normal samples, suggesting an extensive difference in chromatin accessibility between healthy and AMD tissues.

Chromatin accessibility is broadly decreased in the retina and RPE of AMD samples
To explore the impact of AMD, we analyzed the differences in chromatin accessibility between normal and AMD retinas. When comparing the accessibility profiles, we noticed substantial quantitative differences in peak signal between normal and AMD retina samples. For example, in three known regulatory regions of the rhodopsin gene RHO, chromatin accessibility is progressively decreased from normal to early-stage, and then to late-stage AMD (P < 0.05, Fig. 2a). By comparing the signal for each peak in healthy and AMD samples from both macular and peripheral retinas, we observed that 72,689 (92.3%) peaks have reduced chromatin accessibility in AMD (Fig. 2b). These quantitative differences in chromatin accessibility do not result from the process of normalizing ATAC-Seq data because different normalization approaches gave similar results (Extended Data Fig. 1d). Moreover, we separated retina samples into two groups from the macular and peripheral regions. Relative to the peripheral region, we observed more intense global decrease of chromatin accessibility in the macular region of AMD retina (94.5% for macula and 79.9% for periphery with the reduced chromatin accessibility, Extended Data Fig. 1e and 1f).
To extend this observation, we obtained a pair of eyes from a donor whose AMD status was asymmetrical, with the right eye showing early-stage AMD, and the left eye showing late-stage dry AMD (Fig. 2c). By comparing these eyes, we excluded the contribution of potential genetic and environmental differences that might complicate the analysis of epigenetic changes associated with AMD progression. Interestingly, a large number (76.7%) of peaks in the macular retina from the more severely affected eye had decreased intensities relative to the less severely affected eye (Fig. 2d). For the donors whose left and right eyes had been diagnosed at the same disease stage, the chromatin accessibility profiles were highly symmetrical in their macular retinas ( Fig. 2e and Extended Data Fig.   2a). This analysis confirmed that a widespread decrease in chromatin accessibility is associated with AMD progression.
Next, we analyzed changes of RPE chromatin accessibility in AMD. In all RPE samples, a great number (91.6%) of peaks showed the reduced intensity in AMD relative to normal samples (Fig. 3a). This reduction in intensity associated with AMD RPE was observed in both macular and peripheral regions (88.6% for macula and 94.1% for periphery showed reduced ATAC-Seq signal) (Extended Data Fig. 2b and 2c). In the RPE from the patient who showed different stages of AMD between eyes, the intensities of 42,860 (87.1%) peaks were reduced in the more severely affected left eye (Fig. 3b). In contrast, a symmetrical distribution was observed in donors where both eyes were in the same disease stage ( Fig. 3c and Extended Data Fig. 2d). Taken together, our data show a widespread decrease in chromatin accessibility that is observed in both the retina and RPE from AMD patients.

The retina and RPE show decreased chromatin accessibility at different stages of AMD progression
We further set out to identify changes in ATAC-Seq peak intensity that were associated with disease stage in both the retina and RPE. When we compared 5 retinal samples obtained from early-stage AMD to 11 retinal samples obtained from healthy controls, we observed only 3 statistically significant decreases in peak intensity ( Fig. 3d and Extended Data Fig. 3a). However, by comparing 9 late-stage AMD retinal samples to these same 5 early-stage retinal samples, we observed 939 peaks with significantly decreased intensity, suggesting that the chromatin accessibility changes in the retina occur only during late stages of disease.
In contrast, when we compared 4 RPE samples from early-stage AMD to RPE samples from 8 healthy controls, we observed 5,458 significantly decreased peaks, but observed only 2 significantly decreased peaks when these same early-stage samples were compared to 8 RPE samples from late-stage AMD ( Fig. 3d and Extended Data Fig. 3a).
Likewise, when averaging the intensities of significantly decreased peaks at any stage of AMD progression, we found the striking decrease of chromatin accessibility in the RPE occurred at an earlier disease stage than that observed in the retina (Fig. 3e). This observation fits with the widely accepted theory that changes in RPE function trigger AMD 10 , and suggests that epigenetic changes in RPE cells might be a critical factor that regulates disease onset.

AMD-associated changes in gene regulatory networks
We next sought to determine the functional consequence of the differentially accessible regions (DARs) that were observed in normal and AMD samples. To define statistically significant DARs, we used a linear regression model to take into account the potential effects from other confounding factors such as topographical differences (macula vs. periphery), age, gender, and procurement interval. The model estimated the relative contributions from these factors and the effects of disease stage (normal, early, and late AMD) to variations in peak intensity. Our analysis suggested that the peaks in macular retinas are more likely to be reduced than those in peripheral retinas (Extended Data Fig. 3b). Moreover, a longer procurement interval leads to smaller peaks. Given that the procurement interval of AMD samples is slightly shorter than normal samples (Table   1), it is highly unlikely that an altered procurement interval leads to the decreased peak intensity in AMD samples that is observed in this study.  Table 3). For example, the binding motifs of OTX2 and CRX, factors known to play an important role in controlling gene expression in photoreceptors 11,12 , are enriched in AMD retinas. Moreover, OTX2 showed a significantly decreased footprint in DARs for comparison of late-stage AMD to normal samples (Extended Data Fig. 4b).
This pattern confirmed that chromatin accessibility of OTX2 target sites is decreased in retinal samples with AMD disease, suggesting that reduced target sites binding by retina and RPE-specific TFs play a critical role in AMD pathogenesis.

Genes in DARs show altered expression in AMD
We further tested whether the expression levels of genes associated with DARs are more likely to be altered in AMD. We first checked that DAR-associated genes in both the retina and RPE were highly enriched for genes that were selectively expressed in each tissue (Fig. 4b). Moreover, DAR-associated genes in the retina were substantially more likely to regulate retinal layer lamination and photoreceptor survival, while DARassociated genes in the RPE were more likely to regulate the inflammatory response and apoptosis, which are important biological processes in AMD (Fig. 4b).
Using RNA-Seq data obtained from the patient with differential AMD stages between eyes, we observed that ATAC-Seq peak intensity was highly correlated with gene expression in both the retina and RPE ( Fig. 4c and Extended Data Fig. 4c). In the patient with asymmetric AMD progression, DAR-associated genes were significantly more likely to be downregulated in late-stage AMD relative to early-stage AMD (P = 1.1×10 -5 ,

Epigenetic changes generally occur independently of AMD risk-associated genetic variants
To examine whether the observed changes in chromatin accessibility resulted from AMD-associated genetic variants, we compared the distribution of DARs to that of genetic variants linked to AMD susceptibility by GWAS analysis 5 . For each DAR, we tested whether it overlapped with one or more AMD-associated SNPs identified by GWAS. Interestingly, we observed that very few of AMD-associated SNPs were covered by DARs in genomic location. There are < 0.1% of all DARs in the retina, and < 0.2% of DARs in the RPE, overlapped with AMD-associated SNPs (Extended Data Fig. 4e). Even if we extend a 5kb window in each direction of DARs, the proportion of DARs overlapped with AMD-associated SNPs was still low (< 0.3% for the retina and < 0.4% for the RPE). For comparison, we also examined the fractions of non-DARs and nonpeaks that overlapped with AMD-associated SNPs, and found that the overlap was comparable with those of DAR regions in both retina and RPE. These data imply that the observed differences in chromatin accessibility are unlikely due to local AMD-associated genetic variants.

Cigarette smoke treatment leads to global decreases of chromatin accessibility in human iPSC-derived RPE cells
Since cigarette smoking is the strongest environmental risk factor for AMD 13 , we tested whether cigarette smoke treatment of cultured human RPE cells could trigger similar changes in chromatin accessibility from AMD samples. The cultured iPSCderived RPE cells were examined by flow cytometric analysis of RPE-specific markers (Extended Data Fig. 5a). Furthermore, we performed the ATAC-Seq in the cultured RPE cells and confirmed that RPE cells showed a broadly similar pattern of peak distribution to an average ATAC-Seq profile of all healthy samples from RPE tissue (R = 0.83, Extended Data Fig. 5b). For iPSC-derived RPE cells exposed to cigarette smoke extract, ATAC-Seq profile showed a global decrease in chromatin accessibility (Fig. 4e). More importantly, when comparing cigarette smoke treated RPE cells with RPE tissue from AMD patients, we found that the changes in chromatin accessibility are highly correlated ( Fig. 4f; R = 0.36, P < 10 -20 ). The results indicate that cigarette smoke treatment in RPE cells induces a widespread decrease of chromatin accessibility that is much like that seen in AMD.

Increased HDAC11 expression leads to global decreases in chromatin accessibility
Next, we attempted to identify the genes that could induce similar change in chromatin accessibility observed during AMD progression. To do so, we queried a previously published collection of microarray data from AMD samples 14 , along with our own RNA-Seq analysis, to identify differentially expressed genes that are known to regulate chromatin accessibility. Among all HDAC genes, we found that the class IV histone deacetylase HDAC11 had significantly increased expression in the RPE during early disease stages (Extended Data Table 4 and Fig. 5c). Furthermore, cigarette smoke treatment of iPSC-derived RPE cells also increased HDAC11 expression based on a western-blot analysis (Extended Data Fig. 5d and 5e).
We then tested whether overexpression of HDAC11 induces broad decrease in chromatin accessibility. First, we confirmed that HDAC11 was overexpressed in plasmid transfected RPE cells (Extended Data Fig. 5f). We then examined the chromatin accessibility profile in the monolayers of HDAC11-overexpressing RPE cells, and observed a widespread decrease of chromatin accessibility in HDAC11-overexpressing cells relative to cells transfected with empty vector (Fig. 4g). Interestingly, the changes induced by HDAC11 overexpression were strongly consistent with those changes seen in the RPE of AMD patients (Fig. 4h). These results suggest that HDAC11 overexpression may be partially responsible for the global decreases of chromatin accessibility associated with AMD progression.

Discussion
In this study, we report the first comprehensive analysis of changes in chromatin accessibility in AMD. These changes in chromatin accessibility are seen first in the RPE, and then later in the retina. Likewise, we observed greater changes in the macular retina than in the peripheral retina. These data fit with the typical pattern of AMD progression, where changes in the RPE precede the dysfunction and death of macular photoreceptors 10 .
Several lines of evidences suggested that the changes in chromatin accessibility are not simply due to cell death in the retina or RPE. First, a diverse proportion of cell typespecific genes (e.g. rods, cones and Müller glia) were associated with decreased ATAC-Seq peaks (Extended Data Fig. 5g), suggesting that the decrease of chromatin accessibility is not simply due to photoreceptor death. Second, the genes associated with DARs are enriched for specific cellular functions, suggesting that the regions for reduced chromatin accessibility are selective. Third, both cigarette smoke treatment and HDAC11 overexpression mimic the effect of AMD on chromatin accessibility. Altogether, the results suggest that cell death in AMD is unlikely to cause the observed global reduction in chromatin accessibility.
We hypothesize that global changes in chromatin accessibility may be seen in many other diseases. Indeed, large-scale changes in chromatin accessibility have been previously reported in metastatic cancer, where globally increased chromatin accessibility appears to directly drive disease progression 15 . Cancer cells show generally high levels of metabolic activity, and these changes may partially reflect this fact 16,17 . In contrast, neurodegenerative diseases such as AMD are associated with decreased metabolic activity from early stages onward 18,19 . The resulting reduction in cellular levels of acetyl-CoA, an essential cofactor for histone acetylation, that occurs during disease progression may be a common mechanism that contributes to the observed changes in chromatin accessibility 20 . Our data raise the possibility that global, quantitative reduction in chromatin accessibility may also be observed in other retinal dystrophies, and for neurodegenerative diseases in general. Early-stage AMD was defined by the presence of any RPE pigmentary changes and/or large-size drusen (>125μm diameter). Late-stage AMD was defined by areas of geographic atrophy due to loss of the RPE. In this study, we only included dry AMD and excluded wet AMD. For each eye, we separated the retina and RPE and then obtained a punch (6 mm diameter) of retina or RPE tissue in each of macular and peripheral regions.

Human samples
We obtained paired eyes from some donors and one eye from other donors.

ATAC and RNA sequencing
Biopsy punches of fresh retina and RPE tissues were re-suspended in cold PBS.

Mapping and normalization of ATAC-Seq
After removing adaptors using Trimmomatic 21 , 50bp paired-end ATAC-Seq reads were aligned to the human reference genome (GRCh37/hg19) using Bowtie2 with default parameters 22 .After filtering reads from mitochondrial chromosome M and sex chromosome Y, we included properly paired reads with high mapping quality (MAPQ score > 10, qualified reads) through SAMTools for further analysis 23 . Duplicate reads was removed using Picard tools MarkDuplicates program (http://broadinstitute.github.io/picard/). ATAC-Seq peak regions of each sample were called using MACS2 with parameters -nomodel --shift -100 --extsize 200 24 . The blacklisted regions in human were excluded from peak regions (https://www.encodeproject.org/annotations/ENCSR636HFF/). To get a union set of peaks, we merged the ATAC-Seq peaks of which the distance between proximal ends is less than 10 base pairs. We totally identified 308,019 peaks from retina samples and 208,592 from RPE samples. For each of retina and RPE samples, the fragments were counted across each peak region using HTSeq 25 . We further calculated the normalized fragments (C N ) by dividing the raw fragments (C R ) by the library size (S L ) through the formula:

K-means clustering and multidimensional scaling
The circle plot was performed to visualize genomic ATAC-Seq peaks using Circos 26 .
To compare chromatin accessibility of retina and RPE with other tissues, we mapped chromatin accessibility of 125 cells or tissues measured by DHS-Seq to ATAC-Seq peak regions 27 . DHS-Seq data were downloaded from ENCODE project (https://genome.ucsc.edu/ENCODE/index.html). K-means clustering was used to divide ATAC-Seq peaks into tissue-specific and shared groups. To present a two-dimensional distribution of the retina and RPE samples, we performed multidimensional scaling (MDS), in which all pairwise Euclidean distances were calculated as the distance metric.
The distances in MDS represent the similarity of samples. The analyses were carried in R platform.

Comparison of chromatin accessibility and differentially accessible regions (DARs)
An MA plot (log 2 fold change vs. mean average) was used to visualize the change of chromatin accessibility for all peaks. Meanwhile, the density of peaks was showed across the epigenetic changes. For ATAC-Seq peaks, we accessed the significant change of chromatin accessibility between different groups using edgeR 28

Genomic features and function enrichment of DARs
We used ANNOVAR to annotate ATAC-Seq peaks with the nearest genes 29 . ATAC-Seq peaks were assigned to four categories, including promoter-proximal, exonic, intronic, 3'UTR and intergenic peaks. Promoter is defined as the region within 2kb of the reference transcript start site (TSS) from UCSC database (https://genome.ucsc.edu/). The peaks located in 5'UTR were also taken as promoter-proximal peaks. The peaks located at downstream 2kb of transcript end site were merged to 3'UTR peaks. Gene-proximal peaks include the peaks within 2kb of the gene body.
Using the linear model, we obtained the coefficient of disease stage for each ATAC-Seq peak and nearby gene. We then ranked genes based on the corresponding coefficient of disease stage. If multiple ATAC-Seq peaks were mapped to the same gene, we assigned the lowest coefficient to the gene. Function enrichment of accessibility changeranked genes were further performed using GSEA software 30 . Beside of GSEA database, gene sets of photoreceptor degeneration, retina expressed genes, and RPE expressed genes are from published papers 14,31,32 . Similar analysis of GSEA was performed to identify hot spots where are enriched of DAR-associated genes.

Motifs and footprints of transcription factors associated with DARs
To identify DAR-related transcription factors (TFs), we obtained 1,043 position weight matrices of TF motifs from TRANSFAC database 33

Gene expression analysis
Using RNA-Seq, we measured gene expression of retina and RPE from AMD2 left and right eyes. Using Tophat and Cufflinks with default parameters 36 , raw data were mapped to the hg19 genome assembly, and gene expression FPKM (Fragments Per Kilobase Of Exon Per Million) were obtained. In the study, we also used normal, pre-AMD and dry AMD samples from published expression dataset (GEO accession number: GSE29801) 14 .
If multiple ATAC-Seq peaks mapped to the same gene, we select the highest ATAC-Seq peak among the gene-proximal peaks to associate with gene expression. To compare differential chromatin accessibility and differential mRNA expression, we only included genes with the average expression more than 1 FPKM in left and right eyes of AMD2 patient.

Analysis of GWAS SNPs in DARs
To study the association of DARs with SNPs, we downloaded all SNPs (single nucleotide polymorphisms) from the latest genome-wide association studies (GWAS) performed in AMD 5 . According to genomic location, GWAS SNPs were assigned to corresponding ATAC-Seq peaks. We choose the minimum P value of SNPs in each peak region as the significance of ATAC-Seq peak associated with GWAS SNPs. For ATAC-Seq peak without any SNP, P value of peak was set to 1. Similar analyses were performed for DARs, non-DARs and non-peaks. The flank regions with the same size of peaks were chose as non-peaks. We then calculated the proportion of DARs, non-DARs and non-peaks whose P values are less than a fixed threshold.

Cell type-specific genes involving in DARs
We collected seven sets of genes which are specifically expressed in rods, cones, horizontal cells, bipolar cells, Müller glia, amacrine cells and ganglion cells, which were assembled from previously published data in mouse [37][38][39][40][41][42][43][44] . Using information downloaded from MGI database (http://www.informatics.jax.org/), we mapped these genes to their corresponding human orthologues. We then compared cell type-specific genes to genes in DARs, and identified the number of cell type-specific genes in associated with DARs.

Human RPE derived from induced pluripotent stem cells (iPSC) and plasmid transfection
Human induced pluripotent stem cells (iPSC) were cultured and differentiated into RPE monolayers as previously described 45,46  After washing, gels were incubated in GelCode Blue Stain Reagent (Thermo Fisher Scientific TM ) for 1 hour and then destained by washing in water with several changes over 1-2 hours. Gels were scanned and the densitometric analysis was done by Quantity One software (Bio-Rad Laboratories).

Statistical analysis
All statistical analyses were performed in R platform (https://www.R-project.org/).
Fisher's exact test, Student's t-test, one-way ANOVA, Pearson's correlation coefficient were used to assess the significance. Fold change, P-values and FDR (false discovery rate) were calculated in analysis.

Data availability
The data of ATAC-Seq and RNA-Seq data in this study have been deposited in NCBI's Gene Expression Omnibus (GEO) under accession number GSE99287.
Supplementary Information is available in the online version of the paper.