Error in the `ref` base does not trigger a validation error for an intronic variant

[ielis]

Describe the bug

Error in the reference is not picked up for an intronic variant. This happens for the following query:

• genome_build: hg38
• variant_description: NM_006846.4:c.2666+5A>G
• select_transcripts: mane_select

The 5th intronic base in fact corresponds to G (not A, as indicated in the variant description), but this is not picked up by variant validator. The gene in question (SPINK5) is on the forward strand.

The returned JSON document includes the following:

"validation_warnings": [
      "NM_006846.4:c.2666+5A>G automapped to NM_006846.4:c.2666+5="
    ]

However, I am not sure this is very useful.

To Reproduce

Execute the following query:

curl -X GET https://rest.variantvalidator.org/VariantValidator/variantvalidator/hg38/NM_006846.4%3Ac.2666%2B5A%3EG/mane_select?content-type=application%2Fjson -H 'accept: application/json'

Expected behavior

I believe, that the variant validator should complain about the error in the reference.

Screenshots
Not applicable

Desktop (please complete the following information):
Not applicable

Smartphone (please complete the following information):
Not applicable

Additional context
Not applicable

[Peter-J-Freeman]

@ielis. Thanks for bringing this up.

You are correct. VV has seen the issue, hence the NM_006846.4:c.2666+5 but actually it does make more sense to return an error, just as we would had the issue been exonic.

@leicray @John-F-Wagstaff agreed?

@ifokkema, would this change be OK for you, i.e. returning an error rather than handling

[Peter-J-Freeman]

Due to the complexities of intronic variants, we cannot hard fail these issues like we would for exonic variants but we can warn

        "validation_warnings": [
            "NM_006846.4:c.2666+5A>G mapped to NC_000005.10:g.148123965A>G: Variant reference (A) does not agree with reference sequence (G)",
            "NM_006846.4:c.2666+5A>G automapped to NM_006846.4:c.2666+5="
        ],

This might well be the best we can do otherwise we will miss

"validation_warnings": [
            "NM_000088.3:c.642+1GG>G automapped to NM_000088.3:c.642+1_642+2delinsG",
            "NM_000088.3:c.642+1_642+2delGGinsG mapped to NC_000017.11:g.50197947_50197948delCCinsC: Variant reference (CC) does not agree with reference sequence (AC)",
            "TranscriptVersionWarning: A more recent version of the selected reference sequence NM_000088.3 is available for genome build GRCh38 (NM_000088.4)"
        ]

and

"validation_warnings": [
            "OutOfBoundsError: start or end or both are beyond the bounds of transcript record",
            "NM_001128425.1:c.35+1C>T mapped to NC_000001.10:g.45805891G>A: Variant reference (G) does not agree with reference sequence (C)",
            "ExonBoundaryError: Position c.35+1 does not correspond with an exon boundary for transcript NM_001128425.1"
        ]

[John-F-Wagstaff]

Yes, sorry for not replying earlier. The warning is probably better I think.

This is particularly an issue for transcripts with differing intron sequence between the two reference genomes, or main chromosome VS patch/alt sequence mappings. @Peter-J-Freeman what warnings are you intending to make turn up if we get this kind of issue in the final mapping stage? Some kind of "you should use the NM_(NC_)..." warning with the output ref that matches the input variation would be nice, but looks unpleasant to code in correctly.

These kind of differences can be either artefactual, or real. We have seen inputs which where probably caused by spurious variation detection by the user, but either way users in these cases probably need to see more than just a failure message to know what is going on when they submit a variant that is valid only in one genome (i.e. receive a full validation response with a warning and all genomic mappings, to compare, not just an error).

[Peter-J-Freeman]

The warning should be generated based on the genomic reference sequence under consideration. However, I am unsure whether it would currently be triggered by anything but the primary assembly mapping. I will push what I have to a branch and we may need some additional examples for NT/W_(NM_)

As it stands, I think it will only warn for the NC_ auto-selected by evm for the givn genome build

[Peter-J-Freeman]

I stand corrected

import json
import VariantValidator
vval = VariantValidator.Validator()
variant = "NW_011332691.1(NM_012234.6):c.335+1A>G"
genome_build = 'GRCh38'
select_transcripts = 'all'
transcript_set = 'refseq'
validate = vval.validate(variant, genome_build, select_transcripts, transcript_set)
validation = validate.format_as_dict(with_meta=True)
print(json.dumps(validation, sort_keys=True, indent=4, separators=(',', ': ')))

        "validation_warnings": [
            "NW_011332691.1:g.53828T>C: Variant reference (T) does not agree with reference sequence (C)",
            "TranscriptVersionWarning: A more recent version of the selected reference sequence NM_012234.6 is available for genome build GRCh38 (NM_012234.7)",
            "NM_012234.6:c.335+1A>G automapped to NM_012234.6:c.335+1=",
            "NM_012234.6:c.335+1A>G mapped to NW_011332691.1:g.53828T>C: Variant reference (T) does not agree with reference sequence (C)"
        ],

OK, it ends up as a doubled warning, but I am loathe to change it much more as things will start to break. I think currently what we have is the most safe option

So, if the position is incorrect on the primary assembly, it will warn (i.e. no NW/T_ stated)
or if NT/W_ provided then it warns on what was submitted. This is all correct, because all additional mappings will be lifted and adjusted for reference. We only need to check what we are told to check and in 99% of casesa will be the primary chromosome

[ielis]

Hello everyone, thanks a lot for the background. It seems like you consider a warning as a more appropriate response for this situation.

From my side, I am interested in treating an issue like this as an error, to ensure the resources that we curate are free from errors in reference. So, I will have to update our code to check for errors and for warnings.

In principle, this type of work could be simplified by providing a documentation listing the possible errors/warnings reported by Variant Validator API. I am not sure if such a page exists at this point - I searched for a while but I found nothing except from the Responses section of the REST API docs, which lists the 200 for Success.

In any case, thank you very much for your feedback and for the tool!

[Peter-J-Freeman]

@ielis yes, the reason that for intronic variants a warning is appropriate is because there may be multiple formatting issues. If we kill the code, we miss other problems.

VV outputs warnings and errors in the same list, so generally warning handling and filtering will be required anyway.

In principle, this type of work could be simplified by providing a documentation listing the possible errors/warnings reported by Variant Validator API. I am not sure if such a page exists at this point - I searched for a while but I found nothing except from the Responses section of the REST API docs, which lists the 200 for Success.

We totally agree. We really need to spend some time on this. Also, we want to add error codes such as ThisErrorType to all errors and warnigs to make this simpler. Its just finding the time and resource that is an issue sadly

[ifokkema]

@Peter-J-Freeman

@ifokkema, would this change be OK for you, i.e. returning an error rather than handling

Yeah, that should work. I should soon build a test suite for the VV library, so I can test these changes better, and I can have objective checks to verify if something may be broken or not.

@John-F-Wagstaff

This is particularly an issue for transcripts with differing intron sequence between the two reference genomes, or main chromosome VS patch/alt sequence mappings. @Peter-J-Freeman what warnings are you intending to make turn up if we get this kind of issue in the final mapping stage? Some kind of "you should use the NM_(NC_)..." warning with the output ref that matches the input variation would be nice, but looks unpleasant to code in correctly.

FYI: The HVNC meeting consistently runs over time, and we haven't fully discussed this yet. We still need to come up with a definitive format, but it's possible that NC(NM) will refer to Mutalyzer's interpretation and that NM(NC) will refer to VV's interpretation.

@ielis

From my side, I am interested in treating an issue like this as an error, to ensure the resources that we curate are free from errors in reference. So, I will have to update our code to check for errors and for warnings.

In principle, this type of work could be simplified by providing a documentation listing the possible errors/warnings reported by Variant Validator API. I am not sure if such a page exists at this point - I searched for a while but I found nothing except from the Responses section of the REST API docs, which lists the 200 for Success.

You're right, this can be improved. This is why I built my own VV library for my projects that works as an interface between my software (all LOVD-related projects) and VV, making sure there is consistent feedback to my projects while how VV presents the feedback can change. Are you handling things in a similar way?

[Peter-J-Freeman]

NM(NC)

is a typo, we use NC(NM)

@ifokkema @ielis Happy to have a separate issue opened to discuss error handling

[ielis]

Hello @ifokkema
I completely agree that all code for interacting with 3rd party resources, such as with VV from our end, should live at the periphery of a code base, and that is the philosophy that we practice. So, fortunately, there are not too many places that I/we need to update. :)

My point above was, however, focused on finding documentation for VV's endpoints, in particular e.g. JSON schema description of the response body. The description could help maintain the code responsible for interaction with VV's REST API. Basically, I am not sure if after accounting for the matter of this issue I am done, or if there is any other edge case waiting for me.

My strategy is to update the code empirically (the existing logic + the logic for checking the warnings for errors in reference in introns) and wait for further developments. A separate issue for error handling is probably not necessary.

[Peter-J-Freeman]

My point above was, however, focused on finding documentation for VV's endpoints, in particular e.g. JSON schema description of the response body. The description could help maintain the code responsible for interaction with VV's REST API. Basically, I am not sure if after accounting for the matter of this issue I am done, or if there is any other edge case waiting for me.

This is a very fair point. Gonna try get on with it ASAP, get the error codes in place then write the manual

[ifokkema]

Sorry for the late reply; I'm going through older emails now...

@Peter-J-Freeman

NM(NC)

is a typo, we use NC(NM)

No, not a typo; the point is that we have two interpretations of the same syntax: NC(NM). So we need to pull them apart or declare one invalid. Mutalyzer's NC(NM) corresponds to the interpretation of NC(GENE), so there is some support for using Mutalyzer's interpretation for NC(NM). However, since VV's interpretation is a requirement for handling intronic variants properly, we'll still need something for that, and NM(NC) was suggested. However, the committee has still not reached a conclusion. Hopefully, Monday's meeting will create some progress here.

@ielis

My point above was, however, focused on finding documentation for VV's endpoints, in particular e.g. JSON schema description of the response body. The description could help maintain the code responsible for interaction with VV's REST API. Basically, I am not sure if after accounting for the matter of this issue I am done, or if there is any other edge case waiting for me.

Yeah, that makes sense completely. I have the same problem.

[Peter-J-Freeman]

@ifokkema As I understand it, the HGVS convention for describing intronic bases correctly is NC(NM) and is plastered all over the HGVS web pages. Johan also talks about this requirement. So, the Mutalyzer Mutalyzer's NC(NM) corresponds to the interpretation of NC(GENE) is incorrect according to HGVS rules. We need to stop changing the nomenclature to accommodate Mutalyzer!!! If we have to change to NM(NC) then the HGVS nomenclature needs to change because in my opinion, we created the code to adhere to the standard. Am I missing somethign here? @leicray what do you think? To me, mutalyzer shoud be using NM(NC) for Mutalyzer's NC(NM) corresponds to the interpretation of NC(GENE) as it is not the HGVS interpretation. Also, what is meant by NC(GENE)? Gene is NG_ or ENSG. Anything to the left of the colon in :c. is reserved for Reference Sequences.

Software should meet the needs of the standard, the standard should not be changed to meet the needs of a software!