Single-cell analysis reveals context-dependent, cell-level selection of mtDNA

Abstract

Heteroplasmy occurs when wild-type and mutant mitochondrial DNA (mtDNA) molecules co-exist in single cells¹. Heteroplasmy levels change dynamically in development, disease and ageing^2,3, but it is unclear whether these shifts are caused by selection or drift, and whether they occur at the level of cells or intracellularly. Here we investigate heteroplasmy dynamics in dividing cells by combining precise mtDNA base editing (DdCBE)⁴ with a new method, SCI-LITE (single-cell combinatorial indexing leveraged to interrogate targeted expression), which tracks single-cell heteroplasmy with ultra-high throughput. We engineered cells to have synonymous or nonsynonymous complex I mtDNA mutations and found that cell populations in standard culture conditions purge nonsynonymous mtDNA variants, whereas synonymous variants are maintained. This suggests that selection dominates over simple drift in shaping population heteroplasmy. We simultaneously tracked single-cell mtDNA heteroplasmy and ancestry, and found that, although the population heteroplasmy shifts, the heteroplasmy of individual cell lineages remains stable, arguing that selection acts at the level of cell fitness in dividing cells. Using these insights, we show that we can force cells to accumulate high levels of truncating complex I mtDNA heteroplasmy by placing them in environments where loss of biochemical complex I activity has been reported to benefit cell fitness. We conclude that in dividing cells, a given nonsynonymous mtDNA heteroplasmy can be harmful, neutral or even beneficial to cell fitness, but that the 'sign' of the effect is wholly dependent on the environment.

Fig. 1. SCI-LITE enables ultra-high-throughput analysis of targeted transcripts in single cells.

a, Fixed and permeabilized cells are distributed into wells in which targeted transcripts are labelled with well-specific barcodes. The first barcode is added during reverse transcription (RT). The second barcode and the UMI are added in the ligation step after the first round of pooling and splitting. The third and fourth barcodes are added by PCR after the second round of pooling and splitting. Lengths are presented in base pairs (bp). BC, barcode; GOI, gene of interest; Lig, ligation. Schematics of cells in part a were created using BioRender (https://biorender.com). b–e, Barnyard experiment. Reads were assigned to the HeLa or 293T cell line based on the unique sequence of the MT-ND4 transcript (b). Knee plot showing UMI per cell count indicating the number of barcodes corresponding to single cells (c). Barnyard plot showing the number of HeLa and 293T UMIs per cell (d). The two cell lines were mixed at equal ratios at the beginning of the experiment. Cells with alleles assigned to one cell line are considered singlets and coloured in blue (293T) or red (HeLa). Cells with mixed alleles are considered doublets and coloured in orange. Median UMIs detected per cell are presented as a function of raw sequencing reads (e). f, mtDNA and mtRNA depletion with EtBr measured by qPCR of mtDNA, RT–qPCR of mtRNA and SCI-LITE. n = 3 biological replicates. Error bars reflect the mean ± s.d. ****P ≤ 0.0001, ***P ≤ 0.001, *P ≤ 0.1 and not significant (NS) > 0.05 by Student’s unpaired two-tailed t-test. Source Data

Fig. 2. mtDNA base editing leads to a bimodal distribution of heteroplasmy.

a, Schematic overview of the SCI-LITE experiment. 293T cells were transfected with the ONC DdCBE, introducing a nonsense mutation in the MT-ND4 gene. Cells were cultured for 10 h, 24 h and 72 h, harvested and used for SCI-LITE. b, Western blot of ONC DdCBEs, showing expression of FLAG-tagged DdCBE halves (see Supplementary Fig. 2 for uncropped images). Actin was used as a loading control. n = 3 biological replicates are shown. c, Single-cell heteroplasmy in 293T cells interrogated using SCI-LITE. n = 3 biological replicates of edited cells are shown. Lines represent the mean for single biological replicates. Dots represent single cells. NT, not treated. d, Proposed model for heteroplasmy installation using mtDNA base editing. DdCBEs convert cytosine to uracil within double-stranded mtDNA molecules. Replication of one edited mtDNA molecule leads to the formation of one mutated and one wild-type mtDNA molecule. Further editing and replication of wild-type mtDNA molecules leads to formation of 50% mutant and 50% wild-type molecules. Replication of mutated molecules results in 100% mutant mtDNA molecules. e, For a single cell, assuming a fully active editor, after n rounds of mtDNA replication, we would expect that the heteroplasmy of the cell H = 1 − (1/2)ⁿ. Under this simple model, after just a few rounds of mtDNA replication in the face of an active base editor, a single cell will achieve a heteroplasmy level approaching 100%, but only if the mtDNA replication and turnover rate exceeds the cell division rate. Schematics in parts a,d,e were created using BioRender (https://biorender.com). Source Data

Fig. 3. Heteroplasmic shifting operates on nonsynonymous but not silent mtDNA variants.

293T cells were transfected with LHON or SILENT DdCBEs, and sorted based on the expression of the editors reflected by the fluorescence intensity of the eGFP and mCherry reporters. a, Alleles and their frequency introduced by LHON and SILENT DdCBEs. Missense indicates an on-target missense G11696A edit, and Silent indicates an on-target silent C11698T edit. b, Single-cell heteroplasmy in edited cells interrogated using SCI-LITE or by amplicon sequencing of bulk mtDNA and bulk mtRNA. Dots represent single cells. c,d, Heteroplasmy levels assessed by next-generation sequencing of MT-ND4 amplicon in high heteroplasmic cells treated with the LHON (c) or the SILENT (d) editor cultured in media containing glucose (GLU) or galactose (GAL). e,f, Population doublings and viability in LHON-edited and SILENT-edited low heteroplasmic (e) and high heteroplasmic (f) cells cultured in media containing either glucose or galactose. For c–f, n = 3 independent biological replicates. Error bars reflect the mean ± s.d. ****P ≤ 0.0001, ***P ≤ 0.001 and NS > 0.05, by Student’s unpaired two-tailed t-test. g, Joint distribution of missense and silent heteroplasmy changes over time. The graphs are sorted on the x axis by missense heteroplasmy and show the missense and silent heteroplasmy introduced by the LHON DdCBE in single cells measured by SCI-LITE. Each column represents one cell, with the stacked colours representing the percent missense, silent and wild-type heteroplasmy in red, yellow and blue, respectively. The dashed lines are for reference and indicate the midpoints of the x and y axes. h,i, Single-cell heteroplasmy in LHON-edited (h) and SILENT-edited (i) cells. Cells were grown for 5, 10 and 15 days and were subjected to SCI-LITE. The graphs on the left show the binned relative frequency and the graphs on the right show cumulative distributions of MT-ND4 missense and silent heteroplasmy for n = 3 independent biological replicates. The Kolmogorov–Smirnov test was used to calculate D statistics and P values. Source Data

Fig. 4. Selection against nonsynonymous mtDNA variants occurs at the level of cell fitness.

a, Schematic overview of the lineage tracing experiment. 293T cells were transfected with LHON or SILENT DdCBEs and subsequently transduced with a lentiviral library with unique ancestry barcodes so that each heteroplasmic cell expressed a single, unique ancestry barcode. Ancestral lineages were expanded, cells were harvested at day 0 and day 5, and multiplexed SCI-LITE was performed to capture mtDNA heteroplasmy and ancestry barcodes. Schematics in part a were creating using BioRender (https://biorender.com). b, Distribution of differences in heteroplasmy levels between day 5 and day 0, for all ancestral lineages that were detected at both days. c, Heteroplasmy levels in randomly selected ancestral lineages at day 0 and day 5. Each line represents one unique ancestral lineage and visualizes the mean heteroplasmy level at each day. See Extended Data Fig. 8 for all ancestral lineages. d, Single-cell heteroplasmy in LHON-edited cells. The graphs show empirical cumulative distributions of missense heteroplasmy in ancestry lineages found at only one timepoint or at both timepoints. The Kolmogorov–Smirnov test was used to calculate D statistics and P values. e, Relative frequency of binned missense heteroplasmy in ancestral lineages found at day 0 only or at both timepoints. Cells with ancestry barcodes that are found only at day 0 have significantly higher missense heteroplasmy by Kolmogorov–Smirnov test, suggesting cell-level selection against these lineages. Source Data

Fig. 5. A truncating complex I mtDNA mutation can be harmful, neutral or even beneficial to cell fitness depending on the environment.

a, Bulk mtDNA heteroplasmy in 293T cells transfected with ONC DdCBEs. Error bars represent the mean ± s.d. for n = 3 biological replicates. **P < 0.01 and NS > 0.05, by Student’s unpaired two-tailed t-test. b, Nthy-ori cells were edited with an active (ONC) or inactive (mock) DdCBE and 2 × 10⁶ cells were implanted into immunodeficient mice with an equal volume of Matrigel (50 μl). Tumour volume measurements of xenografts are shown. Error bars indicate mean ± s.e.m.; n = 5; *P < 0.05, by Mann–Whitney unpaired two-tailed test. c, Single-cell heteroplasmy in ONC-edited cells interrogated using SCI-LITE for n = 3 biological replicates and 34,666 single cells. See Extended Data Fig. 9b for a strip plot showing single replicates. d,e, CellTrace analysis. K562 cells edited with the ONC DdCBE were stained with a fluorescent dye allowing for tracing of cell divisions. Cells with higher CellTrace stain intensity underwent fewer divisions than cells with lower CellTrace stain intensity (d). Cells were stained with CellTrace and cultured for 4 days in media containing glucose, galactose or oligomycin, then harvested and sorted based on CellTrace intensity (e). Heteroplasmy levels were measured by bulk amplicon sequencing in sorted populations, and these measurements revealed decreased (galactose) or increased (oligomycin) heteroplasmy levels in the proliferating cells. Data for n = 3 biological replicates are shown. Source Data

Extended Data Fig. 1. Heterogeneity primers increase nucleotide diversity during sequencing of SCI-LITE libraries.

a, Schematic of heterogeneity primer design for amplifying SCI-LITE libraries. b, Representative zoomed-in cluster images from a single cycle of Illumina sequencing, for a run with PhiX spike-in (left) and for a run with a library prepared with heterogeneity spacers (right) generated with the Illumina Sequencing Analysis Viewer. The subplots are labeled with the called base and the percentage of clusters called for each base c, Base diversity at each bp of read 1 based on sequencing data from libraries prepared with or without heterogeneity primers.

Extended Data Fig. 2. Expression of DdCBEs correlates to heteroplasmy levels.

a, Western blots of ONC DdCBEs showing expression of FLAG-tagged right-TALE DdCBE and HA-tagged left-TALE DdCBE. Actin was used as a loading control. Data for n = 3 biological replicates are shown. b, Population doubling rate in 293 T cells transfected with DdCBEs, n = 3 biological replicates. Error bars reflect the mean ± s.d. c-d, 293 T cells were transfected with LHON DdCBE and sorted based on the fluorescence intensity of eGFP and mCherry reporters into low, medium and high groups. eGFP and mCherry intensity was used as a proxy for the expression of the editors c, Sorting strategy. Representative gating plots are shown. Graphs were generated using FCS Express software. d, Single-cell heteroplasmy measured by SCI-LITE, n = 3 biological replicates of edited cells are shown, except day 5 low where n = 2 biological replicates. Lines represent the mean for single biological replicates. Dots represent single cells. Source Data

Extended Data Fig. 3. Characterization of heteroplasmic cells.

a, Oxygen consumption rate (OCR) in LHON- and SILENT-edited cells for n = 3 biological replicates. Error bars reflect the mean ± s.e.m. and b, relative OCR values in LHON- and SILENT-edited cells for each of the segments in a that characterize different respiratory parameters. Data for n = 3 independent biological replicates are shown. Error bars reflect the mean ± s.d. *P < 0.05; **P < 0.01; by Student’s unpaired two-tailed t-test. FCCP, carbonyl cyanide-4-(trifluoromethoxy)phenylhydrazone. c, OCR in LHON-edited cells with low and high levels of missense heteroplasmy, n = 3 independent biological replicates. Error bars reflect the mean ± s.e.m. d, Single cell heteroplasmy measured by SCI-LITE in Nthy-ori cells edited with oncocytoma (ONC) editor and sorted for low and high levels of editor expression, n = 3 biological replicates of edited cells are shown. Lines represent the mean for single biological replicates. Dots represent single cells. e, OCR in edited Nthy-ori cells with low and high levels of nonsense heteroplasmy, n = 3 independent biological replicates. Error bars reflect the mean ± s.e.m. f, Population doublings in ONC-edited Nthy-ori cells with low and high levels of nonsense heteroplasmy, n = 3 independent biological replicates. Error bars reflect the mean ± s.d. **P < 0.01; by Student’s unpaired two-tailed t-test. g, Relative mtDNA levels in low and high heteroplasmic 293 T cells edited with LHON, SILENT and ONC DdCBE measured by qPCR, n = 3 independent biological replicates. Error bars reflect the mean ± s.d. ns > 0.05; by Student’s unpaired two-tailed t-test h, Bulk heteroplasmy measurement in cells treated with LHON or SILENT editors showing mean frequency of each allele, n = 3 biological replicates. Error bars reflect the mean ± s.d. Source Data

Extended Data Fig. 4. Cells with high levels of missense heteroplasmy have a growth disadvantage in standard cell culture conditions.

Population doublings per day (a) and cumulative cell counts (b) in LHON- and SILENT-edited high heteroplasmic cells. In b, numbers reflect the fraction of total LHON-edited cells vs SILENT-edited cells. Data for n = 3 independent biological replicates are shown. Error bars reflect the mean ± s.d. *P ≤ 0.05; ***P ≤ 0.001; ns > 0.05; by Student’s unpaired two-tailed t-test.

Extended Data Fig. 5. Kimura distribution fitting.

Density histograms showing the distribution of observed heteroplasmy in SILENT-edited (a) or LHON-edited (b) cells, and their corresponding maximum likelihood (ML) fits to the two-parameter Kimura(p,b) distribution after fixing the initial heteroplasmy parameter, p, to be equal to the experimentally-observed day 0 heteroplasmy. For the LHON-edited cells, the large shift in mean heteroplasmy between day 0 and day 5 cannot be accounted for by fitting only the Kimura drift parameter, b, and yields a poor fit (K-S D = 0.26) in contrast to the superior fit observed for the SILENT-edited cells (K-S D = 0.05). Source Data

Extended Data Fig. 6. Decline in nonsynonymous mtDNA heteroplasmy can be largely explained by modeling selection at the level of cell fitness.

We modeled heteroplasmy dynamics in LHON-edited cells under the assumption that the relative probability of a cell dividing, compared to any other cell, is entirely determined by its missense heteroplasmy level, with no intracellular effects. The model was fit based on the LHON-edited SCI-LITE time course, in which cells were grown for 5, 10 and 15 days. a-b, 2D kernel density estimate plots showing the joint missense and silent heteroplasmy distribution in observed data of the LHON-edited cells (a) and in simulated data based on our modeling (b). c, Distribution of mean heteroplasmy estimates from n = 100 simulated experiments (box plots) compared to the observed mean heteroplasmy based on SCI-LITE data (dotted horizontal lines). The boxes show the interquartile range, the line through the box shows the median, the whiskers show 1.5 * interquartile range. d, The chance of a cell dividing on each simulated day of the experiment was weighted by an inverse sigmoid function of the cell’s missense heteroplasmy. The fitted inverse sigmoid function model corresponds to a threshold at ~56% heteroplasmy. Source Data

Extended Data Fig. 7. Parameter search for computer model of selection acting at level of cell fitness.

a, Fitting the three parameters (location, depth, and pitch) of the inverse sigmoid function resulted in the function shown. b, 3D scatter plot showing the 20,000 models we trained to fit the inverse sigmoid to our data. Each dot represents a model trained with a random combination of the three parameters, and the color indicates the log-transformed smoothed and combined mean squared error (see Methods for details) of the resulting fit to our observed data. The minimum MSE is indicated by a large gold-colored marker. c, A range of multinomial sample sizes per cell for simulating the heteroplasmy was tested, and 1000 was chosen. The boxplots represent a total of n = 20,000 models. The 500 box represents 5071 models, the 1000 represents 4979 models, the 2000 represents 4944 models, and the 4000 represents 5006 models. The boxes show the interquartile range, the line through the box shows the median, the whiskers show 1.5 * interquartile range. d, Plots showing the marginal 2D facets of the 3D plot in panel b. Each hexagon is colored based on the mean smoothed and combined MSE for models in that 2D bin. Red lines indicate the chosen optimal parameter combination.

Extended Data Fig. 8. Observed heteroplasmy in ancestral lineages.

a-b, Heteroplasmy levels in all 720 ancestral lineage with cells detected at both days (a), and a random subsample of 100 of those lineages (b), in LHON-edited and SILENT-edited populations at day 0 and day 5. Each line segment represents a unique ancestral lineage and connects that lineage’s mean heteroplasmy at day 0 with that of day 5.

Extended Data Fig. 9. Environment dictates the sign of mtDNA selection.

a, Single cell heteroplasmy in ONC-edited cells interrogated using SCI-LITE. Graph shows cumulative distributions of MT-ND4 nonsense heteroplasmy for n = 3 biological replicates. Number of cells at each day: Day 3: n = 16,544; Day 6: n = 8,770; Day 10: n = 9,352. Kolmogorov–Smirnov test was used to calculate D statistics and P value. b, Single-cell heteroplasmy in ONC-edited cells interrogated using SCI-LITE. Graph shows single replicates from main Fig. 5c. Data for n = 3 biological replicates of edited cells are shown. Lines represent the mean for single biological replicates. Dots represent single cells. c, Western blots of ONC DdCBE showing expression of FLAG-tagged DdCBE half. Data for n = 3 biological replicates are shown.