Transmission of human mtDNA heteroplasmy in the Genome of the Netherlands families: support for a variable-size bottleneck

Although previous studies have documented a bottleneck in the transmission of mtDNA genomes from mothers to offspring, several aspects remain unclear, including the size and nature of the bottleneck. Here, we analyze the dynamics of mtDNA heteroplasmy transmission in the Genomes of the Netherlands (GoNL) data, which consists of complete mtDNA genome sequences from 228 trios, eight dizygotic (DZ) twin quartets, and 10 monozygotic (MZ) twin quartets. Using a minor allele frequency (MAF) threshold of 2%, we identified 189 heteroplasmies in the trio mothers, of which 59% were transmitted to offspring, and 159 heteroplasmies in the trio offspring, of which 70% were inherited from the mothers. MZ twin pairs exhibited greater similarity in MAF at heteroplasmic sites than DZ twin pairs, suggesting that the heteroplasmy MAF in the oocyte is the major determinant of the heteroplasmy MAF in the offspring. We used a likelihood method to estimate the effective number of mtDNA genomes transmitted to offspring under different bottleneck models; a variable bottleneck size model provided the best fit to the data, with an estimated mean of nine individual mtDNA genomes transmitted. We also found evidence for negative selection during transmission against novel heteroplasmies (in which the minor allele has never been observed in polymorphism data). These novel heteroplasmies are enhanced for tRNA and rRNA genes, and mutations associated with mtDNA diseases frequently occur in these genes. Our results thus suggest that the female germ line is able to recognize and select against deleterious heteroplasmies.Heteroplasmy (intra-individual variation) in mitochondrial DNA (mtDNA) plays an important role in mtDNA-related diseases and has also been implicated in aging and cancer (Greaves et al. 2012; Wallace 2012; Chinnery and Hudson 2013; Lombès et al. 2014). Most mtDNA mutations that cause diseases due to defects in mitochondrial function exist as heteroplasmies and only cause disease symptoms when the frequency of the mutant allele exceeds a particular threshold (Wallace and Chalkia 2013). Below this threshold, individuals are asymptomatic, presumably because there are sufficient functional mitochondria for normal metabolism. Changes in the frequency of pathogenic mutations during the transmission of heteroplasmies from mothers to offspring can thus play an important role in the disease risk of the offspring. However, most of our knowledge concerning the dynamics of heteroplasmy transmission comes from studies of pathogenic mutations (Monnot et al. 2011; Shen et al. 2012; de Laat et al. 2013; Wallace and Chalkia 2013), which in blood have been shown to decrease over time and hence may not accurately reflect the overall level of such pathogenic mutations within an individual (Poulton and Morten 1993; ‘t Hart et al. 1996; Rahman et al. 2001; Rajasimha et al. 2008). Mouse models have also been utilized (Cree et al. 2008; Fan et al. 2008; Freyer et al. 2012; Ross et al. 2013), but to date, there have been only a few studies of normal patterns of heteroplasmy transmission in humans (Sekiguchi et al. 2003; Goto et al. 2011; Sondheimer et al. 2011; Guo et al. 2013; Rebolledo-Jaramillo et al. 2014), including studies of oocytes and placenta (Marchington et al. 1997, 2002; Jacobs et al. 2007), and several questions remain.For example, although it is clear that a bottleneck occurs during the transmission of mtDNA genomes from mothers to offspring, the size of the bottleneck remains a contentious issue. Previous estimates of the effective number of transmitted mtDNA genomes range widely, from eight to 200 (Brown et al. 2001; Guo et al. 2013; Rebolledo-Jaramillo et al. 2014). However, all previous studies have assumed a constant size for the bottleneck across individuals; the effect of allowing the bottleneck size to vary among individuals has not been investigated. Moreover, it has been suggested that mtDNA genomes may not behave as independent entities but instead are organized into discrete units called “nucleoids,” each of which contains 5–10 mtDNA genomes (Jacobs et al. 2000; Cao et al. 2007; Khrapko 2008), although recently it has been suggested that the number may be much smaller, on the order of one mtDNA genome per nucleoid (Kukat et al. 2011). Each nucleoid is thought to be homoplasmic for mtDNA genome sequences; thus, mtDNA heteroplasmy at the cellular level would reflect nucleoids that are homoplasmic for different sequence variants. Nucleoid structures within cells have been studied microscopically and biochemically (Bogenhagen 2012), and nucleoid-based models have been found to provide a better fit to the segregation of heteroplasmic mtDNA genomes in cell lines than do simple bottleneck models in some studies (Cao et al. 2007; Khrapko 2008), but not in others (Cree et al. 2008). However, to date, nucleoid-based models have not been investigated in the transmission of mtDNA heteroplasmy from mothers to offspring.Another issue is the degree to which negative (or purifying) selection may act on deleterious variants during the transmission of mtDNA heteroplasmy. There are conflicting results and views as to whether changes in the frequency of a heteroplasmic mutation from mother to offspring are governed solely by genetic drift, or whether there is an additional role for negative (purifying) selection (Jenuth et al. 1997; Durham et al. 2006; Stewart et al. 2008a,b; Wonnapinij et al. 2008; Wallace and Chalkia 2013; Rebolledo-Jaramillo et al. 2014). Negative selection during heteroplasmy transmission, as evidenced by a decrease in the frequency of presumably deleterious heteroplasmic variants in offspring compared to mothers, must operate on the female germ line and/or early in development after fertilization, and hence differs from negative selection operating on homoplasmic variants that reduce viability or fertility (Holt et al. 2014). The opportunities for, and extent of, such negative selection during heteroplasmy transmission in humans remain largely unknown.Here, we utilize the Genomes of the Netherlands (GoNL) project (Boomsma et al. 2014; Genome of the Netherlands Consortium 2014), consisting of whole-genome sequence data from blood samples from 250 families, to carry out the largest study to date (to our knowledge) of the dynamics of heteroplasmy transmission across the entire mtDNA genome. We utilize the data on changes in minor allele frequency (MAF) from mothers to offspring at heteroplasmic sites to compare different models for the inheritance of mtDNA genomes, and we analyze the data for evidence of negative selection during heteroplasmy transmission.     

The initial peopling of the Americas: A growing number of founding mitochondrial genomes from Beringia

Pan-American mitochondrial DNA (mtDNA) haplogroup C1 has been recently subdivided into three branches, two of which (C1b and C1c) are characterized by ages and geographical distributions that are indicative of an early arrival from Beringia with Paleo-Indians. In contrast, the estimated ages of C1d—the third subset of C1—looked too young to fit the above scenario. To define the origin of this enigmatic C1 branch, we completely sequenced 63 C1d mitochondrial genomes from a wide range of geographically diverse, mixed, and indigenous American populations. The revised phylogeny not only brings the age of C1d within the range of that of its two sister clades, but reveals that there were two C1d founder genomes for Paleo-Indians. Thus, the recognized maternal founding lineages of Native Americans are at least 15, indicating that the overall number of Beringian or Asian founder mitochondrial genomes will probably increase extensively when all Native American haplogroups reach the same level of phylogenetic and genomic resolution as obtained here for C1d.While debate is still ongoing among scientists from several disciplines regarding the number of migratory events, their timing, and entry routes into the Americas (Wallace and Torroni 1992; Torroni et al. 1993; Forster et al. 1996; Kaufman and Golla 2000; Goebel et al. 2003, 2008; Schurr and Sherry 2004; Wang et al. 2007; Waters and Stafford 2007; Dillehay et al. 2008; Gilbert et al. 2008a; O''Rourke and Raff 2010), the general consensus is that modern Native American populations ultimately trace their gene pool to Asian groups who colonized northeast Siberia, including parts of Beringia, prior to the last glacial period. These ancestral population(s) probably retreated into refugial areas during the Last Glacial Maximum (LGM), where their genetic variation was reshaped by drift. Thus, pre-LGM haplotypes of Asian ancestry were differently preserved and lost in Beringian enclaves, but at the same time, novel haplotypes and alleles arose in situ due to new mutations, often becoming predominant because of major founder events (Tamm et al. 2007; Achilli et al. 2008; Bourgeois et al. 2009; Perego et al. 2009; Schroeder et al. 2009). The scenario of a temporally important differentiation stage in Beringia explains the predominance in Native Americans of private alleles and haplogroups such as the autosomal 9-repeat at microsatellite locus D9S1120 (Phillips et al. 2008; Schroeder et al. 2009), the Y chromosome haplogroup Q1a3a-M3 (Bortolini et al. 2003; Karafet et al. 2008; Rasmussen et al. 2010), and the pan-American mtDNA haplogroups A2, B2, C1b, C1c, C1d, D1, and D4h3a (Tamm et al. 2007; Achilli et al. 2008; Fagundes et al. 2008; Perego et al. 2009).In the millennia after the initial Paleo-Indian migrations, other groups from Beringia or eastern Siberia expanded into North America. If the gene pool of the source population(s) had in the meantime partially changed, not only because of drift, but also due to the admixture with population groups newly arrived from regions located west of Beringia, this would have resulted in the entry of additional Asian lineages into North America. This scenario, sometimes invoked to explain the presence of certain mtDNA haplogroups such as A2a, A2b, D2a, D3, and X2a only in populations of northern North America (Torroni et al. 1992; Brown et al. 1998; Schurr and Sherry 2004; Helgason et al. 2006; Achilli et al. 2008; Gilbert et al. 2008b; Perego et al. 2009), has recently received support from nuclear and morphometric data showing that some native groups from northern North America harbor stronger genetic similarities with some eastern Siberian groups than with Native American groups located more in the South (González-José et al. 2008; Bourgeois et al. 2009; Wang et al. 2009; Rasmussen et al. 2010).As for the pan-American mtDNA haplogroups, when analyzed at the highest level of molecular resolution (Bandelt et al. 2003; Tamm et al. 2007; Fagundes et al. 2008; Perego et al. 2009), they all reveal, with the exception of C1d, entry times of 15–18 thousand years ago (kya), which are suggestive of a (quasi) concomitant post-LGM arrival from Beringia with early Paleo-Indians. A similar entry time is also shown for haplogroup X2a, whose restricted geographical distribution in northern North America appears to be due not to a later arrival, but to its entry route through the ice-free corridor (Perego et al. 2009). Despite its continent-wide distribution, C1d was hitherto characterized by an expansion time of only 7.6–9.7 ky (Perego et al. 2009). This major discrepancy has been attributed to a poor and possibly biased representation of complete C1d mtDNA sequences (only 10) in the available data sets (Achilli et al. 2008; Malhi et al. 2010). To clarify the issue of the age of haplogroup C1d and its role as a founding Paleo-Indian lineage, we sequenced and analyzed 63 C1d mtDNAs from populations distributed over the entire geographical range of the haplogroup.     

APOBEC-induced mutations in human cancers are strongly enriched on the lagging DNA strand during replication

APOBEC3A and APOBEC3B, cytidine deaminases of the APOBEC family, are among the main factors causing mutations in human cancers. APOBEC deaminates cytosines in single-stranded DNA (ssDNA). A fraction of the APOBEC-induced mutations occur as clusters (“kataegis”) in single-stranded DNA produced during repair of double-stranded breaks (DSBs). However, the properties of the remaining 87% of nonclustered APOBEC-induced mutations, the source and the genomic distribution of the ssDNA where they occur, are largely unknown. By analyzing genomic and exomic cancer databases, we show that >33% of dispersed APOBEC-induced mutations occur on the lagging strand during DNA replication, thus unraveling the major source of ssDNA targeted by APOBEC in cancer. Although methylated cytosine is generally more mutation-prone than nonmethylated cytosine, we report that methylation reduces the rate of APOBEC-induced mutations by a factor of roughly two. Finally, we show that in cancers with extensive APOBEC-induced mutagenesis, there is almost no increase in mutation rates in late replicating regions (contrary to other cancers). Because late-replicating regions are depleted in exons, this results in a 1.3-fold higher fraction of mutations residing within exons in such cancers. This study provides novel insight into the APOBEC-induced mutagenesis and describes the peculiarity of the mutational processes in cancers with the signature of APOBEC-induced mutations.Carcinogenesis is associated with elevated mutation rates due to abnormal metabolic activities in the cell, disruption of repair systems, or environmental factors such as UV light, radiation, and chemical damage (Roberts and Gordenin 2014a,b). However, some normal protein enzymatic activities can also be a source of DNA damage and mutations. Recently, it was shown that some homologs of APOBEC (apolipoprotein B mRNA-editing enzyme, catalytic polypeptide-like), cytidine deaminases that function as viral protecting agents as well as in RNA editing, may be a major factor causing mutations in human cancers (Nik-Zainal et al. 2012; Burns et al. 2013b; Roberts et al. 2013). Deamination of cytidine residues by APOBEC occurs in single-stranded DNA (ssDNA) (Nowarski et al. 2008; Roberts et al. 2012; Smith et al. 2012). Two members of the APOBEC family, APOBEC3A and APOBEC3B, contribute substantially to mutations in cancers (Burns et al. 2013a,b; Roberts et al. 2013; Chan et al. 2015) by deaminating cytosines in the TpC context (henceforth, the mutated nucleotide is underlined) (Nik-Zainal et al. 2012; Burns et al. 2013a,b; Roberts et al. 2013; Taylor et al. 2013; Roberts and Gordenin 2014b; Chan et al. 2015). The APOBEC cytidine deaminase converts cytosines to uracils, which usually results in C → T or C → G mutations, and much less frequently, in C → A mutations (Taylor et al. 2013). The fact that the APOBEC shows the highest specificity for the TpCpW (where W denotes A or T) context was shown in cancer genomic studies and in experimental systems (Burns et al. 2013a,b; Roberts et al. 2013; Taylor et al. 2013).APOBEC-induced mutations are unevenly distributed along the genome. For example, under experimental conditions in yeasts, 26% of them are located in clusters spanning 6–15 kb (Taylor et al. 2013, 2014). This phenomenon, called kataegis, was described for many cancer types and is believed to be the result of APOBEC-induced mutagenesis (Nik-Zainal et al. 2012; Burns et al. 2013a,b; Roberts et al. 2013). Clustered mutations are frequently strand coordinated, i.e., are comprised of mutations in the TpC context that occur in one of the two strands (Nik-Zainal et al. 2012; Roberts et al. 2012, 2013). Although the majority of clusters carry mutations in one strand, 13% of the clusters exhibit strand switches, e.g., when the 5′ part of the cluster carries TpC coordinated mutations on the forward strand, and the 3′ part, on the reverse strand (corresponding to GpA mutations on the forward strand) (Roberts et al. 2012; Taylor et al. 2013). It was shown in cancers and in yeast experimental models that both coordinated and switching clusters are associated with DNA double-stranded breaks (DSBs) (Nik-Zainal et al. 2012; Roberts et al. 2013; Taylor et al. 2013) as a result of the activity of exonucleases causing long stretches of ssDNA near the DSB, which become a target for APOBEC enzymes (Roberts et al. 2012; Taylor et al. 2013). Alternatively, it was suggested that APOBEC enzymes can induce DSBs (Landry et al. 2011; Burns et al. 2013a). Another putative cause of kataegis is the expansion of ssDNA at the 5′ upstream region of a mismatch during base excision repair (BER) (Taylor et al. 2013; Chen et al. 2014). However, the majority of APOBEC mutations are dispersed (Nik-Zainal et al. 2012; Roberts et al. 2012, 2013; Taylor et al. 2013), and the source of ssDNA that may be a substrate for them in cancer still lacks explicit mechanistic explanation.During replication, DNA exists for some time in a single-stranded state. Although such ssDNA should be protected by the replication protein A (RPA), it may be a substrate for APOBEC-induced deamination, especially under replication stress (Roberts et al. 2012; Roberts and Gordenin 2014a). The lagging strand is single-stranded for a longer period of time than the leading strand due to discontinuous synthesis (Okazaki et al. 1968) and is also enriched in mutations (Reijns et al. 2015). Despite that firing of individual replication origins is stochastic (Rhind et al. 2010), genomic regions vary in mean time of the replication during the S phase (Ryba et al. 2010; Pope et al. 2014) and in their propensity to be replicated unidirectionally, and the preferential fork direction is conserved among human tissues (Baker et al. 2012).Here, we hypothesize that the APOBEC-induced mutagenesis is associated with the lagging strand. By the analysis of large genomic and exomic cancer data sets, we investigate the source of ssDNA targeted by APOBEC in cancer as well as the other APOBEC mutational properties.     

Evolutionary assembly patterns of prokaryotic genomes

Evolutionary innovation must occur in the context of some genomic background, which limits available evolutionary paths. For example, protein evolution by sequence substitution is constrained by epistasis between residues. In prokaryotes, evolutionary innovation frequently happens by macrogenomic events such as horizontal gene transfer (HGT). Previous work has suggested that HGT can be influenced by ancestral genomic content, yet the extent of such gene-level constraints has not yet been systematically characterized. Here, we evaluated the evolutionary impact of such constraints in prokaryotes, using probabilistic ancestral reconstructions from 634 extant prokaryotic genomes and a novel framework for detecting evolutionary constraints on HGT events. We identified 8228 directional dependencies between genes and demonstrated that many such dependencies reflect known functional relationships, including for example, evolutionary dependencies of the photosynthetic enzyme RuBisCO. Modeling all dependencies as a network, we adapted an approach from graph theory to establish chronological precedence in the acquisition of different genomic functions. Specifically, we demonstrated that specific functions tend to be gained sequentially, suggesting that evolution in prokaryotes is governed by functional assembly patterns. Finally, we showed that these dependencies are universal rather than clade-specific and are often sufficient for predicting whether or not a given ancestral genome will acquire specific genes. Combined, our results indicate that evolutionary innovation via HGT is profoundly constrained by epistasis and historical contingency, similar to the evolution of proteins and phenotypic characters, and suggest that the emergence of specific metabolic and pathological phenotypes in prokaryotes can be predictable from current genomes.A fundamental question in evolutionary biology is how present circumstances affect future adaptation and phenotypic change (Gould and Lewontin 1979). Studies of specific proteins, for example, indicate that epistasis between sequence residues limits accessible evolutionary trajectories and thereby renders certain adaptive paths more likely than others (Weinreich et al. 2006; Gong et al. 2013; de Visser and Krug 2014; Harms and Thornton 2014). Similarly, both phenotypic characters (Ord and Summers 2015) and specific genetic adaptations (Conte et al. 2012; Christin et al. 2015) show strong evidence of parallel evolution rather than convergent evolution. That is, a given adaptation is more likely to repeat in closely related organisms than in distantly related ones. This inverse relationship between the repeatability of evolution and taxonomic distance implies a strong effect of lineage-specific contingency on evolution, also potentially mediated by epistasis (Orr 2005).Such observations suggest that genetic adaptation is often highly constrained, and the present state of an evolving system can impact future evolution. Yet, the aforementioned studies are limited to small data sets and specific genetic pathways, and a more principled understanding of the rules by which future evolutionary trajectories are governed by the present state of the system is still lacking. For example, it is not known whether such adaptive constraints are a feature of genome-scale evolution or whether they are limited to finer scales. Moreover, the mechanisms that underlie observed constraints are often completely unknown. Addressing these questions is clearly valuable for obtaining a more complete theory of evolutionary biology, but more pressingly, is essential for tackling a variety of practical concerns, including our ability to combat evolving infectious diseases or engineer complex biological systems.Here, we address this challenge by analyzing horizontal gene transfer (HGT) in prokaryotes. HGT is an ideal system to systematically study genome-wide evolutionary constraints because it involves gene-level innovation, occurs at very high rates relative to sequence substitution (Nowell et al. 2014; Puigbò et al. 2014), and is a principal source of evolutionary novelty in prokaryotes (Gogarten et al. 2002; Jain et al. 2003; Lerat et al. 2005; Puigbò et al. 2014). Clearly, many or most acquired genes are rapidly lost due to fitness costs (van Passel et al. 2008; Baltrus 2013; Soucy et al. 2015), indicating that genes retained in the long term are likely to provide a selective advantage. Moreover, not all genes are equally transferrable (Jain et al. 1999; Sorek et al. 2007; Cohen et al. 2011), and not all species are equally receptive to the same genes (Smillie et al. 2011; Soucy et al. 2015). However, differences in HGT among species have been attributed not only to ecology (Smillie et al. 2011) or to phylogenetic constraints (Popa et al. 2011; Nowell et al. 2014), but also to interactions with the host genome (Jain et al. 1999; Cohen et al. 2011; Popa et al. 2011). Indeed, studies involving single genes or single species support the influence of genome content on the acquisition and retention of transferred genes (Pal et al. 2005; Sorek et al. 2007; Iwasaki and Takagi 2009; Chen et al. 2011; Press et al. 2013; Johnson and Grossman 2014). For example, it has been demonstrated that the presence of specific genes facilitates integration of others into genetic networks (Chen et al. 2011), and genes are more commonly gained in genomes already containing metabolic genes in the same pathway (Pal et al. 2005; Iwasaki and Takagi 2009). However, to date, a systematic, large-scale analysis of such dependencies has not been presented. In this paper, we therefore set out to characterize a comprehensive collection of genome-wide HGT-based dependencies among prokaryotic genes, analyze the obtained set of epistatic interactions, and identify patterns in the evolution of prokaryotic genomes.