A cost-effective sequencing method for genetic studies combining high-depth whole exome and low-depth whole genome.

Whole genome sequencing (WGS) at high-depth (30X) allows the accurate discovery of variants in the coding and non-coding DNA regions and helps elucidate the genetic underpinnings of human health and diseases. Yet, due to the prohibitive cost of high-depth WGS, most large-scale genetic association studies use genotyping arrays or high-depth whole exome sequencing (WES). Here we propose a cost-effective method which we call "Whole Exome Genome Sequencing" (WEGS), that combines low-depth WGS and high-depth WES with up to 8 samples pooled and sequenced simultaneously (multiplexed). We experimentally assess the performance of WEGS with four different depth of coverage and sample multiplexing configurations. We show that the optimal WEGS configurations are 1.7-2.0 times cheaper than standard WES (no-plexing), 1.8-2.1 times cheaper than high-depth WGS, reach similar recall and precision rates in detecting coding variants as WES, and capture more population-specific variants in the rest of the genome that are difficult to recover when using genotype imputation methods. We apply WEGS to 862 patients with peripheral artery disease and show that it directly assesses more known disease-associated variants than a typical genotyping array and thousands of non-imputable variants per disease-associated locus.

Association of novel mutation in TRPV4 with familial nonsyndromic craniosynostosis with complete penetrance and variable expressivity.

OBJECTIVE: The aim of this study was to characterize a novel pathogenic variant in the transient receptor potential vanilloid 4 (TRPV4) gene, causing familial nonsyndromic craniosynostosis (CS) with complete penetrance and variable expressivity. METHODS: Whole-exome sequencing was performed on germline DNA of a family with nonsyndromic CS to a mean depth coverage of 300× per sample, with greater than 98% of the targeted region covered at least 25×. In this study, the authors detected a novel variant, c.496C>A in TRPV4, exclusively in the four affected family members. The variant was modeled using the structure of the TRPV4 protein from Xenopus tropicalis. In vitro assays in HEK293 cells overexpressing wild-type TRPV4 or TRPV4 p.Leu166Met were used to assess the effect of the mutation on channel activity and downstream MAPK signaling. RESULTS: The authors identified a novel, highly penetrant heterozygous variant in TRPV4 (NM_021625.4:c.496C>A) causing nonsyndromic CS in a mother and all three of her children. This variant results in an amino acid change (p.Leu166Met) in the intracellular ankyrin repeat domain distant from the Ca2+-dependent membrane channel domain. In contrast to other TRPV4 mutations in channelopathies, this variant does not interfere with channel activity as identified by in silico modeling and in vitro overexpression assays in HEK293 cells. CONCLUSIONS: Based on these findings, the authors hypothesized that this novel variant causes CS by modulating the binding of allosteric regulatory factors to TRPV4 rather than directly modifying its channel activity. Overall, this study expands the genetic and functional spectrum of TRPV4 channelopathies and is particularly relevant for the genetic counseling of CS patients.

Application of Genomic Sequencing to Refine Patient Stratification for Adjuvant Therapy in Renal Cell Carcinoma.

PURPOSE: Patients with resected localized clear-cell renal cell carcinoma (ccRCC) remain at variable risk of recurrence. Incorporation of biomarkers may refine risk prediction and inform adjuvant treatment decisions. We explored the role of tumor genomics in this setting, leveraging the largest cohort to date of localized ccRCC tissues subjected to targeted gene sequencing. EXPERIMENTAL DESIGN: The somatic mutation status of 12 genes was determined in 943 ccRCC cases from a multinational cohort of patients, and associations to outcomes were examined in a Discovery (n = 469) and Validation (n = 474) framework. RESULTS: Tumors containing a von-Hippel Lindau (VHL) mutation alone were associated with significantly improved outcomes in comparison with tumors containing a VHL plus additional mutations. Within the Discovery cohort, those with VHL+0, VHL+1, VHL+2, and VHL+≥3 tumors had disease-free survival (DFS) rates of 90.8%, 80.1%, 68.2%, and 50.7% respectively, at 5 years. This trend was replicated in the Validation cohort. Notably, these genomically defined groups were independent of tumor mutational burden. Amongst patients eligible for adjuvant therapy, those with a VHL+0 tumor (29%) had a 5-year DFS rate of 79.3% and could, therefore, potentially be spared further treatment. Conversely, patients with VHL+2 and VHL+≥3 tumors (32%) had equivalent DFS rates of 45.6% and 35.3%, respectively, and should be prioritized for adjuvant therapy. CONCLUSIONS: Genomic characterization of ccRCC identified biologically distinct groups of patients with divergent relapse rates. These groups account for the ∼80% of cases with VHL mutations and could be used to personalize adjuvant treatment discussions with patients as well as inform future adjuvant trial design.

PrecisionFDA Truth Challenge V2: Calling variants from short and long reads in difficult-to-map regions.

The precisionFDA Truth Challenge V2 aimed to assess the state of the art of variant calling in challenging genomic regions. Starting with FASTQs, 20 challenge participants applied their variant-calling pipelines and submitted 64 variant call sets for one or more sequencing technologies (Illumina, PacBio HiFi, and Oxford Nanopore Technologies). Submissions were evaluated following best practices for benchmarking small variants with updated Genome in a Bottle benchmark sets and genome stratifications. Challenge submissions included numerous innovative methods, with graph-based and machine learning methods scoring best for short-read and long-read datasets, respectively. With machine learning approaches, combining multiple sequencing technologies performed particularly well. Recent developments in sequencing and variant calling have enabled benchmarking variants in challenging genomic regions, paving the way for the identification of previously unknown clinically relevant variants.

TreeTuner: A pipeline for minimizing redundancy and complexity in large phylogenetic datasets.

Various bioinformatics protocols have been developed for trimming the number of operational taxonomic units (OTUs) in phylogenetic datasets, but they typically require significant manual intervention. Here we present TreeTuner, a semiautomated pipeline that allows both coarse and fine-scale tuning of large protein sequence phylogenetic datasets via the minimization of OTU redundancy. TreeTuner facilitates preliminary investigation of such datasets as well as more rigorous downstream analysis of specific subsets of OTUs. For complete details on the use and execution of this protocol, please refer to Maruyama et al. (2013) and Sibbald et al. (2019).

HLA-B*07:02 and HLA-C*07:02 are associated with trimethoprim-sulfamethoxazole respiratory failure.

We have identified an underrecognized severe adverse drug reaction (ADR) of trimethoprim-sulfamethoxazole (TMP-SMX) associated respiratory failure in previously healthy children and young adults. We investigated potential genetic risk factors associated with TMP-SMX induced respiratory failure in a cohort of seven patients. We explored whole genome sequence among seven patients representing nearly half of all reported cases worldwide and 63 unrelated control individuals in two stages: (1) human leukocyte antigen (HLA) locus variation as several other ADRs have been associated HLA genetic variants and (2) coding variation to catalog and explore potential rare variants contributing to this devastating reaction. All cases were either heterozygous (carriers) or homozygous for the common HLA-B*07:02-HLA-C*07:02 haplotype. Despite the small sample size, this observation is statistically significant both in conservative comparison to maximum reported population frequencies (binomial P = 0.00017 for HLA-B and P = 0.00028 for HLA-C) and to our control population assessed by same HLA genotyping approach (binomial P = 0.000001 for HLA-B and P = 0.000018 for HLA-C). No gene elsewhere in the genome harnessed shared rare case enriched coding variation. Our results suggests that HLA-B*07:02 and HLA-C*07:02 are necessary for a patient to develop respiratory failure due to TMP-SMX.

De novo TRPV4 Leu619Pro variant causes a new channelopathy characterised by giant cell lesions of the jaws and skull, skeletal abnormalities and polyneuropathy.

BACKGROUND: Pathogenic germline variants in Transient Receptor Potential Vanilloid 4 Cation Channel (TRPV4) lead to channelopathies, which are phenotypically diverse and heterogeneous disorders grossly divided in neuromuscular disorders and skeletal dysplasia. We recently reported in sporadic giant cell lesions of the jaws (GCLJs) novel, somatic, heterozygous, gain-of-function mutations in TRPV4, at Met713. METHODS: Here we report two unrelated women with a de novo germline p.Leu619Pro TRPV4 variant and an overlapping systemic disorder affecting all organs individually described in TRPV4 channelopathies. RESULTS: From an early age, both patients had several lesions of the nervous system including progressive polyneuropathy, and multiple aggressive giant cell-rich lesions of the jaws and craniofacial/skull bones, and other skeletal lesions. One patient had a relatively milder disease phenotype possibly due to postzygotic somatic mosaicism. Indeed, the TRPV4 p.Leu619Pro variant was present at a lower frequency (variant allele frequency (VAF)=21.6%) than expected for a heterozygous variant as seen in the other proband, and showed variable regional frequency in the GCLJ (VAF ranging from 42% to 10%). In silico structural analysis suggests that the gain-of-function p.Leu619Pro alters the ion channel activity leading to constitutive ion leakage. CONCLUSION: Our findings define a novel polysystemic syndrome due to germline TRPV4 p.Leu619Pro and further extend the spectrum of TRPV4 channelopathies. They further highlight the convergence of TRPV4 mutations on different organ systems leading to complex phenotypes which are further mitigated by possible post-zygotic mosaicism. Treatment of this disorder is challenging, and surgical intervention of the GCLJ worsens the lesions, suggesting the future use of MEK inhibitors and TRPV4 antagonists as therapeutic modalities for unmet clinical needs.

Rare loss-of-function variants in type I IFN immunity genes are not associated with severe COVID-19.

A recent report found that rare predicted loss-of-function (pLOF) variants across 13 candidate genes in TLR3- and IRF7-dependent type I IFN pathways explain up to 3.5% of severe COVID-19 cases. We performed whole-exome or whole-genome sequencing of 1,864 COVID-19 cases (713 with severe and 1,151 with mild disease) and 15,033 ancestry-matched population controls across 4 independent COVID-19 biobanks. We tested whether rare pLOF variants in these 13 genes were associated with severe COVID-19. We identified only 1 rare pLOF mutation across these genes among 713 cases with severe COVID-19 and observed no enrichment of pLOFs in severe cases compared to population controls or mild COVID-19 cases. We found no evidence of association of rare LOF variants in the 13 candidate genes with severe COVID-19 outcomes.

A Distributed Whole Genome Sequencing Benchmark Study.

Population sequencing often requires collaboration across a distributed network of sequencing centers for the timely processing of thousands of samples. In such massive efforts, it is important that participating scientists can be confident that the accuracy of the sequence data produced is not affected by which center generates the data. A study was conducted across three established sequencing centers, located in Montreal, Toronto, and Vancouver, constituting Canada's Genomics Enterprise (www.cgen.ca). Whole genome sequencing was performed at each center, on three genomic DNA replicates from three well-characterized cell lines. Secondary analysis pipelines employed by each site were applied to sequence data from each of the sites, resulting in three datasets for each of four variables (cell line, replicate, sequencing center, and analysis pipeline), for a total of 81 datasets. These datasets were each assessed according to multiple quality metrics including concordance with benchmark variant truth sets to assess consistent quality across all three conditions for each variable. Three-way concordance analysis of variants across conditions for each variable was performed. Our results showed that the variant concordance between datasets differing only by sequencing center was similar to the concordance for datasets differing only by replicate, using the same analysis pipeline. We also showed that the statistically significant differences between datasets result from the analysis pipeline used, which can be unified and updated as new approaches become available. We conclude that genome sequencing projects can rely on the quality and reproducibility of aggregate data generated across a network of distributed sites.

H3.3 G34W Promotes Growth and Impedes Differentiation of Osteoblast-Like Mesenchymal Progenitors in Giant Cell Tumor of Bone.

Glycine 34-to-tryptophan (G34W) substitutions in H3.3 arise in approximately 90% of giant cell tumor of bone (GCT). Here, we show H3.3 G34W is necessary for tumor formation. By profiling the epigenome, transcriptome, and secreted proteome of patient samples and tumor-derived cells CRISPR-Cas9-edited for H3.3 G34W, we show that H3.3K36me3 loss on mutant H3.3 alters the deposition of the repressive H3K27me3 mark from intergenic to genic regions, beyond areas of H3.3 deposition. This promotes redistribution of other chromatin marks and aberrant transcription, altering cell fate in mesenchymal progenitors and hindering differentiation. Single-cell transcriptomics reveals that H3.3 G34W stromal cells recapitulate a neoplastic trajectory from a SPP1+ osteoblast-like progenitor population toward an ACTA2+ myofibroblast-like population, which secretes extracellular matrix ligands predicted to recruit and activate osteoclasts. Our findings suggest that H3.3 G34W leads to GCT by sustaining a transformed state in osteoblast-like progenitors, which promotes neoplastic growth, pathologic recruitment of giant osteoclasts, and bone destruction. SIGNIFICANCE: This study shows that H3.3 G34W drives GCT tumorigenesis through aberrant epigenetic remodeling, altering differentiation trajectories in mesenchymal progenitors. H3.3 G34W promotes in neoplastic stromal cells an osteoblast-like progenitor state that enables undue interactions with the tumor microenvironment, driving GCT pathogenesis. These epigenetic changes may be amenable to therapeutic targeting in GCT.See related commentary by Licht, p. 1794.This article is highlighted in the In This Issue feature, p. 1775.

Latency and interval therapy affect the evolution in metastatic colorectal cancer.

While comparison of primary tumor and metastases has highlighted genomic heterogeneity in colorectal cancer (CRC), previous studies have focused on a single metastatic site or limited genomic testing. Combining data from whole exome and ultra-deep targeted sequencing, we explored possible evolutionary trajectories beyond the status of these mutations, particularly among patient-matched metastatic tumors. Our findings confirm the persistence of known clinically-relevant mutations (e.g., those of RAS family of oncogenes) in CRC primary and metastases, yet reveal that latency and interval systemic therapy affect the course of evolutionary events within metastatic lesions. Specifically, our analysis of patient-matched primary and multiple metastatic lesions, developed over time, showed a similar genetic composition for liver metastatic tumors, which were 21-months apart. This genetic makeup was different from those identified in lung metastases developed before manifestation of the second liver metastasis. These results underscore the role of latency in the evolutionary path of metastatic CRC and may have implications for future treatment options.

Failure to replicate the association of rare loss-of-function variants in type I IFN immunity genes with severe COVID-19.

A recent report found that rare predicted loss-of-function (pLOF) variants across 13 candidate genes in TLR3- and IRF7-dependent type I IFN pathways explain up to 3.5% of severe COVID-19 cases. We performed whole-exome or whole-genome sequencing of 1,934 COVID-19 cases (713 with severe and 1,221 with mild disease) and 15,251 ancestry-matched population controls across four independent COVID-19 biobanks. We then tested if rare pLOF variants in these 13 genes were associated with severe COVID-19. We identified only one rare pLOF mutation across these genes amongst 713 cases with severe COVID-19 and observed no enrichment of pLOFs in severe cases compared to population controls or mild COVID-19 cases. We find no evidence of association of rare loss-of-function variants in the proposed 13 candidate genes with severe COVID-19 outcomes.

A point mutation in the linker domain of mouse STAT5A is associated with impaired NK-cell regulation.

The transcription factor STAT5 is critical for peripheral NK-cell survival, proliferation, and cytotoxic function. STAT5 refers to two highly related proteins, STAT5A and STAT5B. In this study, we verified the importance of STAT5A isoform for NK cells. We characterized an incidental chemically induced W484G mutation in the Stat5a gene and found that this mutation was associated with a reduction of STAT5A protein expression. Closer examination of NK-cell subsets from Stat5a mutant mice showed marked reductions in NK-cell number and maturation. IL-15 treatment of Stat5a mutant NK cells exhibited defective induction of both STAT5 and mTOR signaling pathways and reduced expression of granzyme B and IFN-γ. Finally, we observed that Stat5a mutant mice revealed more tumor growth upon injection of RMA-S tumor cell line. Overall, our results demonstrate that the W484G mutation in the linker domain of STAT5A is sufficient to compromise STAT5A function in NK-cell homeostasis, responsiveness, and tumoricidal function.

GenPipes: an open-source framework for distributed and scalable genomic analyses.

BACKGROUND: With the decreasing cost of sequencing and the rapid developments in genomics technologies and protocols, the need for validated bioinformatics software that enables efficient large-scale data processing is growing. FINDINGS: Here we present GenPipes, a flexible Python-based framework that facilitates the development and deployment of multi-step workflows optimized for high-performance computing clusters and the cloud. GenPipes already implements 12 validated and scalable pipelines for various genomics applications, including RNA sequencing, chromatin immunoprecipitation sequencing, DNA sequencing, methylation sequencing, Hi-C, capture Hi-C, metagenomics, and Pacific Biosciences long-read assembly. The software is available under a GPLv3 open source license and is continuously updated to follow recent advances in genomics and bioinformatics. The framework has already been configured on several servers, and a Docker image is also available to facilitate additional installations. CONCLUSIONS: GenPipes offers genomics researchers a simple method to analyze different types of data, customizable to their needs and resources, as well as the flexibility to create their own workflows.

Whole Exome Sequencing Reveals the Major Genetic Contributors to Nonsyndromic Tetralogy of Fallot.

RATIONALE: Familial recurrence studies provide strong evidence for a genetic component to the predisposition to sporadic, nonsyndromic Tetralogy of Fallot (TOF), the most common cyanotic congenital heart disease phenotype. Rare genetic variants have been identified as important contributors to the risk of congenital heart disease, but relatively small numbers of TOF cases have been studied to date. OBJECTIVE: We used whole exome sequencing to assess the prevalence of unique, deleterious variants in the largest cohort of nonsyndromic TOF patients reported to date. METHODS AND RESULTS: Eight hundred twenty-nine TOF patients underwent whole exome sequencing. The presence of unique, deleterious variants was determined; defined by their absence in the Genome Aggregation Database and a scaled combined annotation-dependent depletion score of ≥20. The clustering of variants in 2 genes, NOTCH1 and FLT4, surpassed thresholds for genome-wide significance (assigned as P<5×10-8) after correction for multiple comparisons. NOTCH1 was most frequently found to harbor unique, deleterious variants. Thirty-one changes were observed in 37 probands (4.5%; 95% CI, 3.2%-6.1%) and included 7 loss-of-function variants 22 missense variants and 2 in-frame indels. Sanger sequencing of the unaffected parents of 7 cases identified 5 de novo variants. Three NOTCH1 variants (p.G200R, p.C607Y, and p.N1875S) were subjected to functional evaluation, and 2 showed a reduction in Jagged1-induced NOTCH signaling. FLT4 variants were found in 2.4% (95% CI, 1.6%-3.8%) of TOF patients, with 21 patients harboring 22 unique, deleterious variants. The variants identified were distinct to those that cause the congenital lymphoedema syndrome Milroy disease. In addition to NOTCH1, FLT4 and the well-established TOF gene, TBX1, we identified potential association with variants in several other candidates, including RYR1, ZFPM1, CAMTA2, DLX6, and PCM1. CONCLUSIONS: The NOTCH1 locus is the most frequent site of genetic variants predisposing to nonsyndromic TOF, followed by FLT4. Together, variants in these genes are found in almost 7% of TOF patients.

Identification of Elongated Primary Cilia with Impaired Mechanotransduction in Idiopathic Scoliosis Patients.

The primary cilium is an outward projecting antenna-like organelle with an important role in bone mechanotransduction. The capacity to sense mechanical stimuli can affect important cellular and molecular aspects of bone tissue. Idiopathic scoliosis (IS) is a complex pediatric disease of unknown cause, defined by abnormal spinal curvatures. We demonstrate significant elongation of primary cilia in IS patient bone cells. In response to mechanical stimulation, these IS cells differentially express osteogenic factors, mechanosensitive genes, and signaling genes. Considering that numerous ciliary genes are associated with a scoliosis phenotype, among ciliopathies and knockout animal models, we expected IS patients to have an accumulation of rare variants in ciliary genes. Instead, our SKAT-O analysis of whole exomes showed an enrichment among IS patients for rare variants in genes with a role in cellular mechanotransduction. Our data indicates defective cilia in IS bone cells, which may be linked to heterogeneous gene variants pertaining to cellular mechanotransduction.

Accuracy of Programs for the Determination of Human Leukocyte Antigen Alleles from Next-Generation Sequencing Data.

The human leukocyte antigen (HLA) genes code for proteins that play a central role in the function of the immune system by presenting peptide antigens to T cells. As HLA genes show extremely high genetic polymorphism, HLA typing at the allele level is demanding and is based on DNA sequencing. Determination of HLA alleles is warranted as HLA alleles are major genetic risk factors in autoimmune diseases and are matched in transplantation. Here, we compared the accuracy of several published HLA-typing algorithms that are based on next-generation sequencing (NGS) data. As genome sequencing is becoming increasingly common in research, we wanted to test how well HLA alleles can be deduced from genome data produced in studies with objectives other than HLA typing and in platforms not especially designed for HLA typing. The accuracies were assessed using datasets consisting of NGS data produced using an in-house sequencing platform, including the full 4 Mbp HLA segment, from 94 stem cell transplantation patients and exome sequences from 63 samples of the 1000 Genomes collection. In the patient dataset, none of the software gave perfect results for all the samples and genes when programs were used with the default settings. However, we found that ensemble prediction of the results or modifications of the settings could be used to improve accuracy. For the exome-only data, most of the algorithms did not perform very well. The results indicate that the use of these algorithms for accurate HLA allele determination is not straightforward when based on NGS data not especially targeted to the HLA typing and their accurate use requires HLA expertise.

A replication study for association of 53 single nucleotide polymorphisms in ScoliScore test with adolescent idiopathic scoliosis in French-Canadian population.

STUDY DESIGN: A replication association study that used genomic data generated from French-Canadian case and control cohorts. OBJECTIVES: To determine whether the 53 single nucleotide polymorphisms (SNPs) that were previously associated with spinal deformity progression in an American Caucasian cohort are similarly associated in French-Canadian population. SUMMARY OF BACKGROUND DATA: It is widely accepted that genetic factors contribute to adolescent idiopathic scoliosis. The identification of genetic variants associated with the predisposition or progression of curvature could facilitate diagnostic/prognostic tool development. Although 53 SNPs have been associated with spinal curve progression in Caucasian cohorts in the United States, these associations were not replicated in a large Japanese population study, arguing that such a discrepancy could be explained by ethnicity, thus raising the importance of a replication study in an independent Caucasian population of European descent. METHODS: Genomic data were collected from the French-Canadian population, using the Illumina HumanOmni 2.5M BeadChip. Fifty-two SNPs, tested in ScoliScore or in high linkage disequilibrium with SNPs in the test, were selected to assess their association with scoliosis generally, and with spinal curve progression. One SNP in ScoliScore, rs16909285, could not be evaluated in our Genome-Wide association study. RESULTS: None of the SNPs used in ScoliScore were associated with adolescent idiopathic scoliosis curve progression or curve occurrence in French-Canadian population. We evaluated 52 SNPs in severe patients by comparing risk allele frequencies with those in nonsevere patients and with those in control individuals. There was no significant difference between the severe group and the nonsevere group or between the severe group and the control group. CONCLUSION: Although the 52 SNPs studied here were previously associated with curve progression in an American population of European descent, we found no association in French-Canadian patients with adolescent idiopathic scoliosis. This second replication cohort suggests that the lack of association of these SNPs in a Japanese cohort is not due to ethnicity. LEVEL OF EVIDENCE: 4.

Mouse ENU Mutagenesis to Understand Immunity to Infection: Methods, Selected Examples, and Perspectives.

Infectious diseases are responsible for over 25% of deaths globally, but many more individuals are exposed to deadly pathogens. The outcome of infection results from a set of diverse factors including pathogen virulence factors, the environment, and the genetic make-up of the host. The completion of the human reference genome sequence in 2004 along with technological advances have tremendously accelerated and renovated the tools to study the genetic etiology of infectious diseases in humans and its best characterized mammalian model, the mouse. Advancements in mouse genomic resources have accelerated genome-wide functional approaches, such as gene-driven and phenotype-driven mutagenesis, bringing to the fore the use of mouse models that reproduce accurately many aspects of the pathogenesis of human infectious diseases. Treatment with the mutagen N-ethyl-N-nitrosourea (ENU) has become the most popular phenotype-driven approach. Our team and others have employed mouse ENU mutagenesis to identify host genes that directly impact susceptibility to pathogens of global significance. In this review, we first describe the strategies and tools used in mouse genetics to understand immunity to infection with special emphasis on chemical mutagenesis of the mouse germ-line together with current strategies to efficiently identify functional mutations using next generation sequencing. Then, we highlight illustrative examples of genes, proteins, and cellular signatures that have been revealed by ENU screens and have been shown to be involved in susceptibility or resistance to infectious diseases caused by parasites, bacteria, and viruses.

Being Aquifex aeolicus: Untangling a hyperthermophile's checkered past.

Lateral gene transfer (LGT) is an important factor contributing to the evolution of prokaryotic genomes. The Aquificae are a hyperthermophilic bacterial group whose genes show affiliations to many other lineages, including the hyperthermophilic Thermotogae, the Proteobacteria, and the Archaea. Previous phylogenomic analyses focused on Aquifex aeolicus identified Thermotogae and Aquificae either as successive early branches or sisters in a rooted bacterial phylogeny, but many phylogenies and cellular traits have suggested a stronger affiliation with the Epsilonproteobacteria. Different scenarios for the evolution of the Aquificae yield different phylogenetic predictions. Here, we outline these scenarios and consider the fit of the available data, including three sequenced Aquificae genomes, to different sets of predictions. Evidence from phylogenetic profiles and trees suggests that the Epsilonproteobacteria have the strongest affinities with the three Aquificae analyzed. However, this pattern is shown by only a minority of encoded proteins, and the Archaea, many lineages of thermophilic bacteria, and members of genus Clostridium and class Deltaproteobacteria also show strong connections to the Aquificae. The phylogenetic affiliations of different functional subsystems showed strong biases: Most but not all genes implicated in the core translational apparatus tended to group Aquificae with Thermotogae, whereas a wide range of metabolic and cellular processes strongly supported the link between Aquificae and Epsilonproteobacteria. Depending on which sets of genes are privileged, either Thermotogae or Epsilonproteobacteria is the most plausible adjacent lineage to the Aquificae. Both scenarios require massive sharing of genes to explain the history of this enigmatic group, whose history is further complicated by specific affinities of different members of Aquificae to different partner lineages.