Answers.Body: <ol>
  <li><p>Don't know. These could be retired or inconsistent SNPs (e.g. shifted locations, or otherwise found to be different from expected). It might be that the imputation suggests that the SNP should be something other than what is observed. I can't see anything obvious from the <a href="https://www.ncbi.nlm.nih.gov/snp/" rel="nofollow noreferrer">dbSNP</a> information that I've looked at.</p>
  </li>
  <li><p>a) In this specific case for rs116059446, that's a SNP that has more than two variants. One way of normalising VCF files is to rearrange them so that each line only refer to a single polymorphism from the reference sequence, so if there are two non-reference variants, the variant will get two lines in the VCF file.</p>
  <p>b) For <a href="https://www.ncbi.nlm.nih.gov/snp/rs331869" rel="nofollow noreferrer">rs331869</a>, this is another SNP with more than two variants recorded in dbSNP. Given that you're digging through the weeds of a huge dataset, the way of thinking about this needs to change from &quot;Is this common&quot; to &quot;Is this possible&quot;. Problems are much more likely to exist with SNPs that aren't biallelic, because it indicates that something funky is going on with genetics or testing. I'll propose three possibilities for what could be happening here, in what I consider to be descending order of likelihood:</p>
  <ol>
  <li>you have one non-reference allele on one of the sister chromosome copies, and another non-reference allele on another copy</li>
  <li>you have a somatic mutation (i.e. a mutation that happened during your lifetime), and some of your cells have one variant, and others have a different variant.</li>
  <li>the probes are tagging different regions in your genome from what is expected. Probes should be location-unique, but that's impossible to guarantee for everyone, because genetic variation is a thing that happens.</li>
  </ol>
  </li>
  <li><p>In this specific case for <a href="https://www.ncbi.nlm.nih.gov/snp/rs71634747" rel="nofollow noreferrer">rs71634747</a>, it looks like that SNP has been retired, and has merged with <a href="https://www.ncbi.nlm.nih.gov/snp/rs724135" rel="nofollow noreferrer">rs724135</a>. It's possible that this is a multi-nucleotide variant where variants are only observed in tandem. This is an intronic variant, so seeing a 3-nucleotide variation is a little unexpected, but given that you're picking out exceptions from a very large dataset, not completely unexpected.</p>
  </li>
  </ol>
  <p>I think the main take-home from this is that you shouldn't ever expect results to be a perfect representation of reality. It can be useful to look at raw data to help understand what's going on, but there's too much weird stuff going on with biology, chemistry and physics for anything to be 100% accurate.</p>
  <p>If you did exactly the same test at a later date, some of your 23andMe results will be different. If you did imputation at a later date, with a different program version, or different imputation reference dataset, your results are likely to be different.</p>
  <p>Accuracy can be substantially improved by using different technologies to attempt to discover the same thing, and only concentrating on consistently-reported differences. You could try an Infinium Omni2.5 SNPchip, for example, or whole-genome sequencing. Unfortunately, this will lead to another problem: what should be done with the inconsistent results? The answer to <em>that</em> is situation dependent, and is one of the hardest challenges of bioinformatics.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-03-10T13:54:03.100Z
Answers.Text: Thanks for these answers. I'm not expecting data to reflect reality (I've been doing bioinformatics for too long for that ;-) I just expect data to be sane. This is just a basic sanity test, flagging things that looked wrong. I think point 1 and 3 are now tidied up (I don't think the software is behaving correctly, but I can at least understand how the software is behaving, and it's not bonkers.
Answers.UserDisplayName: Dan Bolser
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-03-10T13:57:14.477Z
Answers.Text: I still don't understand point 2. Lets look at your potential sources for this data sanity check failure... 2b.1) makes perfect sense, but I would expect the genotype to be 0/1 and 0/1. That would be easy to interpret. That's not what I see (1/0 and 1/1). 2b.2) This may be true, but the data could never show that (imputation). 2b.3) Same, not possible in this data (imputed). Increasingly I think this call is just wrong. Not biologically wrong, just fails basic sanity checks.
Answers.UserDisplayName: Dan Bolser
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-03-10T14:02:09.223Z
Answers.Text: Like imagine I showed you a vcf that had `3 3333 . A G . . . GT 4/4`... Well, you could say, ah ha, this is a sign of a quadruplex DNA at this location in this individual... well maybe. I'd just say the VCF is broken. It's just wrong / impossible to interpret (imho).
Answers.UserDisplayName: Dan Bolser
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-03-08T00:01:14.717Z
Answers.LastActivityDate: 2022-03-10T03:32:09.807Z
Answers.LastEditDate: 2022-03-10T03:32:09.807Z
Answers.UserDisplayName: gringer
Answers.Body: <p>The answer is to relax and throw away the 0.001% of stuff that doesn't make sense, lower your standards and just get on with finishing the job you set out to do without actually caring about consistency, logic, meaning or formats that are actually defined or interpretable.</p>
  <p>Nobody cares about what these SNPs mean, how to interpret them programmatically, why there are version issues in the output of the imputation server or why it throws away ground truth in the output, leaving the user to struggle with all these issues. I mean... if nobody cares, why should I?</p>
  <p>I've written up as detailed a description of the problems as I can here:
  <a href="https://github.com/richarddurbin/pbwt/issues/52" rel="nofollow noreferrer">https://github.com/richarddurbin/pbwt/issues/52</a> and here: <a href="https://github.com/richarddurbin/pbwt/issues/51" rel="nofollow noreferrer">https://github.com/richarddurbin/pbwt/issues/51</a></p>
  <p>and written a script to throw away all the junk here:
  <a href="https://github.com/Geromics/covcheck/blob/wip/report-v2/Research/debug_imputation_results.py" rel="nofollow noreferrer">https://github.com/Geromics/covcheck/blob/wip/report-v2/Research/debug_imputation_results.py</a></p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2022-06-10T10:56:11.587Z
Answers.LastActivityDate: 2022-06-10T10:56:11.587Z
Answers.UserDisplayName: Dan Bolser
Body: <p>I ran my 23andMe data through the <a href="https://imputation.sanger.ac.uk/" rel="nofollow noreferrer">Sanger Imputation Service</a> that uses <a href="https://alkesgroup.broadinstitute.org/Eagle/" rel="nofollow noreferrer">Eagle</a> v2.4 for phasing, <a href="https://github.com/VertebrateResequencing/pbwt" rel="nofollow noreferrer">PBWT</a> v3.1 for imputation.</p>
  <p>However, some aspects of the results are very confusing to me. Perhaps see #3 (the most confusing) first.</p>
  <p><strong>Point 1</strong></p>
  <ol>
  <li>There are 62,623 rsIDs for chromosome 3 in my 23andMe file and 2,821,894 SNPs (actually, unique rsIDs) after imputation. However, only 60,144 of the original 62,623 rsIDs are in that ~2.8 million. Where did the other 2,479 go? Here are some examples of missing rsIDs:
  <a href="https://www.ncbi.nlm.nih.gov/snp/rs1516332" rel="nofollow noreferrer">rs1516332</a>,
  <a href="https://www.ncbi.nlm.nih.gov/snp/rs12374025" rel="nofollow noreferrer">rs12374025</a>,
  <a href="https://www.ncbi.nlm.nih.gov/snp/rs1170695" rel="nofollow noreferrer">rs1170695</a>,
  <a href="https://www.ncbi.nlm.nih.gov/snp/rs2731343" rel="nofollow noreferrer">rs2731343</a> and
  <a href="https://www.ncbi.nlm.nih.gov/snp/rs63749864" rel="nofollow noreferrer">rs63749864</a>.</li>
  </ol>
  <p><strong>EDIT Point 1</strong></p>
  <p>Looking at <code>rs1516332</code>, the position is 3:105793, which does have a result in the file, but it doesn't have an rsID:</p>
  <pre><code>3    105793    .    C    T    .    PASS    RefPanelAF=0.17878;AN=2;AC=0;INFO=1    GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  </code></pre>
  <p>Looking at <code>rs63749864</code> (3:13564048) it's missing by rsID and position.</p>
  <p><strong>EDIT Point 1 (again)</strong></p>
  <p>As described here: <a href="https://github.com/richarddurbin/pbwt/issues/51" rel="nofollow noreferrer">https://github.com/richarddurbin/pbwt/issues/51</a></p>
  <blockquote>
  <p>&quot;... of the 2,479 'missing' SNPs, 1,806 can be 'found' [by position] leaving 673 'anchors' [or rsIDs called by 23andMe] missing from chromosome 3.&quot;</p>
  </blockquote>
  <p><strong>Point 2</strong></p>
  <ol start="2">
  <li>Some imputed rsIDs appear more than once in the results. How should I interpret that? e.g.</li>
  </ol>
  <pre><code>3   64071   rs116059446 G   T   .   PASS    RefPanelAF=0.000369572;AN=2;AC=0;INFO=1 GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  3   64071   rs116059446 G   A   .   PASS    RefPanelAF=9.2393e-05;AN=2;AC=0;INFO=1  GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  </code></pre>
  <p><strong>EDIT Point 2</strong>
  From the answer below, this is a variant with more than two alleles: <a href="https://www.ncbi.nlm.nih.gov/snp/rs116059446" rel="nofollow noreferrer">rs116059446</a></p>
  <p>Seems the format of this VCF file only allows one alt per line, hence the two alts are split over two lines.</p>
  <p>See my question here: <a href="https://bioinformatics.stackexchange.com/questions/18752/is-it-valid-vcf-not-to-squash-positions-with-more-than-one-alt-allele">Is it valid VCF not to &#39;squash&#39; positions with more than one ALT allele?</a></p>
  <p>And the answer from Durbo here: <a href="https://github.com/richarddurbin/pbwt/issues/52#issuecomment-1062951105" rel="nofollow noreferrer">https://github.com/richarddurbin/pbwt/issues/52#issuecomment-1062951105</a></p>
  <p>In the above case, the genotype is easy to interpret (REF/REF in both lines).</p>
  <p>However, I'm now finding examples that I just can't interpret:</p>
  <pre><code>3   255395  rs331869    C   G   .   PASS    RefPanelAF=0.385864;AN=2;AC=1;INFO=1    GT:ADS:DS:GP    1|0:1,0:1:0,1,0
  3   255395  rs331869    C   T   .   PASS    RefPanelAF=0.980998;AN=2;AC=2;INFO=1    GT:ADS:DS:GP    1|1:1,1:2:0,0,1
  </code></pre>
  <p><em>Is my genotype CG or TT at this position?</em></p>
  <p><strong>Point 3</strong></p>
  <ol start="3">
  <li>Sometimes a rsID appears twice but at different locations:</li>
  </ol>
  <pre><code>3   4942430 rs71634747  G   A   .   PASS    RefPanelAF=0.445673;AN=2;AC=1;INFO=1    GT:ADS:DS:GP    1|0:1,0:1:0,1,0
  3   4942432 rs71634747  C   G   .   PASS    RefPanelAF=0.248907;AN=2;AC=0;INFO=1    GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  </code></pre>
  <p>Unless I'm really confused, this just seems wrong. If it wasn't for this weirdness, I'd seek concrete answers for #1 and #2, but it could be that they are both just a symptom of this potential error?</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2022-03-07T14:22:00.000Z
Comments.Text: Oh, vr-impute@sanger.ac.uk isn't responding to my email. I have asked the above of them.
Comments.UserDisplayName: Dan Bolser
Tags: <vcf><snp><imputation><pbwt>
Title: Confusing result from Sanger Imputation Service (Eagle v2.4 for phasing, PBWT v3.1 for imputation)
ThreadText: # Confusing result from Sanger Imputation Service (Eagle v2.4 for phasing, PBWT v3.1 for imputation)
  
  (Asked by: username_0 on 2022-03-07)
  
  I ran my 23andMe data through the [Sanger Imputation Service](https://imputation.sanger.ac.uk/) that uses [Eagle](https://alkesgroup.broadinstitute.org/Eagle/) v2.4 for phasing, [PBWT](https://github.com/VertebrateResequencing/pbwt) v3.1 for imputation.
  
  However, some aspects of the results are very confusing to me. Perhaps see #3 (the most confusing) first.
  
  **Point 1**
  
  1. There are 62,623 rsIDs for chromosome 3 in my 23andMe file and 2,821,894 SNPs (actually, unique rsIDs) after imputation. However, only 60,144 of the original 62,623 rsIDs are in that ~2.8 million. Where did the other 2,479 go? Here are some examples of missing rsIDs:
     [rs1516332](https://www.ncbi.nlm.nih.gov/snp/rs1516332),
     [rs12374025](https://www.ncbi.nlm.nih.gov/snp/rs12374025),
     [rs1170695](https://www.ncbi.nlm.nih.gov/snp/rs1170695),
     [rs2731343](https://www.ncbi.nlm.nih.gov/snp/rs2731343) and
     [rs63749864](https://www.ncbi.nlm.nih.gov/snp/rs63749864).
  
  **EDIT Point 1**
  
  Looking at `rs1516332`, the position is 3:105793, which does have a result in the file, but it doesn't have an rsID:
  
  ```
  3    105793    .    C    T    .    PASS    RefPanelAF=0.17878;AN=2;AC=0;INFO=1    GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  
  ```
  
  Looking at `rs63749864` (3:13564048) it's missing by rsID and position.
  
  **EDIT Point 1 (again)**
  
  As described here: <https://github.com/richarddurbin/pbwt/issues/51>
  
  > "... of the 2,479 'missing' SNPs, 1,806 can be 'found' [by position] leaving 673 'anchors' [or rsIDs called by 23andMe] missing from chromosome 3."
  
  **Point 2**
  
  2. Some imputed rsIDs appear more than once in the results. How should I interpret that? e.g.
  
  ```
  3   64071   rs116059446 G   T   .   PASS    RefPanelAF=0.000369572;AN=2;AC=0;INFO=1 GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  3   64071   rs116059446 G   A   .   PASS    RefPanelAF=9.2393e-05;AN=2;AC=0;INFO=1  GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  
  ```
  
  **EDIT Point 2**
  From the answer below, this is a variant with more than two alleles: [rs116059446](https://www.ncbi.nlm.nih.gov/snp/rs116059446)
  
  Seems the format of this VCF file only allows one alt per line, hence the two alts are split over two lines.
  
  See my question here: [Is it valid VCF not to 'squash' positions with more than one ALT allele?](https://bioinformatics.stackexchange.com/questions/18752/is-it-valid-vcf-not-to-squash-positions-with-more-than-one-alt-allele)
  
  And the answer from Durbo here: <https://github.com/richarddurbin/pbwt/issues/52#issuecomment-1062951105>
  
  In the above case, the genotype is easy to interpret (REF/REF in both lines).
  
  However, I'm now finding examples that I just can't interpret:
  
  ```
  3   255395  rs331869    C   G   .   PASS    RefPanelAF=0.385864;AN=2;AC=1;INFO=1    GT:ADS:DS:GP    1|0:1,0:1:0,1,0
  3   255395  rs331869    C   T   .   PASS    RefPanelAF=0.980998;AN=2;AC=2;INFO=1    GT:ADS:DS:GP    1|1:1,1:2:0,0,1
  
  ```
  
  *Is my genotype CG or TT at this position?*
  
  **Point 3**
  
  3. Sometimes a rsID appears twice but at different locations:
  
  ```
  3   4942430 rs71634747  G   A   .   PASS    RefPanelAF=0.445673;AN=2;AC=1;INFO=1    GT:ADS:DS:GP    1|0:1,0:1:0,1,0
  3   4942432 rs71634747  C   G   .   PASS    RefPanelAF=0.248907;AN=2;AC=0;INFO=1    GT:ADS:DS:GP    0|0:0,0:0:1,0,0
  
  ```
  
  Unless I'm really confused, this just seems wrong. If it wasn't for this weirdness, I'd seek concrete answers for #1 and #2, but it could be that they are both just a symptom of this potential error?
  
  **Tags:** vcf, snp, imputation, pbwt
  
  > **Comments:**
  > **username_0**: Oh, vr-impute@sanger.ac.uk isn't responding to my email. I have asked the above of them.  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_1 on 2022-03-08. Score: 3
  
  1. Don't know. These could be retired or inconsistent SNPs (e.g. shifted locations, or otherwise found to be different from expected). It might be that the imputation suggests that the SNP should be something other than what is observed. I can't see anything obvious from the [dbSNP](https://www.ncbi.nlm.nih.gov/snp/) information that I've looked at.
  2. a) In this specific case for rs116059446, that's a SNP that has more than two variants. One way of normalising VCF files is to rearrange them so that each line only refer to a single polymorphism from the reference sequence, so if there are two non-reference variants, the variant will get two lines in the VCF file.
  
     b) For [rs331869](https://www.ncbi.nlm.nih.gov/snp/rs331869), this is another SNP with more than two variants recorded in dbSNP. Given that you're digging through the weeds of a huge dataset, the way of thinking about this needs to change from "Is this common" to "Is this possible". Problems are much more likely to exist with SNPs that aren't biallelic, because it indicates that something funky is going on with genetics or testing. I'll propose three possibilities for what could be happening here, in what I consider to be descending order of likelihood:
  
     1. you have one non-reference allele on one of the sister chromosome copies, and another non-reference allele on another copy
     2. you have a somatic mutation (i.e. a mutation that happened during your lifetime), and some of your cells have one variant, and others have a different variant.
     3. the probes are tagging different regions in your genome from what is expected. Probes should be location-unique, but that's impossible to guarantee for everyone, because genetic variation is a thing that happens.
  3. In this specific case for [rs71634747](https://www.ncbi.nlm.nih.gov/snp/rs71634747), it looks like that SNP has been retired, and has merged with [rs724135](https://www.ncbi.nlm.nih.gov/snp/rs724135). It's possible that this is a multi-nucleotide variant where variants are only observed in tandem. This is an intronic variant, so seeing a 3-nucleotide variation is a little unexpected, but given that you're picking out exceptions from a very large dataset, not completely unexpected.
  
  I think the main take-home from this is that you shouldn't ever expect results to be a perfect representation of reality. It can be useful to look at raw data to help understand what's going on, but there's too much weird stuff going on with biology, chemistry and physics for anything to be 100% accurate.
  
  If you did exactly the same test at a later date, some of your 23andMe results will be different. If you did imputation at a later date, with a different program version, or different imputation reference dataset, your results are likely to be different.
  
  Accuracy can be substantially improved by using different technologies to attempt to discover the same thing, and only concentrating on consistently-reported differences. You could try an Infinium Omni2.5 SNPchip, for example, or whole-genome sequencing. Unfortunately, this will lead to another problem: what should be done with the inconsistent results? The answer to *that* is situation dependent, and is one of the hardest challenges of bioinformatics.
  
  > **Comments:**
  > **username_0**: Thanks for these answers. I'm not expecting data to reflect reality (I've been doing bioinformatics for too long for that ;-) I just expect data to be sane. This is just a basic sanity test, flagging things that looked wrong. I think point 1 and 3 are now tidied up (I don't think the software is behaving correctly, but I can at least understand how the software is behaving, and it's not bonkers.  
  > **username_0**: I still don't understand point 2. Lets look at your potential sources for this data sanity check failure... 2b.1) makes perfect sense, but I would expect the genotype to be 0/1 and 0/1. That would be easy to interpret. That's not what I see (1/0 and 1/1). 2b.2) This may be true, but the data could never show that (imputation). 2b.3) Same, not possible in this data (imputed). Increasingly I think this call is just wrong. Not biologically wrong, just fails basic sanity checks.  
  > **username_0**: Like imagine I showed you a vcf that had `3 3333 . A G . . . GT 4/4`... Well, you could say, ah ha, this is a sign of a quadruplex DNA at this location in this individual... well maybe. I'd just say the VCF is broken. It's just wrong / impossible to interpret (imho).  
  
  ---
  
  ### Answer by username_0 on 2022-06-10. Score: 0
  
  The answer is to relax and throw away the 0.001% of stuff that doesn't make sense, lower your standards and just get on with finishing the job you set out to do without actually caring about consistency, logic, meaning or formats that are actually defined or interpretable.
  
  Nobody cares about what these SNPs mean, how to interpret them programmatically, why there are version issues in the output of the imputation server or why it throws away ground truth in the output, leaving the user to struggle with all these issues. I mean... if nobody cares, why should I?
  
  I've written up as detailed a description of the problems as I can here:
  <https://github.com/richarddurbin/pbwt/issues/52> and here: <https://github.com/richarddurbin/pbwt/issues/51>
  
  and written a script to throw away all the junk here:
  <https://github.com/Geromics/covcheck/blob/wip/report-v2/Research/debug_imputation_results.py>
  
  
  ---
  
  