2015-08-21

Variant Normalization in hgvs: Application

Variant normalization was recently implemented in the hgvs Python package [see blog post]. A major motivation for variant normalization is to facilitate text-based comparisons of variant observations using equivalent canonical representations. From a clinical standpoint, this means we would be better able to compare published variant interpretations.

I apply the normalizer in hgvs package to normalize the variants in Clinvitae (http://clinvitae.invitae.com/) and analyze the results.


The Clinvitae collects and stores clinically observed genetic variants from multiple public databases. All the variants can be downloaded in one file.

The Clinvitae has 180,974 variant in total, including both genome level and transcripts level variants. In this analysis, I will focus on transcripts level variants (which type is ‘c.’ and the accession starts with ‘NM’). There are a total of 141,338 raw transcripts level variants in Clinvitae. However, 1,162 (0.82%) could not be parsed correctly and are dropped in the following analysis. 23,119 (16.36%) variants are intronic variants and is not supported to be normalized, because the reference sequence of introns is not defined for RefSeq transcripts. 342 variants are defined on sequences unknown to UTA or that do not have genomic alignments in UTA. 1,249 (0.88%) variants are invalid, including variants that base start position greater than end position etc. 183 variants have wrong reference allele. 179 variants are identity variants. For the remaining set of 115,104 (81.44%) variants that could be properly parsed and validated, duplicated variants were removed, reducing the set to 93,761 distinct variants. These distinct variants comprised the test set and was normalized using the hgvs.normalizer module
.

By default, the hgvs normalizer right (3’) shuffles variants and does not permit shuffling across exon boundaries. However, these options are selectable at runtime. To demonstrate this functionality and highly possible errors in variant reporting, we used  four configurations during normalization:
  1. shuffle to 3’ and allow crossing exon-intron boundary
  2. shuffle to 5’ and allow crossing exon-intron boundary
  3. shuffle to 3’ and not allow crossing exon-intron boundary
  4. shuffle to 5’ and not allow crossing exon-intron boundary

shuffle.png

Among these 93,761 variants, 31 variants span the exon-intron boundary. Among these 31 variants, 22 could be normalized when allowing crossing exon-intron boundary. For the remaining 93,730 variants, 86,114 (91.91%) variants remained the same no matter using whichever configuration of normalization, suggesting that these variants are less likely to have other equivalent forms reported 7,586 (8.09%) variants are normalized by at least one configuration.


For the variants that could be normalized, 1,011 could be normalized by 3’ shuffling and 7,020 could be normalized by 5’ shuffling. This indicates that most variants in Clinvitae have already located at 3’ most, consistent with the HGVS recommendations. 445 variants could be shuffled to either 3’ or 5’ direction.

27 variants are located at 3’ end of exons, so they could be 3’ shuffled only when allowing crossing exon-intron boundary and remained the same when that is not allowed. Two variants located near the 3’ end of exons normalized to different results depending on whether exon boundary crossing was permitted. Similarly, 62 variants locate at 5’ end of exons could be 5’ shuffled when allowing boundary crossing and remained the same when that is not allowed. 40 variants located near 5’ end of an exon and the normalization results are different for crossing exon-intron boundary is allowed or not. Other variants normalization results are not affected by whether allowing crossing exon-intron boundary.

The hgvs normalizer also rewrites variants according to HGVS recommendations depending on sequence context. In the test set, 169 insertions were converted to duplications and 22 delins were converted to inversions.

In summary, some variants in Clinvitae are not normalized as 3’ most as possible according to the HGVS recommendations. And some variants are not correctly described. By using a standard, freely-available normalization process, we will be able to more reliably correlate variant observations with clinical significance.





Implementing Sequence Variant Normalization in the hgvs package


In principle, sequence variants have many equivalent representations. For instance, delGinsA could also be written as delTGinsTA.  In addition, variants may have more than one description depending on the sequence context, such as variants occurring in repeat regions. For example, AGTTTC to AGTTC could be described as 3delT or 4delT or 5delT.


The most basic notion of variant comparison is equivalence of a normalized representation. Variant normalization is essential for comparing and counting variants. Several tools have been developed to perform normalization of variants in VCF files, like vt (https://github.com/atks/vt). Although the HGVS recommendations specify rules for variant normalization, such as using the 3’-most representation of a variant, the hgvs Python package did not support them. Unlike variants in VCF format which is only represented by the reference allele and alternative allele, the variants in HGVS format contain much more information and different types of variants are explicitly expressed in different ways like substitution(>), delins, del, ins, dup, dupN, inv and con etc. So the normalization of HGVS variants is more complicated than the normalization of VCF variants.


As part of my Google Summer of Code project, I implemented normalization of variants in the hgvs package. The normalization procedure contains two parts:
  1. trim the common prefix and suffix of reference allele and alternative allele;
  2. shuffle the variants to as 3’ most (or 5’ most) as possible.


In hgvs normalizer, we utilize a dynamic extension local window to perform the variants normalization. By default, the window size is 3 times length of the maximum length of reference allele and alternative allele. And the reference sequence and alternative sequence are reconstructed based on the UTA database and the variant itself. The trimming and shuffling process are based on the code of vgraph (https://github.com/bioinformed/vgraph). If the shuffling stops before reaching the cutting edge of the window boundary, or reaching the end of sequence, the normalization is finished. When the shuffling reaches the edge of the window, the normalizer will extend the window and continue the shuffling process, until the shuffling stops before reaching to the boundary of the window.







Here is a real example illustrating how one variant is normalized in hgvs:




Here is the summary of the normalizer in hgvs:

  • The normalizer is configurable in shuffling direction (3’ most or 5’ most) and whether allowing shuffling crossing the exon-intron boundary;
HGVS recommendations are to shuffle variants to the 3'-most position on the reference sequence. But for variants in VCF files, they are usually required to be shuffled to the 5’-most position. And currently there is no guidance about whether the shuffling should cross the exon-intron boundary or not in HGVS recommendations.


  • The priority of variant output type is dupN > dup > ins;
Fox example, when the reference is AGTTC and alternative sequence is AGTTTTC, the normalized result is c.4dup2, rather than c.3_4dupT or c.5_6insTT.


  • The normalizer in hgvs supports the normalization of all kinds of variants, except conversions and protein variants (p.).
When the crossing exon-intron boundary is disabled, variants that cross the exon-intron boundary could not be normalized, since that would cause confusing.



Variant normalization will appear in hgvs 0.4.0, due for release in August 2015.