Biophysical insights into OR2T7: Investigation of a potential prognostic marker for glioblastoma.

Glioblastoma multiforme (GBM) is the most aggressive and prevalent form of brain cancer, with an expected survival of 12-15 months following diagnosis. GBM affects the glial cells of the central nervous system, which impairs regular brain function including memory, hearing, and vision. GBM has virtually no long-term survival even with treatment, requiring novel strategies to understand disease progression. Here, we identified a somatic mutation in OR2T7, a G-protein-coupled receptor (GPCR), that correlates with reduced progression-free survival for glioblastoma (log rank p-value = 0.05), suggesting a possible role in tumor progression. The mutation, D125V, occurred in 10% of 396 glioblastoma samples in The Cancer Genome Atlas, but not in any of the 2504 DNA sequences in the 1000 Genomes Project, suggesting that the mutation may have a deleterious functional effect. In addition, transcriptome analysis showed that the p38α mitogen-activated protein kinase (MAPK), c-Fos, c-Jun, and JunB proto-oncogenes, and putative tumor suppressors RhoB and caspase-14 were underexpressed in glioblastoma samples with the D125V mutation (false discovery rate < 0.05). Molecular modeling and molecular dynamics simulations have provided preliminary structural insight and indicate a dynamic helical movement network that is influenced by the membrane-embedded, cytofacial-facing residue 125, demonstrating a possible obstruction of G-protein binding on the cytofacial exposed region. We show that the mutation impacts the "open" GPCR conformation, potentially affecting Gα-subunit binding and associated downstream activity. Overall, our findings suggest that the Val125 mutation in OR2T7 could affect glioblastoma progression by downregulating GPCR-p38 MAPK tumor-suppression pathways and impacting the biophysical characteristics of the structure that facilitates Gα-subunit binding. This study provides the theoretical basis for further experimental investigation required to confirm that the D125V mutation in OR2T7 is not a passenger mutation. With validation, the aforementioned mutation could represent an important prognostic marker and a potential therapeutic target for glioblastoma.

CAGm: a repository of germline microsatellite variations in the 1000 genomes project.

The human genome harbors an abundance of repetitive DNA; however, its function continues to be debated. Microsatellites-a class of short tandem repeat-are established as an important source of genetic variation. Array length variants are common among microsatellites and affect gene expression; but, efforts to understand the role and diversity of microsatellite variation has been hampered by several challenges. Without adequate depth, both long-read and short-read sequencing may not detect the variants present in a sample; additionally, large sample sizes are needed to reveal the degree of population-level polymorphism. To address these challenges we present the Comparative Analysis of Germline Microsatellites (CAGm): a database of germline microsatellites from 2529 individuals in the 1000 genomes project. A key novelty of CAGm is the ability to aggregate microsatellite variation by population, ethnicity (super population) and gender. The database provides advanced searching for microsatellites embedded in genes and functional elements. All data can be downloaded as Microsoft Excel spreadsheets. Two use-case scenarios are presented to demonstrate its utility: a mononucleotide (A) microsatellite at the BAT-26 locus and a dinucleotide (CA) microsatellite in the coding region of FGFRL1. CAGm is freely available at http://www.cagmdb.org/.

Estimating the number of genetic mutations (hits) required for carcinogenesis based on the distribution of somatic mutations.

Individual instances of cancer are primarily a result of a combination of a small number of genetic mutations (hits). Knowing the number of such mutations is a prerequisite for identifying specific combinations of carcinogenic mutations and understanding the etiology of cancer. We present a mathematical model for estimating the number of hits based on the distribution of somatic mutations. The model is fundamentally different from previous approaches, which are based on cancer incidence by age. Our somatic mutation based model is likely to be more robust than age-based models since it does not require knowing or accounting for the highly variable mutation rate, which can vary by over three orders of magnitude. In fact, we find that the number of somatic mutations at diagnosis is weakly correlated with age at cancer diagnosis, most likely due to the extreme variability in mutation rates between individuals. Comparing the distribution of somatic mutations predicted by our model to the actual distribution from 6904 tumor samples we estimate the number of hits required for carcinogenesis for 17 cancer types. We find that different cancer types exhibit distinct somatic mutational profiles corresponding to different numbers of hits. Why might different cancer types require different numbers of hits for carcinogenesis? The answer may provide insight into the unique etiology of different cancer types.

DNA sequencing of anatomy lab cadavers to provide hands-on precision medicine introduction to medical students.

BACKGROUND: Medical treatment informed by Precision Medicine is becoming a standard practice for many diseases, and patients are curious about the consequences of genomic variants in their genome. However, most medical students' understanding of Precision Medicine derives from classroom lectures. This format does little to foster an understanding for the potential and limitations of Precision Medicine. To close this gap, we implemented a hands-on Precision Medicine training program utilizing exome sequencing to prepare a clinical genetic report of cadavers studied in the anatomy lab. The program reinforces Precision Medicine related learning objectives for the Genetics curriculum. METHODS: Pre-embalmed blood samples and embalmed tissue were obtained from cadavers (donors) used in the anatomy lab. DNA was isolated and sequenced and illustrative genetic reports provided to the students. The reports were used to facilitate discussion with students on the implications of pathogenic genomic variants and the potential correlation of these variants in each "donor" with any anatomical anomalies identified during cadaver dissection. RESULTS: In 75% of cases, analysis of whole exome sequencing data identified a variant associated with increased risk for a disease/abnormal condition noted in the donor's cause of death or in the students' anatomical findings. This provided students with real-world examples of the potential relationship between genomic variants and disease risk. Our students also noted that diseases associated with 92% of the pathogenic variants identified were not related to the anatomical findings, demonstrating the limitations of Precision Medicine. CONCLUSION: With this study, we have established protocols and classroom procedures incorporating hands-on Precision Medicine training in the medical student curriculum and a template for other medical educators interested in enhancing their Precision Medicine training program. The program engaged students in discovering variants that were associated with the pathophysiology of the cadaver they were studying, which led to more exposure and understanding of the potential risks and benefits of genomic medicine.

Biophysical comparison of ATP synthesis mechanisms shows a kinetic advantage for the rotary process.

The ATP synthase (F-ATPase) is a highly complex rotary machine that synthesizes ATP, powered by a proton electrochemical gradient. Why did evolution select such an elaborate mechanism over arguably simpler alternating-access processes that can be reversed to perform ATP synthesis? We studied a systematic enumeration of alternative mechanisms, using numerical and theoretical means. When the alternative models are optimized subject to fundamental thermodynamic constraints, they fail to match the kinetic ability of the rotary mechanism over a wide range of conditions, particularly under low-energy conditions. We used a physically interpretable, closed-form solution for the steady-state rate for an arbitrary chemical cycle, which clarifies kinetic effects of complex free-energy landscapes. Our analysis also yields insights into the debated "kinetic equivalence" of ATP synthesis driven by transmembrane pH and potential difference. Overall, our study suggests that the complexity of the F-ATPase may have resulted from positive selection for its kinetic advantage.

Speed of conformational change: comparing explicit and implicit solvent molecular dynamics simulations.

Adequate sampling of conformation space remains challenging in atomistic simulations, especially if the solvent is treated explicitly. Implicit-solvent simulations can speed up conformational sampling significantly. We compare the speed of conformational sampling between two commonly used methods of each class: the explicit-solvent particle mesh Ewald (PME) with TIP3P water model and a popular generalized Born (GB) implicit-solvent model, as implemented in the AMBER package. We systematically investigate small (dihedral angle flips in a protein), large (nucleosome tail collapse and DNA unwrapping), and mixed (folding of a miniprotein) conformational changes, with nominal simulation times ranging from nanoseconds to microseconds depending on system size. The speedups in conformational sampling for GB relative to PME simulations, are highly system- and problem-dependent. Where the simulation temperatures for PME and GB are the same, the corresponding speedups are approximately onefold (small conformational changes), between ∼1- and ∼100-fold (large changes), and approximately sevenfold (mixed case). The effects of temperature on speedup and free-energy landscapes, which may differ substantially between the solvent models, are discussed in detail for the case of miniprotein folding. In addition to speeding up conformational sampling, due to algorithmic differences, the implicit solvent model can be computationally faster for small systems or slower for large systems, depending on the number of solute and solvent atoms. For the conformational changes considered here, the combined speedups are approximately twofold, ∼1- to 60-fold, and ∼50-fold, respectively, in the low solvent viscosity regime afforded by the implicit solvent. For all the systems studied, 1) conformational sampling speedup increases as Langevin collision frequency (effective viscosity) decreases; and 2) conformational sampling speedup is mainly due to reduction in solvent viscosity rather than possible differences in free-energy landscapes between the solvent models.

H++ 3.0: automating pK prediction and the preparation of biomolecular structures for atomistic molecular modeling and simulations.

The accuracy of atomistic biomolecular modeling and simulation studies depend on the accuracy of the input structures. Preparing these structures for an atomistic modeling task, such as molecular dynamics (MD) simulation, can involve the use of a variety of different tools for: correcting errors, adding missing atoms, filling valences with hydrogens, predicting pK values for titratable amino acids, assigning predefined partial charges and radii to all atoms, and generating force field parameter/topology files for MD. Identifying, installing and effectively using the appropriate tools for each of these tasks can be difficult for novice and time-consuming for experienced users. H++ (http://biophysics.cs.vt.edu/) is a free open-source web server that automates the above key steps in the preparation of biomolecular structures for molecular modeling and simulations. H++ also performs extensive error and consistency checking, providing error/warning messages together with the suggested corrections. In addition to numerous minor improvements, the latest version of H++ includes several new capabilities and options: fix erroneous (flipped) side chain conformations for HIS, GLN and ASN, include a ligand in the input structure, process nucleic acid structures and generate a solvent box with specified number of common ions for explicit solvent MD.

Potential immunosuppressive clonal hematopoietic mutations in tumor infiltrating immune cells in breast invasive carcinoma.

A hallmark of cancer is a tumor cell's ability to evade immune destruction. Somatic mutations in tumor cells that prevent immune destruction have been extensively studied. However, somatic mutations in tumor infiltrating immune (TII) cells, to our knowledge, have not been previously studied. Understandably so since normal hematopoiesis prevents the accumulation of somatic mutations in immune cells. However, clonal hematopoiesis does result in the accumulation of somatic mutations in immune cells. These mutations cannot "drive" tumor growth, however, they may "facilitate" it by inhibiting an effective anti-tumor immune response. To identify potential immunosuppressive clonal hematopoietic (CH) mutations in TII cells, we analyzed exome and RNA sequencing data from matched tumor and normal blood samples, and single-cell RNA sequencing data, from breast cancer patients. We selected mutations that were somatic, present in TII cells, clonally expanded, potentially pathogenic, expressed in TII cells, unlikely to be a passenger mutation, and in immune response associated genes. We identified eight potential immunosuppressive CH mutations in TII cells. This work is a first step towards determining if immunosuppressive CH mutations in TII cells can affect the progression of solid tumors. Subsequent experimental confirmation could represent a new paradigm in the etiology of cancer.

Cranial manipulation affects cholinergic pathway gene expression in aged rats.

CONTEXT: Age-dependent dementia is a devastating disorder afflicting a growing older population. Although pharmacological agents improve symptoms of dementia, age-related comorbidities combined with adverse effects often outweigh their clinical benefits. Therefore, nonpharmacological therapies are being investigated as an alternative. In a previous pilot study, aged rats demonstrated improved spatial memory after osteopathic cranial manipulative medicine (OCMM) treatment. OBJECTIVES: In this continuation of the pilot study, we examine the effect of OCMM on gene expression to elicit possible explanations for the improvement in spatial memory. METHODS: OCMM was performed on six of 12 elderly rats every day for 7 days. Rats were then euthanized to obtain the brain tissue, from which RNA samples were extracted. RNA from three treated and three controls were of sufficient quality for sequencing. These samples were sequenced utilizing next-generation sequencing from Illumina NextSeq. The Cufflinks software suite was utilized to assemble transcriptomes and quantify the RNA expression level for each sample. RESULTS: Transcriptome analysis revealed that OCMM significantly affected the expression of 36 genes in the neuronal pathway (false discovery rate [FDR] <0.004). The top five neuronal genes with the largest-fold change were part of the cholinergic neurotransmission mechanism, which is known to affect cognitive function. In addition, 39.9% of 426 significant differentially expressed (SDE) genes (FDR<0.004) have been previously implicated in neurological disorders. Overall, changes in SDE genes combined with their role in central nervous system signaling pathways suggest a connection to previously reported OCMM-induced behavioral and biochemical changes in aged rats. CONCLUSIONS: Results from this pilot study provide sufficient evidence to support a more extensive study with a larger sample size. Further investigation in this direction will provide a better understanding of the molecular mechanisms of OCMM and its potential in clinical applications. With clinical validation, OCMM could represent a much-needed low-risk adjunct treatment for age-related dementia including Alzheimer's disease.

Ethnically biased microsatellites contribute to differential gene expression and glutathione metabolism in Africans and Europeans.

Approximately three percent of the human genome is occupied by microsatellites: a type of short tandem repeat (STR). Microsatellites have well established effects on (a) the genetic structure of diverse human populations and (b) expression of nearby genes. These lines of inquiry have uncovered 3,984 ethnically biased microsatellite loci (EBML) and 28,375 expression STRs (eSTRs), respectively. We hypothesize that a combination of EBML, eSTRs, and gene expression data (RNA-seq) can be used to show that microsatellites contribute to differential gene expression and phenotype in human populations. In fact, our previous study demonstrated a degree of mutual overlap between EBML and eSTRs but fell short of quantifying effects on gene expression. The present work aims to narrow the gap. First, we identify 313 overlapping EBML/eSTRs and recapitulate their mutual overlap. The 313 EBML/eSTRs are then characterized across ethnicity and tissue type. We use RNA-seq data to pursue validation of 49 regions that affect whole blood gene expression; 32 out of 54 affected genes are differentially expressed in Africans and Europeans. We quantify the relative contribution of these 32 genes to differential expression; fold change tends to be less than other differentially expressed genes. Repeat length correlates with expression for 15 of the 32 genes; two are conspicuously involved in glutathione metabolism. Finally, we repurpose a mathematical model of glutathione metabolism to investigate how a single polymorphic microsatellite affects phenotype. We conclude with a testable prediction that microsatellite polymorphisms affect GPX7 expression and oxidative stress in Africans and Europeans.

Statistics and physical origins of pK and ionization state changes upon protein-ligand binding.

This work investigates statistical prevalence and overall physical origins of changes in charge states of receptor proteins upon ligand binding. These changes are explored as a function of the ligand type (small molecule, protein, and nucleic acid), and distance from the binding region. Standard continuum solvent methodology is used to compute, on an equal footing, pK changes upon ligand binding for a total of 5899 ionizable residues in 20 protein-protein, 20 protein-small molecule, and 20 protein-nucleic acid high-resolution complexes. The size of the data set combined with an extensive error and sensitivity analysis allows us to make statistically justified and conservative conclusions: in 60% of all protein-small molecule, 90% of all protein-protein, and 85% of all protein-nucleic acid complexes there exists at least one ionizable residue that changes its charge state upon ligand binding at physiological conditions (pH = 6.5). Considering the most biologically relevant pH range of 4-8, the number of ionizable residues that experience substantial pK changes (DeltapK > 1.0) due to ligand binding is appreciable: on average, 6% of all ionizable residues in protein-small molecule complexes, 9% in protein-protein, and 12% in protein-nucleic acid complexes experience a substantial pK change upon ligand binding. These changes are safely above the statistical false-positive noise level. Most of the changes occur in the immediate binding interface region, where approximately one out of five ionizable residues experiences substantial pK change regardless of the ligand type. However, the physical origins of the change differ between the types: in protein-nucleic acid complexes, the pK values of interface residues are predominantly affected by electrostatic effects, whereas in protein-protein and protein-small molecule complexes, structural changes due to the induced-fit effect play an equally important role. In protein-protein and protein-nucleic acid complexes, there is a statistically significant number of substantial pK perturbations, mostly due to the induced-fit structural changes, in regions far from the binding interface.

Functional bias in molecular evolution rate of Arabidopsis thaliana.

BACKGROUND: Characteristics derived from mutation and other mechanisms that are advantageous for survival are often preserved during evolution by natural selection. Some genes are conserved in many organisms because they are responsible for fundamental biological function, others are conserved for their unique functional characteristics. Therefore one would expect the rate of molecular evolution for individual genes to be dependent on their biological function. Whether this expectation holds for genes duplicated by whole genome duplication is not known. RESULTS: We empirically demonstrate here, using duplicated genes generated from the Arabidopsis thaliana alpha-duplication event, that the rate of molecular evolution of genes duplicated in this event depend on biological function. Using functional clustering based on gene ontology annotation of gene pairs, we show that some duplicated genes, such as defense response genes, are under weaker purifying selection or under stronger diversifying selection than other duplicated genes, such as protein translation genes, as measured by the ratio of nonsynonymous to synonymous divergence (dN/dS). CONCLUSIONS: These results provide empirical evidence indicating that molecular evolution rate for genes duplicated in whole genome duplication, as measured by dN/dS, may depend on biological function, which we characterize using gene ontology annotation. Furthermore, the general approach used here provides a framework for comparative analysis of molecular evolution rate for genes based on their biological function.

An N log N approximation based on the natural organization of biomolecules for speeding up the computation of long range interactions.

Presented here is a method, the hierarchical charge partitioning (HCP) approximation, for speeding up computation of pairwise electrostatic interactions in biomolecular systems. The approximation is based on multiple levels of natural partitioning of biomolecular structures into a hierarchical set of its constituent structural components. The charge distribution in each component is systematically approximated by a small number of point charges, which, for the highest level component, are much fewer than the number of atoms in the component. For short distances from the point of interest, the HCP uses the full set of atomic charges available. For long-distance interactions, the approximate charge distributions with smaller sets of charges are used instead. For a structure consisting of N charges, the computational cost of computing the pairwise interactions via the HCP scales as O(N log N), under assumptions about the structural organization of biomolecular structures generally consistent with reality. A proof-of-concept implementation of the HCP shows that for large structures it can lead to speed-up factors of up to several orders of magnitude relative to the exact pairwise O(N(2)) all-atom computation used as a reference. For structures with more than 2000-3000 atoms the relative accuracy of the HCP (relative root-mean-square force error per atom), approaches the accuracy of the particle mesh Ewald (PME) method with parameter settings typical for biomolecular simulations. When averaged over a set of 600 representative biomolecular structures, the relative accuracies of the two methods are roughly equal. The HCP is also significantly more accurate than the spherical cutoff method. The HCP has been implemented in the freely available nucleic acids builder (NAB) molecular dynamics (MD) package in Amber tools. A 10 ns simulation of a small protein indicates that the HCP based MD simulation is stable, and that it can be faster than the spherical cutoff method. A critical benefit of the HCP approximation is that it is algorithmically very simple, and unlike the PME, the HCP is straightforward to use with implicit solvent models.

Crossing complexity of space-filling curves reveals entanglement of S-phase DNA.

Space-filling curves have been used for decades to study the folding principles of globular proteins, compact polymers, and chromatin. Formally, space-filling curves trace a single circuit through a set of points (x,y,z); informally, they correspond to a polymer melt. Although not quite a melt, the folding principles of Human chromatin are likened to the Hilbert curve: a type of space-filling curve. Hilbert-like curves in general make biologically compelling models of chromatin; in particular, they lack knots which facilitates chromatin folding, unfolding, and easy access to genes. Knot complexity has been intensely studied with the aid of Alexander polynomials; however, the approach does not generalize well to cases of more than one chromosome. Crossing complexity is an understudied alternative better suited for quantifying entanglement between chromosomes. Do Hilbert-like configurations limit crossing complexity between chromosomes? How does crossing complexity for Hilbert-like configurations compare to equilibrium configurations? To address these questions, we extend the Mansfield algorithm to enable sampling of Hilbert-like space filling curves on a simple cubic lattice. We use the extended algorithm to generate equilibrium, intermediate, and Hilbert-like configurational ensembles and compute crossing complexity between curves (chromosomes) in each configurational snapshot. Our main results are twofold: (a) Hilbert-like configurations limit entanglement between chromosomes and (b) Hilbert-like configurations do not limit entanglement in a model of S-phase DNA. Our second result is particularly surprising yet easily rationalized with a geometric argument. We explore ergodicity of the extended algorithm and discuss our results in the context of more sophisticated models of chromatin.

Identifying multi-hit carcinogenic gene combinations: Scaling up a weighted set cover algorithm using compressed binary matrix representation on a GPU.

Despite decades of research, effective treatments for most cancers remain elusive. One reason is that different instances of cancer result from different combinations of multiple genetic mutations (hits). Therefore, treatments that may be effective in some cases are not effective in others. We previously developed an algorithm for identifying combinations of carcinogenic genes with mutations (multi-hit combinations), which could suggest a likely cause for individual instances of cancer. Most cancers are estimated to require three or more hits. However, the computational complexity of the algorithm scales exponentially with the number of hits, making it impractical for identifying combinations of more than two hits. To identify combinations of greater than two hits, we used a compressed binary matrix representation, and optimized the algorithm for parallel execution on an NVIDIA V100 graphics processing unit (GPU). With these enhancements, the optimized GPU implementation was on average an estimated 12,144 times faster than the original integer matrix based CPU implementation, for the 3-hit algorithm, allowing us to identify 3-hit combinations. The 3-hit combinations identified using a training set were able to differentiate between tumor and normal samples in a separate test set with 90% overall sensitivity and 93% overall specificity. We illustrate how the distribution of mutations in tumor and normal samples in the multi-hit gene combinations can suggest potential driver mutations for further investigation. With experimental validation, these combinations may provide insight into the etiology of cancer and a rational basis for targeted combination therapy.

Why Computed Protein Folding Landscapes Are Sensitive to the Water Model.

We investigate the effect of solvent models on the computed thermodynamics of protein folding. Atomistic folding simulations of a fast-folding mini-protein, CLN025, were employed to compare two commonly used explicit solvent water models, TIP3P and TIP4P/Ew, and one implicit solvent (AMBER generalized Born) model. Although all three solvent models correctly identify the same native folded state (RMSD = 1.5 ± 0.1 Å relative to the experimental structure), the corresponding free energy landscapes vary drastically between water models: almost an order-of-magnitude difference is seen in the predicted fraction of the unfolded state between the two explicit solvent models, with even larger differences between the implicit and the explicit models. Quantitative arguments are presented for why the sensitivity is expected to hold for other proteins, as well as for other conformational transitions involving large changes in solvent exposed areas such as protein-ligand binding. Comparing protein-solvent and solvent-solvent contributions to the folding energy between different water models, water-water electrostatic interactions are identified as the largest contributor to the differences in the predicted folding energy, which helps explain the strong sensitivity of the folding landscape to subtle details of the water model. For the two explicit solvent models, differences in water model parameters also result in the average number of water molecules surrounding the protein being noticeably different. Water models that poorly reproduce certain bulk properties of liquid water such as self-diffusion are likely to misrepresent water-water interactions; we argue that within a pairwise additive energy function this error cannot, in general, be compensated by an adjustment to the solute-solute and solute-solvent parts of the energy.

Abundance of ethnically biased microsatellites in human gene regions.

Microsatellites-a type of short tandem repeat (STR)-have been used for decades as putatively neutral markers to study the genetic structure of diverse human populations. However, recent studies have demonstrated that some microsatellites contribute to gene expression, cis heritability, and phenotype. As a corollary, some microsatellites may contribute to differential gene expression and RNA/protein structure stability in distinct human populations. To test this hypothesis, we investigate genotype frequencies, functional relevance, and adaptive potential of microsatellites in five super-populations (ethnicities) drawn from the 1000 Genomes Project. We discover 3,984 ethnically-biased microsatellite loci (EBML); for each EBML at least one ethnicity has genotype frequencies statistically different from the remaining four. South Asian, East Asian, European, and American EBML show significant overlap; on the contrary, the set of African EBML is mostly unique. We cross-reference the 3,984 EBML with 2,060 previously identified expression STRs (eSTRs); repeats known to affect gene expression (64 total) are over-represented. The most significant pathway enrichments are those associated with the matrisome: a broad collection of genes encoding the extracellular matrix and its associated proteins. At least 14 of the EBML have established links to human disease. Analysis of the 3,984 EBML with respect to known selective sweep regions in the genome shows that allelic variation in some of them is likely associated with adaptive evolution.

Author Correction: Differentiating between cancer and normal tissue samples using multi-hit combinations of genetic mutations.

An amendment to this paper has been published and can be accessed via a link at the top of the paper.

Differentiating between cancer and normal tissue samples using multi-hit combinations of genetic mutations.

Cancer is known to result from a combination of a small number of genetic defects. However, the specific combinations of mutations responsible for the vast majority of cancers have not been identified. Current computational approaches focus on identifying driver genes and mutations. Although individually these mutations can increase the risk of cancer they do not result in cancer without additional mutations. We present a fundamentally different approach for identifying the cause of individual instances of cancer: we search for combinations of genes with carcinogenic mutations (multi-hit combinations) instead of individual driver genes or mutations. We developed an algorithm that identified a set of multi-hit combinations that differentiate between tumor and normal tissue samples with 91% sensitivity (95% Confidence Interval (CI) = 89-92%) and 93% specificity (95% CI = 91-94%) on average for seventeen cancer types. We then present an approach based on mutational profile that can be used to distinguish between driver and passenger mutations within these genes. These combinations, with experimental validation, can aid in better diagnosis, provide insights into the etiology of cancer, and provide a rational basis for designing targeted combination therapies.