The correlation of copy number variations with longevity in a genome-wide association study of Han Chinese

The numbers of CNVs increase in the genomes of long-lived individuals

To identify the CNVs associated with longevity in Chinese population, we recruited 1950 long-lived individuals called “cases”, and 2057 middle-aged Chinese as “controls. These individuals were further separated into four cohorts: 1000 long-lived and 1215 middle-aged who lived in the north of China, and the rest of 950 long-lived and 842 middle-aged who resided in the south. The two Northern cohorts were also called “discovery cohorts”, and the two Southern cohorts were referred as “replicate cohorts”.

Using the Principal Component Analysis (PCA) method, no extra sub-cluster was found between case and control comparison (Supplementary Figure 1). Our data obtained from the Northern cohorts revealed a significantly increased number of CNVs among long-lived individuals compared to the controls (6.94 ± 0.35 vs 5.96 ± 0.27, p = 0.027). The similar results were found in the Southern replicate cohorts (5.96 ± 0.27 vs 4.60 ± 0.17; p = 0.001). Total identified CNVs were summarized in detail in Table 1.

In all participating subjects including both cases and the controls, we identified 10046 deletions and 6932 duplications in total. We found that the numbers of CNVs increased significantly in older ages (Spearman rho = 0.386, p = 0.002; Figure 1), particularly, the deletions increased much more than the duplications (Spearman rho = 0.356, p = 0.004). On average, the centenarians (101.51 ± 0.07 years of age) contained 4.15 ± 0.20 deletions and 2.30 ± 0.10 duplications, showing a significant increase (p = 0.001) compared to the middle-aged (48.22 ± 0.16 years of age), who had 3.25 ± 0.16 deletions and 2.14 ± 0.07 duplications. The Spearman correlation was also calculated to show the significance between ages and total added lengths of CNVs (Spearman rho = 0.31, p = 0.017); long-lived people had 508.54 ± 19.8 kb of total CNVs, which was significantly longer than 453.69 ± 17.25 kb in middle-aged controls (p = 0.024). In summary, we have shown that the CNV numbers, especially the deletion numbers, increased significantly in the genomes of long-lived people.

Figure 1. CNV burden in different age groups. (a) The numbers of CNVs in different age groups. (b) The added lengths of CNVs in different age groups. Triangles represent the numbers of the CNVs. The areas in each histogram represent the percentage of different lengths. * p<0.05 by t-test.

Identification of CNVRs associated with longevity

We identified the CNVRs associated with longevity, according to a previously published method [17]. These regions included 62 deletions and 5 duplications identified from the Northern cohorts, 21 deletions and 6 duplications in the Southern cohorts, and 99 deletions and 27 duplications in the combined samples (Figure 2). Among them, we identified eleven CNVs (2 duplications and 9 deletions) with case frequency > 1%, p value < 3.97×10^-4 and length > 10k in long-lived individuals from the north as well as the north and south combined samples (Table 2, Supplementary Table 1). These CNVRs have been previously deposited in the database of genomic variants (DGV, http://dgv.tcag.ca/dgv/app/home) through other studies unrelated to longevity (Supplementary Figure 2). In order to confirm the accuracy of genotyping and CNV detecting methods, top six CNVRs (with p-value less than 10^-5) were chose to do the quantitative PCR. (Supplementary Figure 3) and all these six CNVRs were validated accurate.

Figure 2. The distribution of CNVs among 4007 individuals. The first circle indicates the positions of chromosomal bands; the second is a histogram representing the frequencies of CNVs in long-lived individuals (red: deletions; blue: duplications; height: CNV frequencies); the third is also a histogram showing the CNV frequencies in the long-lived (orange) and middle-aged (orange) individuals. Frequency: orange and green, >1%; grey, <1%; red, the validated CNVRs by qPCR. The heatmap shows the p-values of CNVRs: orange for the long-lived, green for the middle-aged, and the color gradience towards dark indicates decreasing p-values. The text presents the names of the CNVRs identified in this study.

Functional analysis of eleven CNVRs

To determine if our identified CNVRs can affect biological function, we analyzed these loci for encoded genes. We found several genes within or near 11 CNVRs including Zinc Finger Protein 716 (ZNF716), Zinc Finger Protein 208 (ZNF208), Zinc Finger Protein 257 (ZNF257), Zic Family Member 5 (ZIC5), Zic Family Member 2 (ZIC2), Pleckstrin Homology Like Domain Family A Member 2 (PHLDA2), Solute Carrier Family 22 Member 18 (SLC22A18), Family With Sequence Similarity 53 Member A (FAM53A), Peptidyl-TRNA Hydrolase 1 Homolog (PTRH1), SH2 Domain Containing 3C (SH2D3C), Torsin Family 2 Member A (TOR2A), Tetratricopeptide Repeat Domain 16 (TTC16), Laminin Subunit Alpha 5 (LAMA5), Rho Guanine Nucleotide Exchange Factor 10 (ARHGEF10), Phosphatidylinositol Glycan Anchor Biosynthesis Class P (PIGP), Tetratricopeptide Repeat Domain 3 (TTC3), Tetratricopeptide Repeat Domain 6 (TTC6), Forkhead Box A1 (FOXA1) and Nuclear Receptor Subfamily 2 Group F Member 2 (NR2F2) (Table 2).

Some of these genes have previously been linked to longevity. For example, it has been shown that the rs4925386 in LAMA5 gene in 20q13.33 belongs to laminin alpha family. The age-stratified analyses showed that the rs4925386-T allele was positively associated with longevity (p = 0.001) [18]. The CNVR in 9q34.1 encodes genes SH2D3C and TOR2A. SH2D3C is a signaling adapter protein required for T cell activation, and a decrease in SH2D3C might cause B cell dysfunction and impaired immune function, which could affect lifespan [19]. In addition, the CNVR in 14q21.1 encodes the gene, FOXA1; upregulation of FOXA1 has been shown to decrease diet-restriction-induced longevity in C. elegans [20]. Additionally, NR2F2, a gene about 14kb away from 15q26.2 (p = 8.43×10^-5, OR = 4.57) showed higher expression in older samples and was related to vascular development during oxidative stress-induced cellular senescence; vascular-related disease could be used as a marker for aging people [21].

To further understand whether these CNVRs affect specific pathways regulating aging process, The Database for Annotation, Visualization and Integrated Discovery (DAVID), Gene Ontology (GO), and Functional Enrichment analysis tool (FunRich) [22–24] were used for the enriched pathway analyses (Table 3). Six pathways were found to be enriched using FunRich including FOXA1 and FOXA transcription factor networks, which were most directly related to the regulation of longevity.

To determine if the difference in CNVR could be found between the nonagenarians (90-99) and centenarians (>100), we conducted the CNVR analyses in these two groups separately We detected four new CNVRs that were unique for the nonagenarians, and six new CNVRs only in the centenarians (Supplementary Table 2).

The variants in the CNVRs and eQTL

To determine if SNPs could alter gene expression in CNVRs, we investigated our 11 CNVRs using the Haploreg database [25]. Four of them were found to contain multiple SNPs, defining ≥10 eQTLs in tissue database that might affect gene function (Supplementary Table 4). These CNVRs are the deletions in 21q22.13, 20q13.33, and 11p15.4, and a duplication in 19p12. The eQTLs in the 21q22.13 deletion region could alter the expressions of the PIGP and TTC3 genes primarily in brain tissue. According to the annotation database Genecards (http://www.genecards.org/), PIGP resides in the protein complex that catalyzes the transfer of N-acetylglucosamine from UDP-GlcNAc to phosphatidylinositol, the first step of the glycosylphosphatidylinositol (GPI) biosynthesis. TTC3 is an ubiquitin-protein ligase that mediates the ubiquitination and subsequent degradation of phosphorylated AKT proteins. Several SNPs in the 20q13.33 deletion might affect the LAMA5 and CABLES2 gene expression. LAMA5 has been shown previously to affect cellular adhesion, migration, and organization. Two SNPs, rs1661052 and rs450244, in the 11p15.4 deletion region, showed strong cis-eQTL characteristics and might affect the SLC22A18 expression and its function in transporting organic cations. Finally, the SNPs located in the 19p12 duplication region can affect the expression of the ZNF gene family member such as ZNF208 and ZNF 257.

Commonly shared CNVRs between Han Chinese and other ethnicities

To determine if our identified CNVRs in Han Chinese were also found in other ethnicities, we compared our data to a Danish and an American study. The Danish study investigated 603 nonagenarians or long-lived individuals of 90.0–102.5 years of age or 96.9 on average [6], while the U.S. ‘wellderly’ healthy aging cohort included 1,354 individuals from 80 to 105 years old or 84.2 on average [26].

Among the 11 CNVRs, the duplication in 7p11.2 (Chr7: 57787663- 57839005) was interesting. The Danish study identified a similar duplication in 7p11.2 (Chr7: 57208666 - 57882950) with a population frequency of 1%. Another similar duplication was also found in the U. S. ‘wellderly’ with a population frequency of 0.2% and a 19% overlap (Table 4). To further determine the distribution of this duplication in different ethnicities, we examined the genomes in the 1000 Genomes Project, which contained randomly picked adults older than 18 years old from healthy populations of multiple ethnics. A similar duplication (Chr7: 57729553 - 57807369) was found in the phase 3 database with low frequencies or no distribution across all populations: 0.1% in European, 0.8% in African, and none in East Asian, South Asian, and mixed American (an admixed population of Americans). In contrast, the frequency of this duplication was significantly higher in long-lived individuals: 15.89%, 1%, and 0.2% in the Han Chinese, Danish, and American population, respectively. Thus, these studies strongly supported that the duplication in 7p11.2 associated with longevity.

Our three other CNVRs were found to be comparable with the CNVRs identified in the U.S. or Danish study (Table 4). The deletion in 20q13.33 (Chr20: 60872280 - 60909518) partially overlapped with the Danish CNVR (Chr20: 60872280 - 60909518). Its population frequencies were 1.07%, 0.33%, and 0 in Han Chinese, Danish, and the 1000 genomes database, respectively. The other two CNVRs consisted of one duplication in 19p12 (Chr19: 22145147 – 22231026) and one deletion in 8p23 (Chr8: 1789006 – 1796964). Both shared 9% and 12% overlapping regions with a similar duplication and deletion in the ‘wellderly’, respectively. The population frequencies in Han Chinese, Ad Mixed American, and the 1000 genomes database for the duplication CNVR were 2.36%, 0.1%, and 0, and 1.02%, 0.59%, and 0 for the deletion CNVR, correlatively These results also indicated that the CNV frequencies were higher in long-lived population than those represented by the 1000 Genomes Project, which composed of randomly picked adults.

The population studied

Our study included 1950 long-lived cases and 2057 middle-age controls from the Chinese Longitudinal Healthy Longevity Survey (CLHLS), which represented 22 provinces of China, following the method described previously [5,41; http://web5.pku.edu.cn/ageing/html/datadownload.html).

These people were divided, based on their geographical locations, into two population cohorts, representing the discovery and replicate samples. The former was composed of 1000 long-lived cases and 1215 middle-age controls from 11 provinces in the north of China, while the latter consisted of 950 long-lived and 842 middle-aged individuals from 11 provinces in the south of China.

The whole genome genotyping

Genotyping was carried out following the standard protocols of the Illumina HumanOmniZhongHua-8 BeadChips (Illumina Inc.). A beadchip contains 900,015 SNPs markers with a mean and median spacing of 3.3 kb and1.7kb, respectively. Quality control (QC) was conducted using plink 1.06 and call rate > 0.95; 339 samples from the replicate was excluded. Principal component analysis (PCA) was also subsequently performed using plink [42].

CNV detection

CNVs were detected using the software PennCNV 1.04 [43]. In the software, a hidden Markov model (HMM) determined the CNVs with the total signal intensity Log R ratio (LRR) and B allele frequency (BAF) for each SNP marker. The LRR and BAF values can be downloaded from the Illumina Genomestudio2011. The population frequency of the B allele (PFB) was calculated using the BAF values collected from the discovery or replicate study. The GC file was obtained from UCSC genome bioinformatics (http://hgdownload.cse.ucsc.edu/goldenPath/hg19/database/gc5Base.txt.gz). The final CNVs were obtained after excluding those that met these conditions: LRR SD (standard deviation) > 0.3, GC wave factor > 0.05, CNV number >100, CNV length < 10kb, consecutive SNPs < 10, and confidence <10.

CNVs encoded gene annotation and pathway analysis

The centenarian specific CNVs were annotated using human genome reference version 19 (hg19) within the range of 100kb. Genome ontology (GO) enrichment was performed using Go Ontology (http://geneontology.org/) [23]. Functional enrichment and interaction network analysis of genes and proteins were performed using The Database for Annotation, Visualization and Integrated Discovery (DAVID) (https://david.ncifcrf.gov/home.jsp) [24] and Functional Enrichment analysis tool (FunRich) 3.1.3 [25].

Determination of the CNVRs associated with longevity

We compared the long-lived to the middle-aged individuals to identify CNVs associated with longevity using the CNVR searching software, ParseCNV20 [17]. The CNVs were called based on p values calculated by Fisher’s exact test. The SNPs in the CNVRs were also compared among the old and middle-aged groups with Fisher’s test.

The identified CNVRs were then annotated based on the hg19 reference genome. Those CNVs overlapping with telomeres were discarded.

We then compared our identified CNVs with a Danish study (personal communications with Dr. Nygaard) and a U.S. study (personal communications with Dr. Ali Torkamani).

The validations of CNVRs

Quantitative polymerase chain reaction (qPCR) was used to validate CNVRs. Considering the long lived DNA are valuable, limited samples could be used for the validation. In addition, we still need samples for proceeding the next stage of GWAS experiments, so the top six CNVRs (with p-value less than 10^-5) were chose to do the quantitative PCR and confirm the accuracy of genotyping and CNV detecting methods. Primers sets were designed for each of selected CNVRs (Supplementary Table 3). The qPCR conditions followed the ABI Step One (Life technologies Inc.) in 96-well plates: 95°C for 10 min to denature DNA samples, followed by 40 cycles (95°C for 15 s and then 58°C for 1 min). Each sample was repeated three times. The relative copy number (RCN) was calculated using CT values (2^–∆∆Ct). RCN was then adjusted by setting the no-copy-number-change as “1”. Larger than “1.2” was regarded as duplication, while smaller than “0.8” indicated deletion.

Expression Quantitative Trait Loci (eQTL) analysis

Haploreg V4.1, a web-based tool for annotating noncoding genomic regions for variations was used to examine the effects of SNPs on gene expression and determine eQTLs (http://archive.broadinstitute.org/mammals/haploreg/haploreg.php). These eQTLs from Haploreg were also identified from the GTEx Project (https://www.gtexportal.org/home/) and Geuvadis Data Browser (https://www.ebi.ac.uk/Tools/geuvadis-das/).

Statistical analysis of CNV burden

To determine if the CNV numbers or total added lengths of CNVs were statistically different between the long-lived and middle-aged cohorts, a Bonferroni correction independent students’ t test was applied. The Spearman correlation was calculated between CNV numbers and CNV lengths. Since genetic association studies are often confronted by the problem of complex population structure that can result in false positive results, in this study, we used an adjusted p-value of 3.97×10^-4 (0.05/126), based on multiple testing, as a cutoff to identify statistically significant in long-lived specific CNVRs. The t test analyses were performed using SPSS Statistics 19.0 software (SPSS, Inc.) and p<0.05 was regarded as significant. Data correlation analyses were performed using a R software from the R Project for Statistical Computing (3.1.4). Figure 1 and supplement figures were constructed with Origin 8.5, a graphing software (Origin Lab Inc.), and Figure 2 was by Circos 0.69, a software package for visualizing genomic data and information in circular layout (http://circos.ca/).