{'fraction_non_alphanumeric': 0.03343736067766384, 'fraction_numerical': 0.03670679149947986, 'mean_word_length': 5.0358744394618835, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 5, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Failure to detect mutations in U2AF1 due to changes in the GRCh38 reference sequence

The U2AF1 gene is a core part of mRNA splicing machinery and frequently contains somatic mutations that contribute to oncogenesis in MDS, AML, and other cancers. A change introduced in the GRCh38 version of the human reference build prevents mutations in this gene from being detected by many variant calling pipelines. We describe the problem in detail and show that a modified GRCh38 reference build with unchanged coordinates can be used to ameliorate the issue. This reference is available at https://zenodo.org/record/4684553 (doi:10.5281/zenodo.4684553)

The U2AF1 gene encodes a core component of the messenger RNA splicing machinery that is frequently mutated in MDS and other cancers [1][2][3] . Though predominantly associated with hematopoietic cancers (73%), mutations are also recurrent in lung tumors (6.5%) and have been reported in 24 other tumor types. Specifically, mutations at residues S34 and Q157 have been shown to promote exon skipping and are confirmed driver mutations contributing to cancer pathogenesis [4][5][6][7] .
As part of the National Myelodysplastic Syndrome (MDS) Natural History Study (NCT02775383), we used a targeted gene panel to sequence bone marrow samples from 120 patients either diagnosed with MDS or suspected to have MDS 8 . Of these patients, 38 were eventually confirmed to have MDS or MPN. Initial analyses looking at sequencing quality metrics revealed coverage levels and mutation frequencies that closely matched expectations, with one exception: mutations in the U2AF1 spliceosome gene are typically observed in nearly 10% of MDS or Myeloproliferative Neoplasm (MPN) patients, but only 2 of the 38 (5.2%) MDS or MPN patients in this group had such mutations, both at the Q157 hotspot 1,9 . While this deviation was not significant, it was concerning that both mutations had far lower sequence coverage than would be expected, as did the rest of the U2AF1 gene. Manual inspection with IGV revealed that this poor coverage extended to a 150 kilobase region of the genome where the majority of reads had a mapping quality of zero ( Fig. 1) 10 . This region spans the CBS, U2AF1, FRGCA, and CRYAA genes.
Further investigation revealed that in the GRCh38 reference build, content added to the p-arm of chromosome 21 (chr21:6427259-6580181) contained sequence that nearly exactly replicated the sequence of the U2AF1 locus (chr21:43035875-43187577). The same issue does not exist in prior reference builds GRCh36 or GRCh37. After consultation with members of the Genome Reference Consortium (GRC), it was determined that a BAC clone (https://www.ncbi.nlm.nih.gov/nuccore/FP236240.8) was incorrectly added to the reference genome, creating this duplicate sequence. This resulted in the alignment algorithm, BWA-MEM, splitting the reads among these two loci and assigning mapping quality scores of zero to reads that were multi-mapped 11 . This lowered the overall coverage substantially (Fig. 1,   Supplementary Fig. 1), and as such alignments are typically excluded or down-weighted during variant calling, explained the paucity of mutations we observed in U2AF1, especially at the S34F position. To address this issue, we created a modified version of GRCh38 that maintains the coordinate system, but replaces the new duplicate sequence on chromosome 21p with "N" characters. We realigned the data to this reference, and observed a substantial increase in coverage and mapping qualities across the affected region. Over exons of U2AF1, the coverage rose from a median of 0.3x to a median of 1,195x. This enabled the discovery of an additional U2AF1 mutation (S34F) in this cohort (Supplementary Table 1).
To validate this finding in an orthogonal data set, we examined data from Acute Myeloid Leukemia (AML) patients sequenced for The Cancer Genome Atlas paper 12 . The data in this study were originally aligned to genome build GRCh36, and 6 mutations in the U2AF1 gene were reported from exome sequencing. We downloaded exome sequence data that had been realigned (using BWA-MEM) by the core Genomic Data Commons pipelines against the GRCh38.d1.vd1 reference 11,13 . The GDC also provided somatic variant files from four different algorithms, none of which reported mutations in the U2AF1 gene in any sample. Running our own somatic variant calling pipeline (https://git.io/JYlpq) on these alignments also did not reveal any U2AF1 mutations.
Review of the alignments from these TCGA cases revealed the same coverage and mapping quality issues on chromosome 21. One sample lacked any U2AF1 reads in the available exome data due to unrelated sequencing problems (the TCGA consortium only identified the mutation using orthogonal assays). After realigning the data from the remaining 5 cases to our modified version of GRCh38, the same somatic pipeline as above identified all 5 evaluable mutations in U2AF1 (Supplemental Table 1).
The reference genome is essential to modern cancer genomics, but problems with these assemblies have the potential to cause both false positives and negatives. In this study, we describe the latter, where changes introduced in the GRCh38 human reference build cause mutations in a key driver gene to be missed with standard analysis approaches. GRCh38 was first released in late 2013, but widespread adoption lagged somewhat, so some studies involving myeloid malignancies (TCGA, Beat AML) have been spared this issue because they used older versions of the reference. Nonetheless, we expect that large cancer databases, including the NCI's Genomic Data Commons, may be missing many U2AF1 mutations due to this artifact. This also has clinical implications, because U2AF1 mutations have strong associations with prognosis and clinical trials of splice-modulating drugs are being planned or are underway 5,14,15 .
We have reported these findings to the Genome Reference Consortium (https://www.ncbi.nlm.nih.gov/grc/human/issues/HG-2544) but as there are no patches or new releases currently scheduled for the human reference genome, the problem remains unresolved in the current release of the GRCh38 reference (GRCh38.p13). In order to ameliorate the problem, we have generated a corrected reference that is coordinate-compatible, and the differences in the resulting alignments are almost entirely in the masked and U2AF1 regions.
This means it can be used interchangeably with the standard GRCh38 reference in applications where detection of U2AF1 mutations is critical, including sequencing of hematological cancers or studies of spliceosome dysfunction. The reference fasta file is available at https://zenodo.org/record/4684553 (doi:10.5281/zenodo.4684553).