Whole-exome sequencing revealed exonic profiles of SAM strains

No novel exonic mutations commonly detected among SAMP strains

The Ogg1 p.R304W mutation was previously observed in all of the SAMP strains, but this same mutation was also detected in NZB/N, NFS/N, SJL/J, and NOD/ShiLtJ strains [23]. The Ogg1 gene encodes the enzyme 8-oxoguanine DNA glycosylase, by which oxidatively modified bases are repaired [24, 25]. The methyl-CpG-binding domain protein, encoded by the Mbd4 gene, is also a DNA repair enzyme that is responsible for removing mismatched thymine or uracil within methylated CpG sites [26]. Similar to Ogg1 p.R304W, Mbd4 p.D129N was previously found in normal mice strains including 129P2/OlaHsd, 129S1/SvImJ, 129S5SvEvBrd, DBA/2J, LP/J, NOD/ShiLtJ, and NZO/HlLtJ [21]. It is interesting that all of the SAMP strains as well as the NOD/ShiLtJ strain share these genes that are involved in DNA repair, i.e., Ogg1 and Mbd4. NOD/ShiLtJ is a mouse model of type 1 diabetes, showing a short life span [27, 28]. Nevertheless, we should be careful to conclude that the combination of these mutations regulates the accelerated-senescence phenotype of SAMP, because the short life span of NOD/ShiLtJ is generally attributed to diabetes caused by insulitis.

Unique deleterious mutations identified in each substrain

We also detected 52 novel deleterious mutations in SAMR strains (Additional file 1: Table S11). These results are not surprising, because it has been reported that SAMR strains exhibit several diseases such as non-thymic lymphoma, histiocytic sarcoma, and ovarian cysts [29], although the SAMR strains have been used as control groups against the SAMP strains. Novel deleterious mutations including Fbxl13 p.S734N, Sh3bp5l p.R217W, Tnrc6a p.A278V, and Zkscan2 p.C232X were detected among all of the SAMP strains in addition to being found in several SAMR strains (Additional file 1: Table S12). These mutations may be associated with susceptibility to diseases in SAMP strains as well as in SAMR ones.

Prxp.R167C mutation in SAMP3, SAMP10, and SAMP11 strains

Interestingly, we detected deleterious mutations in several genes that had been earlier reported to be associated with severe genetic disorders in both humans and mice. For example, the Prx p.R167C mutation (SIFT score: 0.01, PolyPhen-2 score: 0.998) was detected in 3 SAMP strains (SAMP3/SlcIdr, SAMP10/TaSlc, and SAMP11/SlcIder; Table 5). The Prx gene encodes periaxin, a protein required for the maintenance of myelin [30]. Because myelin is necessary for the conduction of high-frequency and high-velocity nerve impulses by saltatory conduction, its defects lead to severe neuropathy. In humans, nonsense mutations in periaxin cause an autosomal recessive form of CMT4F (Dejerin-Sottas disease), which is one of the severe hereditary motor and sensory neuropathies [31, 32]. It was also reported that periaxin-knockout mice exhibit peripheral demyelination, mechanical allodynia, and thermal hyperalgesia [33]. The Prx p.R167C mutation is located within the nuclear localization signal (NLS), which is necessary for this protein to be imported into the nucleus from the cytoplasm (Figure 1) [34]. Localization of periaxin in the nucleus is observed in murine embryonic Schwann cells only for a limited period of time [35]. Because the arginine at position167 is highly conserved among mammalian species, this mutation in NLS might disturb the transfer of periaxin into the nucleus, thereby adversely affecting the normal differentiation of Schwann cells.

Ldb3p.R473W mutation in all of SAMP strains except for SAMP8

In the Ldb3 gene, encoding LIM domain-binding protein 3, the p.R467W mutation (SIFT score: 0.02, PolyPhen-2 score: 0.968) was detected in all of the SAMP strains except for SAMP8/TaSlc (Table 5). Ldb3 is a component of the sarcomere Z disk protein complex expressed in cardiac and skeletal muscles, and it is connected to calsarcin-1 and α-actinin [36]. Mutations in the Ldb3 gene are responsible for myofibrillar myopathy and dilated cardiomyopathy in humans [37, 38]. In addition, LDB3 exon 4 is aberrantly spliced in myotonic dystrophy type 1 [39]. Pathological changes in skeletal and cardiac muscles of SAMP strains, however, have not been fully analyzed.

Gja3p.S405P mutation in SAMP3, SAMP6, SAMP10, and SAMP11 strains

We detected the Gja3 p.S405P mutation (SIFT score: 0.09, PolyPhen-2 score: 0.917) in 4 SAMP strains (SAMP3/SlcIdr, SAMP6/TaSlc, SAMP10/TaSlc and SAMP11/SlcIdr; Table 5). Gap junction protein alpha 3, encoded by Gja3, is specifically expressed in the plasma membrane of lens fiber cells to form gap junctions [40]. Gap junctions directly connect the cytoplasm of adjacent cells, and allow various molecules and ions to pass freely between cells, functioning for the maintenance of osmotic and metabolic balance in the avascular lens. A large number of studies have reported the association of mutations of the GJA3 gene with cataract in humans [41, 42].

Il4ra p.S540fs and Zdhhc12p.R112C mutations in the SAMP6 strain

We next analyzed deleterious SNVs unique to a single SAMP strain, i.e., not shared with other SAMP strains. We detected 7 deleterious mutations specific to the SAMP6/TaSlc strain, which has been used as a mouse model of osteoporosis or ulcerative colitis [43, 44]. Among these 7 mutations, we focused on the Il4ra p.S540fs frameshift mutation (Figure 2, Table 6). The Il4ra p.S540fs substitution immediately generates a stop codon at this position. Osteoporosis is caused by bone resorption in excess of bone formation. The differentiation of osteoclasts is promoted by RANKL, a membrane-bound cytokine, as well as by inflammatory cytokines such as TNF-α, IL-1 and IL-6 [45–48]. These proteins are mainly expressed in type 1 T-helper lymphocytes (Th1 cells). IL-4, a Th2 cytokine suppresses the formation of Th1 cells to keep the proper balance between Th1/Th2 cytokines. A functional loss of the Il4ra protein would lead to activation of Th1 cells. Therefore, both osteoporosis and ulcerative colitis in SAMP6 might be explained by activation of Th1 cytokines, which would be induced by the Il4ra p.S540fs frameshift mutation.

The Zdhhc12 p.R112C mutation (SIFT score: 0.000, PolyPhen-2 score: 0.999) was also detected uniquely in SAMP6 (Table 6). Zinc-finger DHHC domain-containing protein 12, encoded by the Zdhhc12 gene, has a predicted DHHC cysteine-rich palmitoyl acyltransferase domain [49]. Several gene mutations in the Zdhhc family have been implicated in human diseases and abnormal phenotypes of mice. Remarkably, Zdhhc13-truncated mutant mice develop alopecia, osteoporosis, and systemic amyloidosis [50]; and the osteoporotic phenotype can be explained by the finding that protein palmitoylation regulates osteoblast differentiation through bone morphogenesis protein (BMP)-induced Osterix expression [51]. Thus we speculate that Zdhhc12 p.R112C mutation might contribute to the osteoporotic phenotype in SAMP6.

Aifm3p.K582N mutation specific to SAMP8

Five detected deleterious mutations were unique to the SAMP8/TaSlc strain, which show deficits in learning and memory, emotional disorder, and abnormal circadian rhythm at early ages (Table 7) [52, 53]. It is remarkable that the K582N mutation in the Aifm3 gene (SIFT score: 0.01, PolyPhen-2 score: 0.879), encoding apoptosis-inducing factor mitochondrion-associated protein 3, was detected in SAMP8/TaSlc. Although the function of Aifm3 has not been fully elucidated, it has been reported that Aifm3 shares 35% homology with Aifm1 and that overexpression of Aifm3 induces apoptosis in HEK 293 cells [54]. Because the lysine at 582 in Aifm3 is highly conserved among mammalian species, this p.K582N mutation therefore may alter the function of Aifm3, contributing to the mitochondrial dysfunction in SAMP8 mice.

Whole-exome sequencing identified new candidate mutations responsible for age-related phenotypes in SAMP strains

In the present study, we identified the entire spectrum of the SNVs in 6 SAMP and 3 SAMR strains by whole-exome sequencing. We summarized the candidate mutations regulating various pathogenic phenotypes in SAMP strains in Figure 3. Our study has clarified that several disease-causing mutations were common among multiple SAMP strains. Two of these mutations, Ogg1 p.R304W and Mbd4 p.D129N, were common among all SAMP strains and would be involved in the susceptibility to diseases via defects in DNA repair. In all SAMP strains except SAMP8/TaSlc, we detected a p.R473W missense mutation in the Ldb3 gene, which has been associated with myofibrillar myopathy. In 3 SAMP strains (SAMP3/SlcIdr, SAMP10/TaSlc, and SAMP11/SlcIdr), we identified a p.R167C missense mutation in the Prx gene, which has been linked to hereditary motor and sensory neuropathy. In 4 SAMP strains (SAMP3/SlcIdr, SAMP6/TaSlc, SAMP10/TaSlc, and SAMP11/SlcIdr), we detected a p.S405P missense mutation in the Gja3 gene, which is a cause of cataract. As the candidate gene mutations responsible for strain-specific phenotypes, we detected 24 deleterious mutations specific to a single SAMP strain, including the Il4ra p.S540fs frameshift mutation in SAMP6/TaSlc, which is used as a model for osteoporosis, and the Aifm3 p.K582N mutation in SAMP8/TaSlc mice, which display deficits in learning and memory and mitochondrial dysfunction.

We detected Ogg1 p.R304W and Mbd4 p.D129N deleterious mutations, which were common to all of the SAMP strains, but absent in the SAMR and AKR/J strains; although these 2 mutations were also detected in other mouse strains. It was already investigated as to whether a defect in Ogg1 protein would affect the life spans in SAMP strains. Mori et al. reported that hybrid mice with the homozygous mutation in Ogg1 p.R304W exhibited a complete loss of the glycosylase activity as well as a higher level of 8-oxoguanine in their hepatic nuclear DNA [23]. However, the average life span of the SAMP1×B10.BR hybrid was not different among the mice homozygous, heterozygous or nullzygous (B10.BL allele) for the SAMP1 allele. Moreover, NZB/N, NFS/N, SJL/J, and NOD/ShiLtJ also have the Ogg1 p.R304W mutation. These results suggest that Ogg1 p.R304W alone is not sufficient to cause accelerated senescence and a short life span. We assume that the combination of Ogg1 p.R304W and Mbd4 p.D129N causes accelerated senescence. Both mutations were detected in the NOD/ShiLtJ strain, which is a type 1 diabetes model [27]. Although NOD/ShiLtJ mice may live for only 6 to 8 months due to diabetes under normal food and water conditions [28], we cannot predicate these mutations to be essential for the accelerated-senescence phenotype of SAMP because the cause of death is different between SAMP and NOD/ShiLtJ strains. Nevertheless, mouse strains that possess Ogg1 p.R304W mutation are known for their pathologic phenotypes: NZB/N for autoimmune hemolytic anemia; SJL/J for reticulum cell sarcomas, in addition to NOD/ShiLtJ for type 1 diabetes [55, 56]. Somatic mutations have been implicated in various diseases, and the accumulation of such mutations is one of the most accepted theories to explain aging. The Ogg1 p.R304W mutation might partly contribute to the phenotypes of these mouse strains as well as to the accelerated-senescence phenotypes of SAMP strains.

In several SAMP strains, missense mutations were detected in the Prx, Ldb3, and Gja3 genes, which mutations have been found in various human degenerative diseases. The pathogenesis of myofibrillar myopathy and peripheral neurodegeneration has not been fully analyzed in SAMP strains. Age-related muscle atrophy and a decline in peripheral neuronal function are assumed to be a common phenomenon that probably occurs in the course of the senescence process [57, 58]. Nevertheless, genetic susceptibility to degeneration of skeletal muscle and peripheral neurons may be different among SAMP strains. Prx p.R167C and Ldb3 p.R473W mutations possibly contributed to the degenerative phenotypes of 3 of the SAMP strains in the course of the senescence process.

In the present study, the Gja3 p.S405P mutation was detected in 4 SAMP strains (SAMP3, SAMP6, SAMP10, and SAMP11), among which only the SAMP3 strain is reported to develop cataract [59]. As a lack of reports of cataract in SAMP6, SAMP10, and SAMP11 does not indicate the actual lack of cataract, careful ophthalmologic examinations for cataracts in these 3 other SAMP strains may reveal a pathogenic association. Alternatively, because it is suggested that the pathogenic mechanism underlying the development of cataract in SAMP strains is different from that of murine hereditary cataract, which is generally regulated by single-gene mutations [60], the Gja3 mutation alone may not be sufficient to cause cataract. The SAMP3 strain may have additional mutations besides the Gja3 p.S405P mutation that are responsible for cataract.

The Il4ra p.S540fs frameshift mutation can explain the osteoporosis observed in the SAMP6 strain from the viewpoint of osteo-immunology. It is known that IL-4 signaling inhibits osteoclast differentiation by suppressing Th1 cytokines such as RANKL, TNF-α, and IL-1. In fact, Il4 gene knockout mice are sensitive to RANKL-induced bone resorption [61]. A defect in Il4ra might thus enhance osteoclast differentiation due to dysregulation of Th1 cytokines. The Il4ra p.S540fs frameshift mutation can also explain the ulcerative colitis found in the SAMP6 strain. Although the true cause of ulcerative colitis remains unknown, abnormalities of the immune system are possibly related to its pathogenesis. Particularly, a high level of TNF-α was proposed to play an important role in disease progression [62]. In SAMP6 mice, it is expected that up-regulation of TNF-α expression in the colon would occur due to activation of Th1 cells. Thus, both of these pathogenic phenotypes, osteoporosis and ulcerative colitis, may be ascribable to the defect in Il4ra in SAMP6.

We also detected the Aifm3 p.K582N mutation in the SAMP8/TaSlc mice, which display deficits in learning and memory. High oxidative stress derived from brain mitochondrial dysfunction is thought to be one of the causes of age-related neurodegeneration in SAMP8 animals. Actually, decreased activities of NADH-cytochrome c reductase are observed even in 4-week-old SAMP8 mice, suggesting crucial defects in maintenance of the respiratory chain [63]. Aifm3 is likely to be related to mitochondrial maintenance, because it induces apoptosis in vitro and has an oxidoreductase domain, as is the case for Aifm1 [54], which plays roles in maintenance of the mitochondrial respiratory chain [64]. However, the actual roles of Aifm3 in apoptosis in the senescence process and the actual substrates of the oxidoreductase remain unknown. Further investigations are necessary to examine whether the Aifm3 p.K582N contributes to deficits in learning and memory via dysfunction of brain mitochondria.

Overall, it seems that the combinations of different disease-causing mutations specific to each strain cause various degenerative diseases, which combinations are a cause of short life spans of SAMP strains as far as focusing on the mutations of the coding regions is concerned. Actually, it was reported earlier that the life spans of SAMP strains are susceptible to environmental conditions [65]. These observations may be ascribable to the multifactorial nature of the short life span of SAMP unlike other progeroid mice whose life span is regulated by single gene mutation. de Magalhaes JP et al. also reported that the Gompertz mortality curve of the SAMP was not different from that of the SAMR prior to age 1 year despite the difference in age when 50% of mice died, suggesting that the life spans of the SAMP strains may not be related to aging per se[66]. Nevertheless, we think that it is premature to conclude that SAMP strains are degenerative disease models rather than accelerated-senescence models because in vitro studies have shown that primary-cultured cells from several SAMP strains show accelerated senescence and higher oxidative stress and mitochondrial dysfunction than the SAMR1 strain [67–69].

Limitation of the present study

Whole-exome sequencing using 50-bp single-end reads on the SOLiD4 platform is able to detect only single or 2-base nucleotide variations and insertion/deletion. Because accelerated senescence and the various pathogenic phenotypes may not be explained completely by the nucleotide substitutions in the coding regions, we cannot ignore the possibility that other types of genetic variations are also involved in common accelerated-senescence phenotypes of SAMP strain. Fairfield et al. succeeded in identifying causative mutations in several ENU-induced mutants by exome sequencing, but failed to do so in several spontaneous disease models [70]. They suggested that mutations responsible for spontaneous disease models might reside in the non-coding regions. Actually, it has been proven that most of the non-coding regions have some biochemical functions [71].

Carter et al. reported a 15-bp insertion mutation in the Fgf1 gene in SAMP10 [7], suggesting the involvement of a small structural variation in an exon of this gene. A long-read sequencing platform, which can generate over 200-bp fragments, would be required to detect them. It has been suggested that not only small structural variations in exons, but also large genomic structural variations such as copy number variations and gene translocations, contribute to the complex traits of humans [72]. Furthermore, because complementary RNA probes are designed based on reference genome sequences, we were limited to find variants in comparison with the reference sequence. De novo assembly by whole genome sequencing or mate-pair library sequencing and comparative genomic hybridization (CGH) array analysis should be performed to detect these sequence variations.

In present study, we focused on only novel deleterious mutations that could be predicted by SIFT and PolyPhen-2. Although these bioinformatics tools are useful to narrow down the candidate mutations, a recent study indicated that SIFT and PolyPhen-2 show 63 and 79% correct prediction rates, respectively [73]. In the future, functional analyses should be conducted to confirm whether the mutations that were predicted to be deleterious in the present study really affect the functions of these genes.

DNA extraction

Genomic DNA was extracted from the livers of 11 mouse strains, i.e., SAMP1/SkuSlc, SAMP3/SlcIdr, SAMP6/TaSlc, SAMP8/TaSlc, SAMP10/TaSlc, SAMP11/SlcIdr, SAMR1/SlcIdr, SAMR1/TaSlc, SAMR3B/SlcIdr, AKR/J and C57BL/6J strains. RNase treatment was performed to obtain a high-quality DNA library. All experimental procedures using laboratory animals were approved by the Animal Care and Use Committee of the Tokyo Metropolitan Institute of Gerontology, the Institute for Developmental Research of the Aichi Human Service Center, and by Shinshu University School of Medicine.

Targeted capture and next-generation sequencing

Target enrichment was performed by use of a SureSelect^XT Mouse All Exon kit (Agilent Technologies, Santa Clara, California, US) optimized for the ABI SOLiD system and 3 μg of genomic DNA according to the manufacturer’s protocol. The kit is designed to enrich for 221,784 exons within 24,306 genes covering a total of 49.6 Mb genomic sequences. DNA was sheared by acoustic fragmentation (Covaris, Woburn, Massachusetts, US) and purified with an Agencourt AMPure XP kit (Beckman Coulter, Brea, California, US). The quality of the fragmentation and purification was assessed with an Agilent 2100 Bioanalyzer. The fragment ends were repaired and adaptors were ligated to the fragments (Agilent). The modified DNA library was purified by using the Agencourt AMPure XP kit, and amplified by PCR and captured by hybridization to biotinylated RNA library baits (Agilent). Captured DNA was purified with streptavidin-coated magnetic Dynal beads (Life Technologies, Carlsbad, California, US) and amplified with Barcoding Primer. The prepared exome library was pooled and subjected to emulsion PCR and sequenced on the SOLiD4 (Life Technologies) as single-end 50-bp reads. For each sample, 1 quad of a SOLiD sequencing slide was used.

Read mapping and variant analysis

Sequence reads were mapped to the reference mouse genome (UCSC mm9, NCBI build 37) by using Bioscope software version 1.3 (Life Technologies), which utilizes an iterative mapping approach. After removal of low-quality and duplicate reads, single nucleotide variants (SNVs) were detected with Avadis NGS software version1.3 (Strand Life Sciences, Bangalore, Karnataka, India). Avadis NGS performs SNV identification via an adapted version of the MAQ algorithm, which calculates the probability that the consensus genotype is incorrect by using a Bayesian statistical model with mapping quality, base quality and ploidy taken into consideration. We established criteria for SNV detection as a read coverage ≥ 2, and other parameters were set as default values. Detected SNVs were annotated for extracting non-synonymous and homozygous SNVs by using the Avadis NGS with UCSC transcript annotation. Moreover, we extracted novel SNVs by comparison with NCBI dbSNP build 128 and SNV data of 17 inbred strains of laboratory mice obtained by whole-genome sequencing. We compared filtered SNVs among all strains to explore the mutations that were commonly present among the SAMP strains but absent in the SAMR strains, AKR/J strain, and C57BL/6J strain. The unique mutations of each strain were also selected.

Interpretation of novel missense SNVs

To predict whether the candidate SNVs would have deleterious effects or not, we used 2 software programs, i.e., Sorting Intolerant from Tolerant amino acid substitutions (SIFT; J. Craig Venter Institute, San Diego, California, US, http://sift.jcvi.org/) and Polymorphism Phenotyping v2 (PolyPhen-2; Harvard University, Cambridge, Massachusetts, US, http://genetics.bwh.harvard.edu/pph2/). SIFT uses sequence homology to predict amino acid substitutions that will affect protein function, thus contributing to a disease [74]. SIFT predicts substitutions with a score less than 0.05 as being “deleterious“ (Range: 0 to 1). PolyPhen-2 takes into account the physicochemical characteristics of the wild-type and mutated amino acid residue and the consequence of the amino acid change for the structural properties of the protein in addition to evolutional conservation [75]. PolyPhen-2 generates a different scale of reported scores, with the corresponding predictions being “probably damaging” with a score larger than 0.85, “possibly damaging” with a score between 0.85 and 0.15,” and “benign” with a score less than 0.15. Because PolyPhen-2 considers only human protein sequences, the mouse SNVs were investigated in the context of human protein sequences.

Mutation validation

Validating the candidate SNVs was performed by using the standard Sanger sequencing approach. Primers were designed to surround candidate SNVs by using Primer 3 version 4.0, and custom DNA oligos were ordered (Life Technologies; Operon Biotechnologies, Tokyo, Japan). Primer sequences are shown in Additional file 1: Tables S1-S2. PCR reactions were carried out in 10-μl reaction mixtures containing a 0.5 μM concentration of each primer, 0.2 mM dNTPs, 0.25U Ex Taq DNA Polymerase Hot-Start Version, 1.0 μl 10×Ex Taq Buffer (Takara Bio, Shiga, Japan), and 1 μl of extracted DNA. The amplification conditions were 1 cycle at 96°C for 5 min of denaturation, 40 cycles of 94°C for 30 s, 55-68°C for 45 s of annealing in proportion to the Tm value of each primer, and extension at 72°C for 45 s, followed by a final extension at 72°C for 10 min. PCR products were purified by using a MultiScreen_HTS PCR 96-Well Plate (Millipore, Billerica, Massachusetts, US) for sequences. DNA templates were subjected to the sequencing reactions by using a BigDye Terminator version 3.1 Cycle Sequencing Kit (Life Technologies). The sequencing reaction solution contained 4 μl BigDye Terminator v3.1, 0.32 μM M13 forward primer, 1.75 μl 5×Sequence Buffer, and 2.0 μl PCR product in a final volume of 10 μl. PCR conditions were 1 cycle at 94°C of denaturation, 25 cycles of 94°C for 10 s, 50°C for 15 s and 3 min at 60°C, followed by cleaning of the reaction products by ethanol precipitation. The capillary electrophoresis sequencing was performed by using an ABI Prism 3130xl Genetic Analyzer (Life Technologies), and sequence data were analyzed with Sequencher version 4.2.2 (Gene Codes, Ann Arbor, Michigan, US).

Gene Ontology enrichment analysis

Gene Ontology enrichment analysis (GO analysis) was performed by using WebGestalt (http://bioinfo.vanderbilt.edu/webgestalt/) [76]. The obtained p-values were adjusted by Benjamini-Hochberg multiple testing, and the significant level was established at p<0.05.

Multiple alignment

Multiple sequence alignment was performed by using the Clustal Omega program on the UniProt website (http://www.uniprot.org/) [77].

Data access

Exome data were deposited in DDBJ Sequence Read Archive (BioProject Accession Number: PRJDB37).