Understanding the genetic and epigenetic regulation of mitochondrial DNA (mtDNA) is essential for elucidating mechanisms of aging and disease. Long-read sequencing can span the entire mitochondrial genome and directly capture base-modification signals, yet analytical tools for such data remain limited. We developed Himito, a graph-based toolkit for analyzing mitochondrial genome using long reads. Himito filters reads originating from nuclear mitochondrial insertions (NUMTs), constructs a sequence graph to represent mtDNA diversity, assembles primary haplotypes, calls variants, and analyzes 5-methylcytosine (5mC) modifications within a unified framework. Benchmarking on high-quality reference datasets shows Himito achieves superior performance in assembly and variant calling compared with existing tools. Applied to the All of Us (AoU) v8 dataset, Himito identified pathogenic mtDNA variants, revealed population-scale haplogroup diversity, and uncovered age-related genetic and epigenetic patterns. These results demonstrate that long-read sequencing, combined with graph-based analysis, enables integrated characterization of mitochondrial genomic and epigenomic variation. Himito is available at https://github.com/broadinstitute/Himito.

Introduction Sepsis is a highly morbid condition characterized by multi-organ dysfunction resulting from dysregulated inflammation in response to acute infection. Mitochondrial dysfunction may contribute to sepsis pathogenesis, but quantifying mitochondrial dysfunction remains challenging. Objective To assess the extent to which circulating markers of mitochondrial dysfunction are increased in septic shock, and their relationship to severity and mortality. Methods We performed both full-scan and targeted (known markers of genetic mitochondrial disease) metabolomics on plasma to determine markers of mitochondrial dysfunction which distinguish subjects with septic shock (n = 42) from cardiogenic shock without infection (n = 19), bacteremia without sepsis (n = 18), and ambulatory controls (n = 19) – the latter three being conditions in which mitochondrial function, proxied by peripheral oxygen consumption, is presumed intact. Results Nine metabolites were significantly increased in septic shock compared to all three comparator groups. This list includes N -formyl- l -methionine (f-Met), a marker of dysregulated mitochondrial protein translation, and N -lactoyl-phenylalanine (lac-Phe), representative of the N -lactoyl-amino acids (lac-AAs), which are elevated in plasma of patients with monogenic mitochondrial disease. Compared to lactate, the clinical biomarker used to define septic shock, there was greater separation between survivors and non-survivors of septic shock for both f-Met and the lac-AAs measured within 24 h of ICU admission. Additionally, tryptophan was the one metabolite significantly decreased in septic shock compared to all other groups, while its breakdown product kynurenate was one of the 9 significantly increased. Conclusion Future studies which validate the measurement of lac-AAs and f-Met in conjunction with lactate could define a sepsis subtype characterized by mitochondrial dysfunction.

Heteroplasmy occurs when wild-type and mutant mitochondrial DNA (mtDNA) molecules co-exist in single cells 1 . 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) 4 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.

Single cell ATAC-seq (scATAC-seq) enables the mapping of regulatory elements in fine-grained cell types. Despite this advance, analysis of the resulting data is challenging, and large scale scATAC-seq data are difficult to obtain and expensive to generate. This motivates a method to leverage information from previously generated large scale scATAC-seq or scRNA-seq data to guide our analysis of new scATAC-seq datasets. We analyze scATAC-seq data using latent Dirichlet allocation (LDA), a Bayesian algorithm that was developed to model text corpora, summarizing documents as mixtures of topics defined based on the words that distinguish the documents. When applied to scATAC-seq, LDA treats cells as documents and their accessible sites as words, identifying "topics" based on the cell type-specific accessible sites in those cells. Previous work used uniform symmetric priors in LDA, but we hypothesized that nonuniform matrix priors generated from LDA models trained on existing data sets may enable improved detection of cell types in new data sets, especially if they have relatively few cells. In this work, we test this hypothesis in scATAC-seq data from whole C. elegans nematodes and SHARE-seq data from mouse skin cells. We show that nonsymmetric matrix priors for LDA improve our ability to capture cell type information from small scATAC-seq datasets.

Mitochondrial DNA (mtDNA) is a maternally inherited, high-copy-number genome required for oxidative phosphorylation 1 . Heteroplasmy refers to the presence of a mixture of mtDNA alleles in an individual and has been associated with disease and ageing. Mechanisms underlying common variation in human heteroplasmy, and the influence of the nuclear genome on this variation, remain insufficiently explored. Here we quantify mtDNA copy number (mtCN) and heteroplasmy using blood-derived whole-genome sequences from 274,832 individuals and perform genome-wide association studies to identify associated nuclear loci. Following blood cell composition correction, we find that mtCN declines linearly with age and is associated with variants at 92 nuclear loci. We observe that nearly everyone harbours heteroplasmic mtDNA variants obeying two principles: (1) heteroplasmic single nucleotide variants tend to arise somatically and accumulate sharply after the age of 70 years, whereas (2) heteroplasmic indels are maternally inherited as mixtures with relative levels associated with 42 nuclear loci involved in mtDNA replication, maintenance and novel pathways. These loci may act by conferring a replicative advantage to certain mtDNA alleles. As an illustrative example, we identify a length variant carried by more than 50% of humans at position chrM:302 within a G-quadruplex previously proposed to mediate mtDNA transcription/replication switching 2 , 3 . We find that this variant exerts cis -acting genetic control over mtDNA abundance and is itself associated in- trans with nuclear loci encoding machinery for this regulatory switch. Our study suggests that common variation in the nuclear genome can shape variation in mtCN and heteroplasmy dynamics across the human population.

There is widespread interest in identifying interventions that extend healthy lifespan. Chronic continuous hypoxia delays the onset of replicative senescence in cultured cells and extends lifespan in yeast, nematodes, and fruit flies. Here, we asked whether chronic continuous hypoxia is beneficial in mammalian aging. We utilized the Ercc1 Δ/- mouse model of accelerated aging given that these mice are born developmentally normal but exhibit anatomic, physiological, and biochemical features of aging across multiple organs. Importantly, they exhibit a shortened lifespan that is extended by dietary restriction, the most potent aging intervention across many organisms. We report that chronic continuous 11% oxygen commenced at 4 weeks of age extends lifespan by 50% and delays the onset of neurological debility in Ercc1 Δ/- mice. Chronic continuous hypoxia did not impact food intake and did not significantly affect markers of DNA damage or senescence, suggesting that hypoxia did not simply alleviate the proximal effects of the Ercc1 mutation, but rather acted downstream via unknown mechanisms. To the best of our knowledge, this is the first study to demonstrate that "oxygen restriction" can extend lifespan in a mammalian model of aging.

Recently developed single-cell technologies allow researchers to characterize cell states at ever greater resolution and scale. Caenorhabditis elegans is a particularly tractable system for studying development, and recent single-cell RNA-seq studies characterized the gene expression patterns for nearly every cell type in the embryo and at the second larval stage (L2). Gene expression patterns give insight about gene function and into the biochemical state of different cell types; recent advances in other single-cell genomics technologies can now also characterize the regulatory context of the genome that gives rise to these gene expression levels at a single-cell resolution. To explore the regulatory DNA of individual cell types in C. elegans, we collected single-cell chromatin accessibility data using the sci-ATAC-seq assay in L2 larvae to match the available single-cell RNA-seq data set. By using a novel implementation of the latent Dirichlet allocation algorithm, we identify 37 clusters of cells that correspond to different cell types in the worm, providing new maps of putative cell type-specific gene regulatory sites, with promise for better understanding of cellular differentiation and gene regulation.

The human epigenome has been experimentally characterized by measurements of protein binding, chromatin acessibility, methylation, and histone modification in hundreds of cell types. The result is a huge compendium of data, consisting of thousands of measurements for every basepair in the human genome. These data are difficult to make sense of, not only for humans, but also for computational methods that aim to detect genes and other functional elements, predict gene expression, characterize polymorphisms, etc. To address this challenge, we propose a deep neural network tensor factorization method, Avocado, that compresses epigenomic data into a dense, information-rich representation of the human genome. We use data from the Roadmap Epigenomics Consortium to demonstrate that this learned representation of the genome is broadly useful: first, by imputing epigenomic data more accurately than previous methods, and second, by showing that machine learning models that exploit this representation outperform those trained directly on epigenomic data on a variety of genomics tasks. These tasks include predicting gene expression, promoter-enhancer interactions, replication timing, and an element of 3D chromatin architecture. Our findings suggest the broad utility of Avocado’s learned latent representation for computational genomics and epigenomics.

The mammalian mitochondrial proteome is under dual genomic control, with 99% of proteins encoded by the nuclear genome and 13 originating from the mitochondrial DNA (mtDNA). We previously developed MitoCarta, a catalogue of over 1000 genes encoding themammalianmitochondrial proteome. This catalogue was compiled using a Bayesian integration of multiple sequence features and experimental datasets, notably protein mass spectrometry of mitochondria isolated from fourteen murine tissues. Here, we introduce MitoCarta3.0. Beginning with the MitoCarta2.0 inventory, we performed manual review to remove 100 genes and introduce 78 additional genes, arriving at an updated inventory of 1136 human genes. We now include manually curated annotations of sub-mitochondrial localization (matrix, inner membrane, intermembrane space, outer membrane) as well as assignment to 149 hierarchical 'MitoPathways' spanning seven broad functional categories relevant to mitochondria. MitoCarta3.0, including sub-mitochondrial localization and MitoPathway annotations, is freely available at http://www. broadinstitute.org/ mitocarta and should serve as a continued community resource for mitochondrial biology and medicine.

The human epigenome has been experimentally characterized by thousands of measurements for every basepair in the human genome. We propose a deep neural network tensor factorization method, Avocado, that compresses this epigenomic data into a dense, information-rich representation. We use this learned representation to impute epigenomic data more accurately than previous methods, and we show that machine learning models that exploit this representation outperform those trained directly on epigenomic data on a variety of genomics tasks. These tasks include predicting gene expression, promoter-enhancer interactions, replication timing, and an element of 3D chromatin architecture.

In the past decade, the use of high-throughput sequencing assays has allowed researchers to experimentally acquire thousands of functional measurements for each basepair in the human genome. Despite their value, these measurements are only a small fraction of the potential experiments that could be performed while also being too numerous to easily visualize or compute on. In a recent pair of publications [1,2], we address both of these challenges with a deep neural network tensor factorization method, Avocado, that compresses these measurements into dense, information-rich representations. We demonstrate that these learned representations can be used to impute, with high accuracy, the output of tens of thousands of functional experiments that have not yet been performed. Further, we show that, on a variety of genomics tasks, machine learning models that leverage these learned representations outperform those trained directly on the functional measurements. The code is publicly available at https://github.com/jmschrei/avocado.

Additional file 3 Model performances by assay. The performance of ChromImpute, PREDICTD, and Avocado on the six performance measures shown in Table 1 calculated for each assay.

AbstractRecently developed single cell technologies allow researchers to characterize cell states at ever greater resolution and scale.C. elegansis a particularly tractable system for studying development, and recent single cell RNA-seq studies characterized the gene expression patterns for nearly every cell type in the embryo and at the second larval stage (L2). Gene expression patterns are useful for learning about gene function and give insight into the biochemical state of different cell types; however, in order to understand these cell types, we must also determine how these gene expression levels are regulated. We present the first single cell ATAC-seq study inC. elegans. We collected data in L2 larvae to match the available single cell RNA-seq data set, and we identify tissue-specific chromatin accessibility patterns that align well with existing data, including the L2 single cell RNA-seq results. Using a novel implementation of the latent Dirichlet allocation algorithm, we leverage the single-cell resolution of the sci-ATAC-seq data to identify accessible loci at the level of individual cell types, providing new maps of putative cell type-specific gene regulatory sites, with promise for better understanding of cellular differentiation and gene regulation in the worm.

Here, we present Scribe (https://github.com/aristoteleo/Scribe-py), a toolkit for detecting and visualizing causal regulatory interactions between genes and explore the potential for single-cell experiments to power network reconstruction. Scribe employs restricted directed information to determine causality by estimating the strength of information transferred from a potential regulator to its downstream target. We apply Scribe and other leading approaches for causal network reconstruction to several types of single-cell measurements and show that there is a dramatic drop in performance for "pseudotime"-ordered single-cell data compared with true time-series data. We demonstrate that performing causal inference requires temporal coupling between measurements. We show that methods such as "RNA velocity" restore some degree of coupling through an analysis of chromaffin cell fate commitment. These analyses highlight a shortcoming in experimental and computational methods for analyzing gene regulation at single-cell resolution and suggest ways of overcoming it.

Avocado’s model has seven structural hyperparameters: the number of latent factors representing cell types, assay types, and the three scales of genomic positions, as well as two parameters (number of layers and number of nodes per layer) for the deep neural network. We optimized these hyperparameters via random search. The search considered the following grid of values: cell type factors ∈ (16, 32, 64, 128, 256), assay factors ∈ (16, 32, 64, 128, 256), 25 bp resolution genome factors ∈ (5, 10, 15, 20, 25), 250 bp resolution genome factors ∈ (10, 20, 30, 40, 50), 5 kbp resolution genome factors ∈ (15, 30, 45, 60, 75), number of layers in the neural network ∈ (0, 1, 2, 3, 4), and number of neurons in the neural network ∈ (128, 256, 512, 1024, 2048). Note that setting the number of layers to 0 corresponds to training a linear regression model on top of the learned factors. These ranges were selected based on experimental results from Durham et al. [1], suggesting that 100 latent factors for each of the three axes performed well. In this grid we trained 1,000 models out of a possible ∼61,000. Each model was trained on the ENCODE Pilot Regions, which are comprised of 44 regions of 0.5-2 Mb length that jointly make up approximately 1% of the full genome. The data were split into a training set of 764 tracks, a validation set of 100 tracks, and a test set of 150 tracks. We selected the final set of hyperparameters based on performance on the validation set, as measured by mean-squared error (MSE). The different hyperparameter settings displayed a wide variance in performance, with most performing better than ChromImpute and many performing better than PREDICTD on the validation set (Additional file 1: Figure S1). Once the hyperparameters were set, the model was then retrained on both the training and validation sets and tested on the held-out test set. Note that the training, validation, and test sets used here correspond to the same splits used for the PREDICTD approach. The resulting model had a MSE of 0.1130 on the test set, which represents an 18.5% improvement over ChromImpute (MSE 0.1387) and a 4.9% improvement over PREDICTD (MSE 0.1188). We next investigated the effect that each hyperparameter had on the overall predictive performance of Avocado. To do this, we considered each hyperparameter individually and, for each value that the hyperparameter could take, we plotted the MSE of each model that used that value (Additional file 1: Figure S2). The clearest trend was that the performance of the model increased as the size of the neural network increased, both in terms of the number of layers and the number of neurons per layer. In contrast, the number of latent factors did not show a clear trend of improvement over any of the three axes. To attempt to better understand where the allocation of parameters was most beneficial, we considered performance when compared to the total number of parameters in the neural network and when compared to the total number of parameters in the embedding matrices (Additional file 1: Figure S3). We see that the validation set error decreases steadily with an increase in the number of network parameters until leveling off around 10 parameters. In particular, having no hidden layer, i.e., learning a linear regression on top of the tensor factorization, leads to very poor models. However, adding more than two layers does not yield much gain. When considering the number of parameters at each genomic position in the tensor factorization, we see no similar trend of increased complexity leading to increased performance. We focus on the number of parameters per genomic position rather than the total number of parameters in the model

The set of promoter-enhancer interactions used to evaluate the TargetFinder model [1] has been recently shown to contain biases related to the pairwise nature of the task [2]. The bias arises for two reasons. First, the data set includes features derived from the window between the promoter and the enhancer, and these features are highly correlated between examples whose windows overlap. This correlation leads to a leakage of information when regions of the genome are in windows of examples in both the training and the test set. Fortunately, this issue can be easily corrected by simply removing the problematic features. The second issue is that when the data set was constructed, an equal number of positive and negative interactions were sampled at each genomic distance. Consequently, many promoters occur repeatedly and only in the context of a negative interaction. When promoter-enhancer pairs are randomly assigned to both the training and test sets, as is the case with the TargetFinder model, then a sufficiently complicated model can simply memorize these repeated promoters as never interacting. These issues are described more thoroughly by Xi and Beer [2]. To construct a data set without these biases, we choose the simple approach of filtering out interactions such that each promoter occurs only once in each cell type. While we do not also enforce that enhancers can only occur once, we greedily select pairs where the enhancer has not yet been part of an example. This approach yields a data set with 27,048 interactions across all four cell types in chromosomes 1 through 22, where each interaction corresponds to a unique promoter in its cell type. Among these interactions, nearly all (26,707) have unique enhancers as well; 158 enhancers are seen twice, 7 are seen 3 times, and 1 is seen four times. After this filtering step, IMR90 has 4,702 pairs, of which 82 are positive interactions; GM12878 has 7,881 pairs, of which 181 are positive interactions; HeLa-S3 has 7,060 pairs, of which 121 are positive interactions; and K562 has 7,405 pairs, of which 145 are positive interactions. Promoters and enhancers

To develop a catalog of regulatory sites in two major model organisms, Drosophila melanogaster and Caenorhabditis elegans, the modERN (model organism Encyclopedia of Regulatory Networks) consortium has systematically assayed the binding sites of transcription factors (TFs). Combined with data produced by our predecessor, modENCODE (Model Organism ENCyclopedia Of DNA Elements), we now have data for 262 TFs identifying 1.23 M sites in the fly genome and 217 TFs identifying 0.67 M sites in the worm genome. Because sites from different TFs are often overlapping and tightly clustered, they fall into 91,011 and 59,150 regions in the fly and worm, respectively, and these binding sites span as little as 8.7 and 5.8 Mb in the two organisms. Clusters with large numbers of sites (so-called high occupancy target, or HOT regions) predominantly associate with broadly expressed genes, whereas clusters containing sites from just a few factors are associated with genes expressed in tissue-specific patterns. All of the strains expressing GFP-tagged TFs are available at the stock centers, and the chromatin immunoprecipitation sequencing data are available through the ENCODE Data Coordinating Center and also through a simple interface (http://epic.gs.washington.edu/modERN/) that facilitates rapid accessibility of processed data sets. These data will facilitate a vast number of scientific inquiries into the function of individual TFs in key developmental, metabolic, and defense and homeostatic regulatory pathways, as well as provide a broader perspective on how individual TFs work together in local networks and globally across the life spans of these two key model organisms.

BACKGROUND:Enhancers play an important role in morphological evolution and speciation by controlling the spatiotemporal expression of genes. Previous efforts to understand the evolution of enhancers in primates have typically studied many enhancers at low resolution, or single enhancers at high resolution. Although comparative genomic studies reveal large-scale turnover of enhancers, a specific understanding of the molecular steps by which mammalian or primate enhancers evolve remains elusive.RESULTS:We identified candidate hominoid-specific liver enhancers from H3K27ac ChIP-seq data. After locating orthologs in 11 primates spanning around 40 million years, we synthesized all orthologs as well as computational reconstructions of 9 ancestral sequences for 348 active tiles of 233 putative enhancers. We concurrently tested all sequences for regulatory activity with STARR-seq in HepG2 cells. We observe groups of enhancer tiles with coherent trajectories, most of which can be potentially explained by a single gain or loss-of-activity event per tile. We quantify the correlation between the number of mutations along a branch and the magnitude of change in functional activity. Finally, we identify 84 mutations that correlate with functional changes; these are enriched for cytosine deamination events within CpGs.CONCLUSIONS:We characterized the evolutionary-functional trajectories of hundreds of liver enhancers throughout the primate phylogeny. We observe subsets of regulatory sequences that appear to have gained or lost activity. We use these data to quantify the relationship between sequence and functional divergence, and to identify CpG deamination as a potentially important force in driving changes in enhancer activity during primate evolution.

The Encyclopedia of DNA Elements (ENCODE) and the Roadmap Epigenomics Project seek to characterize the epigenome in diverse cell types using assays that identify, for example, genomic regions with modified histones or accessible chromatin. These efforts have produced thousands of datasets but cannot possibly measure each epigenomic factor in all cell types. To address this, we present a method, PaRallel Epigenomics Data Imputation with Cloud-based Tensor Decomposition (PREDICTD), to computationally impute missing experiments. PREDICTD leverages an elegant model called “tensor decomposition” to impute many experiments simultaneously. Compared with the current state-of-the-art method, ChromImpute, PREDICTD produces lower overall mean squared error, and combining the two methods yields further improvement. We show that PREDICTD data captures enhancer activity at noncoding human accelerated regions. PREDICTD provides reference imputed data and open-source software for investigating new cell types, and demonstrates the utility of tensor decomposition and cloud computing, both promising technologies for bioinformatics.

Single-cell transcriptome sequencing now routinely samples thousands of cells, potentially providing enough data to reconstruct causal gene regulatory networks from observational data. Here, we present Scribe, a toolkit for detecting and visualizing causal regulatory interactions between genes and explore the potential for single-cell experiments to power network reconstruction. Scribe employs Restricted Directed Information to determine causality by estimating the strength of information transferred from a potential regulator to its downstream target. We apply Scribe and other leading approaches for causal network reconstruction to several types of single-cell measurements and show that there is a dramatic drop in performance for "pseudotime” ordered single-cell data compared to true time series data. We demonstrate that performing causal inference requires temporal coupling between measurements. We show that methods such as “RNA velocity” restore some degree of coupling through an analysis of chromaffin cell fate commitment. These analyses therefore highlight an important shortcoming in experimental and computational methods for analyzing gene regulation at single-cell resolution and point the way towards overcoming it.