Systems and methods for analyzing sequencing data

Notes

US Patent App. 17/476,113

Recent Publications

Gagoski D, Rube T, Rastogi C, Melo LA, Li X, Voleti R, Shah NH, Bussemaker HJ. Accurate sequence-to-affinity models for SH2 domains from multi-round peptide binding assays coupled with free-energy regression. bioRxiv. 2025;34(11). doi: https://doi.org/10.1002/pro.70317

Publisher's Version

Publisher's Version
Street LA, Majumdar M, Rothamel KL, Gogus E, Brar G, Yeo GW, Rube T, Jovanovic M. SPRINT-MS: A high-throughput platform for identifying protein-protein interactions using pooled IP-MS and sparse signal recovery. bioRxiv. Submitted Submitted.
Wolin E, Guo JK, Blanco MR, Goronzy IN, Gorhe D, Dong W, Perez AA, Keskin A, Valenzuela E, Abdou AA, et al. Resource SPIDR enables multiplexed mapping of RNA-protein interactions and uncovers a mechanism for selective translational suppression upon cell stress. Cell. 2025;188.

Publisher's Version

Publisher's Version
Zhang L, Li G, Zhang Y, Cheng Y, Roberts N, Glenn S, DeZwaan-McCabe D, Rube H, Manthey J, Coleman G, et al. Boosting genome editing efficiency in human cells and plants with novel LbCas12a variants.. Europe PMC. 2023.

Publisher's Version

Publisher's Version
Feng S, Rastogi C, Loker RE, Glassford WJ, Rube T, Bussemaker HJ, Mann RS. Transcription factor paralogs orchestrate alternative gene regulatory networks by context-dependent cooperation with multiple cofactors. bioRxiv. 2022.

Publisher's Version

Publisher's Version

In eukaryotes, members of large transcription factor families often exhibit similar DNA binding properties in vitro, yet initiate paralog-specific gene regulatory networks in vivo. The serially homologous first (T1) and third (T3) thoracic legs of Drosophila, which result from alternative gene regulatory networks specified by the Hox proteins Scr and Ubx, respectively, offer a unique opportunity to address this paradox in vivo. Genome-wide analyses using epitope-tagged alleles of both Hox loci in the T1 and T3 leg imaginal discs, which are the precursors to the adult appendages and ventral body regions, show that ∼8% of Hox binding is paralog-specific. Binding specificity is mediated by interactions with distinct cofactors in different domains: the known Hox cofactor Exd acts in the proximal domain and is necessary for Scr to bind many of its paralog-specific targets, while in the distal leg domain, we identified the homeodomain protein Distal-less (Dll) as a novel Hox cofactor that enhances Scr binding to a different subset of genomic loci. Reporter genes confirm the in vivo roles of Scr+Dll and suggest that ∼1/3 of paralog-specific Hox binding in enhancers is functional. Together, these findings provide a genome-wide view of how Hox paralogs, and perhaps paralogs of other transcription factor families, orchestrate alternative downstream gene networks and suggest the importance of multiple, context-specific cofactors.
Rube T, Rastogi C, Feng S, Kribelbauer JF, Li A, Becerra B, Melo LA, Do BV, Li X, Adam HH, et al. Prediction of protein-ligand binding affinity from sequencing data with interpretable machine learning. Nat. Biotechnol. 2022.

Publisher's Version

Publisher's Version

Protein-ligand interactions are increasingly profiled at high throughput using affinity selection and massively parallel sequencing. However, these assays do not provide the biophysical parameters that most rigorously quantify molecular interactions. Here we describe a flexible machine learning method, called ProBound, that accurately defines sequence recognition in terms of equilibrium binding constants or kinetic rates. This is achieved using a multi-layered maximum-likelihood framework that models both the molecular interactions and the data generation process. We show that ProBound quantifies transcription factor (TF) behavior with models that predict binding affinity over a range exceeding that of previous resources; captures the impact of DNA modifications and conformational flexibility of multi-TF complexes; and infers specificity directly from in vivo data such as ChIP-seq without peak calling. When coupled with an assay called KD-seq, it determines the absolute affinity of protein-ligand interactions. We also apply ProBound to profile the kinetics of kinase-substrate interactions. ProBound opens new avenues for decoding biological networks and rationally engineering protein-ligand interactions.
Basu A, Bobrovnikov DG, Qureshi Z, Kayikcioglu T, Ngo TTM, Ranjan A, Eustermann S, Cieza B, Morgan MT, Hejna M, et al. Measuring DNA mechanics on the genome scale. Nature. 2021;589(7842):462–467.

Publisher's Version

Publisher's Version

Mechanical deformations of DNA such as bending are ubiquitous and have been implicated in diverse cellular functions1. However, the lack of high-throughput tools to measure the mechanical properties of DNA has limited our understanding of how DNA mechanics influence chromatin transactions across the genome. Here we develop 'loop-seq'-a high-throughput assay to measure the propensity for DNA looping-and determine the intrinsic cyclizabilities of 270,806 50-base-pair DNA fragments that span Saccharomyces cerevisiae chromosome V, other genomic regions, and random sequences. We found sequence-encoded regions of unusually low bendability within nucleosome-depleted regions upstream of transcription start sites (TSSs). Low bendability of linker DNA inhibits nucleosome sliding into the linker by the chromatin remodeller INO80, which explains how INO80 can define nucleosome-depleted regions in the absence of other factors2. Chromosome-wide, nucleosomes were characterized by high DNA bendability near dyads and low bendability near linkers. This contrast increases for deeper gene-body nucleosomes but disappears after random substitution of synonymous codons, which suggests that the evolution of codon choice has been influenced by DNA mechanics around gene-body nucleosomes. Furthermore, we show that local DNA mechanics affect transcription through TSS-proximal nucleosomes. Overall, this genome-scale map of DNA mechanics indicates a 'mechanical code' with broad functional implications.
Zhang L, Zuris JA, Viswanathan R, Edelstein JN, Turk R, Thommandru B, Rube T, Glenn SE, Collingwood MA, Bode NM, et al. AsCas12a ultra nuclease facilitates the rapid generation of therapeutic cell medicines. Nat. Commun. 2021;12(1):4500.

Publisher's Version

Publisher's Version

Though AsCas12a fills a crucial gap in the current genome editing toolbox, it exhibits relatively poor editing efficiency, restricting its overall utility. Here we isolate an engineered variant, “AsCas12a Ultra”, that increased editing efficiency to nearly 100% at all sites examined in HSPCs, iPSCs, T cells, and NK cells. We show that AsCas12a Ultra maintains high on-target specificity thereby mitigating the risk for off-target editing and making it ideal for complex therapeutic genome editing applications. We achieved simultaneous targeting of three clinically relevant genes in T cells at >90% efficiency and demonstrated transgene knock-in efficiencies of up to 60%. We demonstrate site-specific knock-in of a CAR in NK cells, which afforded enhanced anti-tumor NK cell recognition, potentially enabling the next generation of allogeneic cell-based therapies in oncology. AsCas12a Ultra is an advanced CRISPR nuclease with significant advantages in basic research and in the production of gene edited cell medicines.
Bussemaker HH, Rastogi C, Rube HT. Systems and methods for analyzing sequencing data. 2021.
Zhang L, Rube T, Vakulskas CA, Behlke MA, Bussemaker HJ, Pufall MA. Systematic in vitro profiling of off-target affinity, cleavage and efficiency for CRISPR enzymes. Nucleic Acids Res. 2020;48(9):5037–5053.

Publisher's Version

Publisher's Version

CRISPR RNA-guided endonucleases (RGEs) cut or direct activities to specific genomic loci, yet each has off-target activities that are often unpredictable. We developed a pair of simple in vitro assays to systematically measure the DNA-binding specificity (Spec-seq), catalytic activity specificity (SEAM-seq) and cleavage efficiency of RGEs. By separately quantifying binding and cleavage specificity, Spec/SEAM-seq provides detailed mechanistic insight into off-target activity. Feature-based models generated from Spec/SEAM-seq data for SpCas9 were consistent with previous reports of its in vitro and in vivo specificity, validating the approach. Spec/SEAM-seq is also useful for profiling less-well characterized RGEs. Application to an engineered SpCas9, HiFi-SpCas9, indicated that its enhanced target discrimination can be attributed to cleavage rather than binding specificity. The ortholog ScCas9, on the other hand, derives specificity from binding to an extended PAM. The decreased off-target activity of AsCas12a (Cpf1) appears to be primarily driven by DNA-binding specificity. Finally, we performed the first characterization of CasX specificity, revealing an all-or-nothing mechanism where mismatches can be bound, but not cleaved. Together, these applications establish Spec/SEAM-seq as an accessible method to rapidly and reliably evaluate the specificity of RGEs, Cas::gRNA pairs, and gain insight into the mechanism and thermodynamics of target discrimination.
Kribelbauer JF, Loker RE, Feng S, Rastogi C, Abe N, Rube T, Bussemaker HJ, Mann RS. Context-dependent gene regulation by homeodomain transcription factor complexes revealed by shape-readout deficient proteins. Mol. Cell. 2020;78(1):152–167.e11.

Publisher's Version

Publisher's Version

Eukaryotic transcription factors (TFs) form complexes with various partner proteins to recognize their genomic target sites. Yet, how the DNA sequence determines which TF complex forms at any given site is poorly understood. Here, we demonstrate that high-throughput in vitro DNA binding assays coupled with unbiased computational analysis provide unprecedented insight into how different DNA sequences select distinct compositions and configurations of homeodomain TF complexes. Using inferred knowledge about minor groove width readout, we design targeted protein mutations that destabilize homeodomain binding both in vitro and in vivo in a complex-specific manner. By performing parallel systematic evolution of ligands by exponential enrichment sequencing (SELEX-seq), chromatin immunoprecipitation sequencing (ChIP-seq), RNA sequencing (RNA-seq), and Hi-C assays, we not only classify the majority of in vivo binding events in terms of complex composition but also infer complex-specific functions by perturbing the gene regulatory network controlled by a single complex.
Rastogi C, Rube T, Kribelbauer JF, Crocker J, Loker RE, Martini GD, Laptenko O, Freed-Pastor WA, Prives C, Stern DL, et al. Accurate and sensitive quantification of protein-DNA binding affinity. Proc. Natl. Acad. Sci. U. S. A. 2018;115(16):E3692–E3701.

Publisher's Version

Publisher's Version

Transcription factors (TFs) control gene expression by binding to genomic DNA in a sequence-specific manner. Mutations in TF binding sites are increasingly found to be associated with human disease, yet we currently lack robust methods to predict these sites. Here, we developed a versatile maximum likelihood framework named No Read Left Behind (NRLB) that infers a biophysical model of protein-DNA recognition across the full affinity range from a library of in vitro selected DNA binding sites. NRLB predicts human Max homodimer binding in near-perfect agreement with existing low-throughput measurements. It can capture the specificity of the p53 tetramer and distinguish multiple binding modes within a single sample. Additionally, we confirm that newly identified low-affinity enhancer binding sites are functional in vivo, and that their contribution to gene expression matches their predicted affinity. Our results establish a powerful paradigm for identifying protein binding sites and interpreting gene regulatory sequences in eukaryotic genomes.
Rube T, Rastogi C, Kribelbauer JF, Bussemaker HJ. A unified approach for quantifying and interpreting DNA shape readout by transcription factors. Mol. Syst. Biol. 2018;14(2):e7902.

Publisher's Version

Publisher's Version

Transcription factors (TFs) interpret DNA sequence by probing the chemical and structural properties of the nucleotide polymer. DNA shape is thought to enable a parsimonious representation of dependencies between nucleotide positions. Here, we propose a unified mathematical representation of the DNA sequence dependence of shape and TF binding, respectively, which simplifies and enhances analysis of shape readout. First, we demonstrate that linear models based on mononucleotide features alone account for 60-70% of the variance in minor groove width, roll, helix twist, and propeller twist. This explains why simple scoring matrices that ignore all dependencies between nucleotide positions can partially account for DNA shape readout by a TF Adding dinucleotide features as sequence-to-shape predictors to our model, we can almost perfectly explain the shape parameters. Building on this observation, we developed a post hoc analysis method that can be used to analyze any mechanism-agnostic protein-DNA binding model in terms of shape readout. Our insights provide an alternative strategy for using DNA shape information to enhance our understanding of how cis-regulatory codes are interpreted by the cellular machinery.
Zhang L, Martini GD, Rube T, Kribelbauer JF, Rastogi C, FitzPatrick VD, Houtman JC, Bussemaker HJ, Pufall MA. SelexGLM differentiates androgen and glucocorticoid receptor DNA-binding preference over an extended binding site. Genome Res. 2018;28(1):111–121.

Publisher's Version

Publisher's Version

The DNA-binding interfaces of the androgen (AR) and glucocorticoid (GR) receptors are virtually identical, yet these transcription factors share only about a third of their genomic binding sites and regulate similarly distinct sets of target genes. To address this paradox, we determined the intrinsic specificities of the AR and GR DNA-binding domains using a refined version of SELEX-seq. We developed an algorithm, SelexGLM, that quantifies binding specificity over a large (31-bp) binding site by iteratively fitting a feature-based generalized linear model to SELEX probe counts. This analysis revealed that the DNA-binding preferences of AR and GR homodimers differ significantly, both within and outside the 15-bp core binding site. The relative preference between the two factors can be tuned over a wide range by changing the DNA sequence, with AR more sensitive to sequence changes than GR. The specificity of AR extends to the regions flanking the core 15-bp site, where isothermal calorimetry measurements reveal that affinity is augmented by enthalpy-driven readout of poly(A) sequences associated with narrowed minor groove width. We conclude that the increased specificity of AR is correlated with more enthalpy-driven binding than GR. The binding models help explain differences in AR and GR genomic binding and provide a biophysical rationale for how promiscuous binding by GR allows functional substitution for AR in some castration-resistant prostate cancers.
Kim M, Kreig A, Lee C-Y, Rube T, Calvert J, Song JS, Myong S. Quantitative analysis and prediction of G-quadruplex forming sequences in double-stranded DNA. Nucleic Acids Res. 2016;44(10):4807–4817.

Publisher's Version

Publisher's Version

G-quadruplex (GQ) is a four-stranded DNA structure that can be formed in guanine-rich sequences. GQ structures have been proposed to regulate diverse biological processes including transcription, replication, translation and telomere maintenance. Recent studies have demonstrated the existence of GQ DNA in live mammalian cells and a significant number of potential GQ forming sequences in the human genome. We present a systematic and quantitative analysis of GQ folding propensity on a large set of 438 GQ forming sequences in double-stranded DNA by integrating fluorescence measurement, single-molecule imaging and computational modeling. We find that short minimum loop length and the thymine base are two main factors that lead to high GQ folding propensity. Linear and Gaussian process regression models further validate that the GQ folding potential can be predicted with high accuracy based on the loop length distribution and the nucleotide content of the loop sequences. Our study provides important new parameters that can inform the evaluation and classification of putative GQ sequences in the human genome.
Bell RJA, Rube T, Xavier-Magalhães A, Costa BM, Mancini A, Song JS, Costello JF. Understanding TERT promoter mutations: A common path to immortality. Mol. Cancer Res. 2016;14(4):315–323.

Publisher's Version

Publisher's Version

Telomerase (TERT) activation is a fundamental step in tumorigenesis. By maintaining telomere length, telomerase relieves a main barrier on cellular lifespan, enabling limitless proliferation driven by oncogenes. The recently discovered, highly recurrent mutations in the promoter of TERT are found in over 50 cancer types, and are the most common mutation in many cancers. Transcriptional activation of TERT, via promoter mutation or other mechanisms, is the rate-limiting step in production of active telomerase. Although TERT is expressed in stem cells, it is naturally silenced upon differentiation. Thus, the presence of TERT promoter mutations may shed light on whether a particular tumor arose from a stem cell or more differentiated cell type. It is becoming clear that TERT mutations occur early during cellular transformation, and activate the TERT promoter by recruiting transcription factors that do not normally regulate TERT gene expression. This review highlights the fundamental and widespread role of TERT promoter mutations in tumorigenesis, including recent progress on their mechanism of transcriptional activation. These somatic promoter mutations, along with germline variation in the TERT locus also appear to have significant value as biomarkers of patient outcome. Understanding the precise molecular mechanism of TERT activation by promoter mutation and germline variation may inspire novel cancer cell-specific targeted therapies for a large number of cancer patients.
Rube T, Lee W, Hejna M, Chen H, Yasui DH, Hess JF, LaSalle JM, Song JS, Gong Q. Sequence features accurately predict genome-wide MeCP2 binding in vivo. Nat. Commun. 2016;7(1):11025.

Publisher's Version

Publisher's Version

Methyl-CpG binding protein 2 (MeCP2) is critical for proper brain development and expressed at near-histone levels in neurons, but the mechanism of its genomic localization remains poorly understood. Using high-resolution MeCP2-binding data, we show that DNA sequence features alone can predict binding with 88% accuracy. Integrating MeCP2 binding and DNA methylation in a probabilistic graphical model, we demonstrate that previously reported genome-wide association with methylation is in part due to MeCP2's affinity to GC-rich chromatin, a result replicated using published data. Furthermore, MeCP2 co-localizes with nucleosomes. Finally, MeCP2 binding downstream of promoters correlates with increased expression in Mecp2-deficient neurons.
Jin H, Rube T, Song JS. Categorical spectral analysis of periodicity in nucleosomal DNA. Nucleic Acids Res. 2016;44(5):2047–2057.

Publisher's Version

Publisher's Version

DNA helical twist imposes geometric constraints on the location of histone-DNA interaction sites along nucleosomal DNA. Certain 10.5-bp periodic nucleotides in phase with these geometric constraints have been suggested to facilitate nucleosome positioning. However, the extent of nucleotide periodicity in nucleosomal DNA and its significance in directing nucleosome positioning still remain unclear. We clarify these issues by applying categorical spectral analysis to high-resolution nucleosome maps in two yeast species. We find that only a small fraction of nucleosomal sequences contain significant 10.5-bp periodicity. We further develop a spectral decomposition method to show that the previously observed periodicity in aligned nucleosomal sequences mainly results from proper phasing among nucleosomal sequences, and not from a preponderant occurrence of periodicity within individual sequences. Importantly, we show that this phasing may arise from the histones' proclivity for putting preferred nucleotides at some of the evenly spaced histone-DNA contact points with respect to the dyad axis. We demonstrate that 10.5-bp periodicity, when present, significantly facilitates rotational, but not translational, nucleosome positioning. Finally, although periodicity only moderately affects nucleosome occupancy genome wide, reduced periodicity is an evolutionarily conserved signature of nucleosome-depleted regions around transcription start/termination sites.
Bell RJA, Rube T, Kreig A, Mancini A, Fouse SD, Nagarajan RP, Choi S, Hong C, He D, Pekmezci M, et al. Cancer. The transcription factor GABP selectively binds and activates the mutant TERT promoter in cancer. Science. 2015;348(6238):1036–1039.

Publisher's Version

Publisher's Version

Reactivation of telomerase reverse transcriptase (TERT) expression enables cells to overcome replicative senescence and escape apoptosis, which are fundamental steps in the initiation of human cancer. Multiple cancer types, including up to 83% of glioblastomas (GBMs), harbor highly recurrent TERT promoter mutations of unknown function but specific to two nucleotide positions. We identified the functional consequence of these mutations in GBMs to be recruitment of the multimeric GA-binding protein (GABP) transcription factor specifically to the mutant promoter. Allelic recruitment of GABP is consistently observed across four cancer types, highlighting a shared mechanism underlying TERT reactivation. Tandem flanking native E26 transformation-specific motifs critically cooperate with these mutations to activate TERT, probably by facilitating GABP heterotetramer binding. GABP thus directly links TERT promoter mutations to aberrant expression in multiple cancers.
Rube T, Song JS. Quantifying the role of steric constraints in nucleosome positioning. Nucleic Acids Res. 2014;42(4):2147–2158.

Publisher's Version

Publisher's Version

Statistical positioning, the localization of nucleosomes packed against a fixed barrier, is conjectured to explain the array of well-positioned nucleosomes at the 5' end of genes, but the extent and precise implications of statistical positioning in vivo are unclear. We examine this hypothesis quantitatively and generalize the idea to include moving barriers as well as nucleosomes actively packed against a barrier. Early experiments noted a similarity between the nucleosome profile aligned and averaged across genes and that predicted by statistical positioning; however, we demonstrate that aligning random nucleosomes also generates the same profile, calling the previous interpretation into question. New rigorous results reformulate statistical positioning as predictions on the variance structure of nucleosome locations in individual genes. In particular, a quantity termed the variance gradient, describing the change in variance between adjacent nucleosomes, is tested against recent high-throughput nucleosome sequencing data. Constant variance gradients provide support for generalized statistical positioning in ∼ 50% of long genes. Genes that deviate from predictions have high nucleosome turnover and cell-to-cell gene expression variability. The observed variance gradient suggests an effective nucleosome size of 158 bp, instead of the commonly perceived 147 bp. Our analyses thus clarify the role of statistical positioning in vivo.
Kallosh R, Linde A, Olive KA, Rube T. Chaotic inflation and supersymmetry breaking. Physical Review D. 2011;84(8):083519.

Publisher's Version

Publisher's Version

We investigate the recently proposed class of chaotic inflation models in supergravity with an arbitrary inflaton potential V(\phi). These models are extended to include matter fields in the visible sector and we employ a mechanism of supersymmetry breaking based on a particular phenomenological version of the KKLT mechanism (the KL model). We describe specific features of reheating in this class of models and show how one can solve the cosmological moduli and gravitino problems in this context.
Rube T, Wacker JG. The simplicity of perfect atoms: degeneracies in supersymmetric hydrogen. Journal of mathematical physics. 2011;52(6):062102.

Publisher's Version

Publisher's Version

Supersymmetric QED hydrogen-like bound states are remarkably similar to nonsupersymmetric hydrogen, including an accidental degeneracy of the fine structure and is broken by the Lamb shift. This article classifies the states, calculates the leading order spectrum, and illustrates the results in several limits. The relation to other nonrelativistic bound states is explored.
Kallosh R, Linde A, Rube T. General inflaton potentials in supergravity. Physical Review D. 2011.

Publisher's Version

Publisher's Version

We describe a way to construct supergravity models with an arbitrary inflaton potential V(\phi) and show that all other scalar fields in this class of models can be stabilized at the inflationary trajectory by a proper choice of the Kähler potential.
Behbahani SR, Jankowiak M, Rube T, Wacker JG. Nearly supersymmetric dark atoms. Advances in High Energy Physics. 2011;2011.

Publisher's Version

Publisher's Version

Theories of dark matter that support bound states are an intriguing possibility for the identity of the missing mass of the Universe. This article proposes a class of models of supersymmetric composite dark matter where the interactions with the Standard Model communicate supersymmetry breaking to the dark sector. In these models, supersymmetry breaking can be treated as a perturbation on the spectrum of bound states. Using a general formalism, the spectrum with leading supersymmetry effects is computed without specifying the details of the binding dynamics. The interactions of the composite states with the Standard Model are computed, and several benchmark models are described. General features of nonrelativistic supersymmetric bound states are emphasized.
Kallosh R, Lee CH, Rube T. N= 8 Supergravity 4-point amplitudes. Journal of High Energy Physics. 2009;2009(02):050.

Publisher's Version

Publisher's Version

We present the explicit expressions in Script N= 8 supergravity for the bosonic 4-particle tree and 1-loop amplitudes including vectors and scalars. We also present the candidate 4-point UV divergences in a form of helicity amplitudes, corresponding to 3-loop manifestly Script N= 8 supersymmetric and Lorentz covariant counterterm. This may shed some light on the 3-loop finiteness of Script N= 8 SG and on a conjectured higher loop finiteness. We perform a supersymmetric deformation to complex momentum of the 4-point generating function including higher-loop counterterms and the 1-loop UV finite amplitudes. Using the explicit form of the scalar part of the 3-loop counterterm and of the 1-loop UV finite scalar 4-point amplitudes we find that they both have an unbroken E (7, 7) symmetry. We derive from E (7, 7) symmetry the low-energy theorem for the 1-loop n-point amplitudes.

Rube Lab

Mathematical Biology | Genomics | Biophysics

Systems and methods for analyzing sequencing data

Notes