2
$\begingroup$

I recently used Dante to obtain a WGS and have downloaded several files produced during this process. There is a file that ends '.vcf' which I did think was a VCF file but I'm not sure now after a previous question. Part of that file is included below. Its about 1.1GB in size.

It appears to me to follow VCF format at least to a degree but it is missing the rsid's. The field where an rsid would usually be in a VCF file is replaced with a period (.).

I would like to backfill (annotate?) this file with known rsid's (where defined). There are nearly 4 million of them and online tools don't seem ideal, I would prefer to be able to do it locally, if that is possible. I am not looking for software to do this, but a data file which specifies which rsid is assigned to a particular chr:pos:allele. That may or may not be possible - perhaps I misunderstand and its far more complex than that. But from what I can tell, for a particular reference genome, and a particular build of dbsnp, there is a simple relationship as described.

I want to do this, for instance, to be able to use various tools, like SNPedia. Dante appears to use GRCh37, vs SNPedia using mainly GRCh38, which I can see is an issue, and will probably require me to buy a GRCh38 conversion from Dante (or do it myself). But please ignore that particular issue right now.

I don't quite understand what the dependencies are. I've seen comments that suggest that reference genome is not important to generating rsid's from VCF files and that only dbsnp build number needs to match the environment that the VCF file was generated in, but that does not seem to make much sense to me. From what I can see, if you are converting chr:pos information to rsid's, you need to be using the right reference genome. Does dbsnp build really matter? Should you not use the latest build, if it is continually being added to - surely rsid's are not reused? And then a post on the ANNOVAR site seems to suggest that no one really knows what rsid's are and that to different people they are different things, further muddying the waters.

So in terms of question - is there a file, downloadable from somewhere, that contains mappings from chr:pos to rsid, and if so, how does one select the correct file for a particular dataset (IE generated with a particular reference genome, with a particular dbsnp, to be used with references that have used a particular reference/dbsnp build like GRCh38 / dbsnp build 151)?

This file (header first), downloaded from Dante, is named .filtered.snp.vcf

    ##fileformat=VCFv4.2
    ##FILTER=<ID=PASS,Description="All filters passed">
    ##FORMAT=<ID=AD,Number=R,Type=Integer,Description="Allelic depths (counting only informative reads out of the total reads) for the ref and alt alleles in the order listed">
    ##FORMAT=<ID=AF,Number=A,Type=Float,Description="Allele fractions for alt alleles in the order listed">
    ##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Approximate read depth (reads with MQ=255 or with bad mates are filtered)">
    ##FORMAT=<ID=F1R2,Number=R,Type=Integer,Description="Count of reads in F1R2 pair orientation supporting each allele">
    ##FORMAT=<ID=F2R1,Number=R,Type=Integer,Description="Count of reads in F2R1 pair orientation supporting each allele">
    ##FORMAT=<ID=GP,Number=G,Type=Float,Description="Phred-scaled posterior probabilities for genotypes as defined in the VCF specification">
    ##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
    ##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
    ##FORMAT=<ID=MB,Number=4,Type=Integer,Description="Per-sample component statistics to detect mate bias">
    ##FORMAT=<ID=PL,Number=G,Type=Integer,Description="Normalized, Phred-scaled likelihoods for genotypes as defined in the VCF specification">
    ##FORMAT=<ID=PRI,Number=G,Type=Float,Description="Phred-scaled prior probabilities for genotypes">
    ##FORMAT=<ID=PS,Number=1,Type=Integer,Description="Physical phasing ID information, where each unique ID within a given sample (but not across samples) connects records within a phasing group">
    ##FORMAT=<ID=SB,Number=4,Type=Integer,Description="Per-sample component statistics which comprise the Fisher's Exact Test to detect strand bias">
    ##FORMAT=<ID=SQ,Number=A,Type=Float,Description="Somatic quality">
    ##DRAGENCommandLine=<ID=HashTableBuild,Version="SW: 01.003.044.3.9.5, HashTableVersion: 8",CommandLineOptions="dragen --build-hash-table true --output-directory grch37 --ht-reference grch37.fasta">
    ##DRAGENCommandLine=<ID=dragen,Version="SW: 05.121.645.4.0.3, HW: 05.121.645",Date="Thu Mar 02 09:58:00 UTC 2023",CommandLineOptions="-f -c output.cfg -r /references/grch37">
    ##INFO=<ID=AC,Number=A,Type=Integer,Description="Allele count in genotypes, for each ALT allele, in the same order as listed">
    ##INFO=<ID=AF,Number=A,Type=Float,Description="Allele Frequency, for each ALT allele, in the same order as listed">
    ##INFO=<ID=AN,Number=1,Type=Integer,Description="Total number of alleles in called genotypes">
    ##INFO=<ID=FS,Number=1,Type=Float,Description="Phred-scaled p-value using Fisher's exact test to detect strand bias">
    ##INFO=<ID=QD,Number=1,Type=Float,Description="Variant Confidence/Quality by Depth">
    ##INFO=<ID=SOR,Number=1,Type=Float,Description="Symmetric Odds Ratio of 2x2 contingency table to detect strand bias">
    ##INFO=<ID=DP,Number=1,Type=Integer,Description="Approximate read depth (informative and non-informative); some reads may have been filtered based on mapq etc.">
    ##INFO=<ID=FractionInformativeReads,Number=1,Type=Float,Description="The fraction of informative reads out of the total reads">
    ##INFO=<ID=MQ,Number=1,Type=Float,Description="RMS Mapping Quality">
    ##INFO=<ID=MQRankSum,Number=1,Type=Float,Description="Z-score From Wilcoxon rank sum test of Alt vs. Ref read mapping qualities">
    ##INFO=<ID=ReadPosRankSum,Number=1,Type=Float,Description="Z-score from Wilcoxon rank sum test of Alt vs. Ref read position bias">
    ##FILTER=<ID=DRAGENSnpHardQUAL,Description="Set if true:QUAL < 3">
    ##FILTER=<ID=DRAGENIndelHardQUAL,Description="Set if true:QUAL < 3">
    ##FILTER=<ID=LowDepth,Description="Set if true:DP <= 1">
    ##FILTER=<ID=PloidyConflict,Description="Genotype call from variant caller not consistent with chromosome ploidy">
    ##FILTER=<ID=base_quality,Description="Site filtered because median base quality of alt reads at this locus does not meet threshold">
    ##FILTER=<ID=filtered_reads,Description="Site filtered because too large a fraction of reads have been filtered out">
    ##FILTER=<ID=fragment_length,Description="Site filtered because absolute difference between the median fragment length of alt reads and median fragment length of ref reads at this locus exceeds threshold">
    ##FILTER=<ID=low_af,Description="Allele frequency does not meet threshold">
    ##FILTER=<ID=low_depth,Description="Site filtered because the read depth is too low">
    ##FILTER=<ID=low_frac_info_reads,Description="Site filtered because the fraction of informative reads is below threshold">
    ##FILTER=<ID=low_normal_depth,Description="Site filtered because the normal sample read depth is too low">
    ##FILTER=<ID=long_indel,Description="Site filtered because the indel length is too long">
    ##FILTER=<ID=mapping_quality,Description="Site filtered because median mapping quality of alt reads at this locus does not meet threshold">
    ##FILTER=<ID=multiallelic,Description="Site filtered because more than two alt alleles pass tumor LOD">
    ##FILTER=<ID=non_homref_normal,Description="Site filtered because the normal sample genotype is not homozygous reference">
    ##FILTER=<ID=no_reliable_supporting_read,Description="Site filtered because no reliable supporting somatic read exists">
    ##FILTER=<ID=panel_of_normals,Description="Seen in at least one sample in the panel of normals vcf">
    ##FILTER=<ID=read_position,Description="Site filtered because median of distances between start/end of read and this locus is below threshold">
    ##FILTER=<ID=RMxNRepeatRegion,Description="Site filtered because all or part of the variant allele is a repeat of the reference">
    ##FILTER=<ID=str_contraction,Description="Site filtered due to suspected PCR error where the alt allele is one repeat unit less than the reference">
    ##FILTER=<ID=too_few_supporting_reads,Description="Site filtered because there are too few supporting reads in the tumor sample">
    ##FILTER=<ID=weak_evidence,Description="Somatic variant score does not meet threshold">
    ##ML=PackageVersion=4,ModelVersion=4,ModelCRC=2921457924
    ##contig=<ID=1,length=249250621>
    ##contig=<ID=2,length=243199373>
    ##contig=<ID=3,length=198022430>
    ##contig=<ID=4,length=191154276>
    ##contig=<ID=5,length=180915260>
    ##contig=<ID=6,length=171115067>
    ##contig=<ID=7,length=159138663>
    ##contig=<ID=8,length=146364022>
    ##contig=<ID=9,length=141213431>
    ##contig=<ID=10,length=135534747>
    ##contig=<ID=11,length=135006516>
    ##contig=<ID=12,length=133851895>
    ##contig=<ID=13,length=115169878>
    ##contig=<ID=14,length=107349540>
    ##contig=<ID=15,length=102531392>
    ##contig=<ID=16,length=90354753>
    ##contig=<ID=17,length=81195210>
    ##contig=<ID=18,length=78077248>
    ##contig=<ID=19,length=59128983>
    ##contig=<ID=20,length=63025520>
    ##contig=<ID=21,length=48129895>
    ##contig=<ID=22,length=51304566>
    ##contig=<ID=X,length=155270560>
    ##contig=<ID=Y,length=59373566>
    ##contig=<ID=MT,length=16569>
    ##contig=<ID=GL000207.1,length=4262>
    ##contig=<ID=GL000226.1,length=15008>
    ##contig=<ID=GL000229.1,length=19913>
    ##contig=<ID=GL000231.1,length=27386>
    ##contig=<ID=GL000210.1,length=27682>
    ##contig=<ID=GL000239.1,length=33824>
    ##contig=<ID=GL000235.1,length=34474>
    ##contig=<ID=GL000201.1,length=36148>
    ##contig=<ID=GL000247.1,length=36422>
    ##contig=<ID=GL000245.1,length=36651>
    ##contig=<ID=GL000197.1,length=37175>
    ##contig=<ID=GL000203.1,length=37498>
    ##contig=<ID=GL000246.1,length=38154>
    ##contig=<ID=GL000249.1,length=38502>
    ##contig=<ID=GL000196.1,length=38914>
    ##contig=<ID=GL000248.1,length=39786>
    ##contig=<ID=GL000244.1,length=39929>
    ##contig=<ID=GL000238.1,length=39939>
    ##contig=<ID=GL000202.1,length=40103>
    ##contig=<ID=GL000234.1,length=40531>
    ##contig=<ID=GL000232.1,length=40652>
    ##contig=<ID=GL000206.1,length=41001>
    ##contig=<ID=GL000240.1,length=41933>
    ##contig=<ID=GL000236.1,length=41934>
    ##contig=<ID=GL000241.1,length=42152>
    ##contig=<ID=GL000243.1,length=43341>
    ##contig=<ID=GL000242.1,length=43523>
    ##contig=<ID=GL000230.1,length=43691>
    ##contig=<ID=GL000237.1,length=45867>
    ##contig=<ID=GL000233.1,length=45941>
    ##contig=<ID=GL000204.1,length=81310>
    ##contig=<ID=GL000198.1,length=90085>
    ##contig=<ID=GL000208.1,length=92689>
    ##contig=<ID=GL000191.1,length=106433>
    ##contig=<ID=GL000227.1,length=128374>
    ##contig=<ID=GL000228.1,length=129120>
    ##contig=<ID=GL000214.1,length=137718>
    ##contig=<ID=GL000221.1,length=155397>
    ##contig=<ID=GL000209.1,length=159169>
    ##contig=<ID=GL000218.1,length=161147>
    ##contig=<ID=GL000220.1,length=161802>
    ##contig=<ID=GL000213.1,length=164239>
    ##contig=<ID=GL000211.1,length=166566>
    ##contig=<ID=GL000199.1,length=169874>
    ##contig=<ID=GL000217.1,length=172149>
    ##contig=<ID=GL000216.1,length=172294>
    ##contig=<ID=GL000215.1,length=172545>
    ##contig=<ID=GL000205.1,length=174588>
    ##contig=<ID=GL000219.1,length=179198>
    ##contig=<ID=GL000224.1,length=179693>
    ##contig=<ID=GL000223.1,length=180455>
    ##contig=<ID=GL000195.1,length=182896>
    ##contig=<ID=GL000212.1,length=186858>
    ##contig=<ID=GL000222.1,length=186861>
    ##contig=<ID=GL000200.1,length=187035>
    ##contig=<ID=GL000193.1,length=189789>
    ##contig=<ID=GL000194.1,length=191469>
    ##contig=<ID=GL000225.1,length=211173>
    ##contig=<ID=GL000192.1,length=547496>
    ##contig=<ID=NC_007605,length=171823>
    ##contig=<ID=hs37d5,length=35477943>
    ##reference=file:///references/grch37/reference.bin
    ##bcftools_filterVersion=1.9+htslib-1.9
    ##bcftools_filterCommand=filter -i 'TYPE="snp" ' /processing/sample.raw.vcf.gz; Date=Thu Mar  2 10:28:45 2023
    #CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  GFX0xxxxxx
    1   12783   .   G   A   6.16    PASS    AC=1;AF=0.5;AN=2;DP=37;FS=0;MQ=11.6;MQRankSum=-2.137;QD=0.45;ReadPosRankSum=2.878;SOR=1.949;FractionInformativeReads=0.973  GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:10,26:0.722:36:4,14:6,12:4:38,0,1:6.1614,2.5231,7.0198:0,34.77,37.77:10,0,26,0:4,6,14,12
    1   13868   .   A   G   4.69    PASS    AC=1;AF=0.5;AN=2;DP=61;FS=2.275;MQ=16.39;MQRankSum=2.368;QD=0.36;ReadPosRankSum=-1.078;SOR=1.136;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:36,25:0.41:61:16,15:20,10:4:37,0,8:4.6918,2.1112,13.42:0,34.77,37.78:15,21,8,17:23,13,12,13
    1   13980   .   T   C   23.14   PASS    AC=2;AF=1;AN=2;DP=2;FS=0;MQ=5.7;MQRankSum=0;QD=6.5;ReadPosRankSum=0;SOR=2.303;FractionInformativeReads=1    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    1/1:0,2:1:2:0,1:0,1:19:61,24,0:23.143,20.887,0.056833:0,34.77,37.77:0,0,2,0:0,0,1,1
    1   14464   .   A   T   14.11   PASS    AC=1;AF=0.5;AN=2;DP=48;FS=7.243;MQ=22.38;MQRankSum=3.14;QD=0.86;ReadPosRankSum=2.527;SOR=0.267;FractionInformativeReads=1   GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:11,37:0.771:48:8,16:3,21:7:77,6,0:14.113,0.98271,7.8591:0,34.77,37.77:3,8,19,18:4,7,20,17
    1   14907   .   A   G   3.67    PASS    AC=1;AF=0.5;AN=2;DP=77;FS=0.893;MQ=26.78;MQRankSum=-0.058;QD=0.41;ReadPosRankSum=5.956;SOR=0.532;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:48,29:0.377:77:31,13:17,16:4:36,0,22:3.6675,2.453,27.555:0,34.77,37.77:27,21,15,14:28,20,15,14
    1   14930   .   A   G   3.26    PASS    AC=1;AF=0.5;AN=2;DP=82;FS=3.051;MQ=27.04;MQRankSum=-0.363;QD=0.37;ReadPosRankSum=-0.277;SOR=0.754;FractionInformativeReads=1    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:51,31:0.378:82:30,16:21,15:3:35,0,21:3.2627,2.7897,26.606:0,34.77,37.77:28,23,14,17:29,22,17,14
    1   15211   .   T   G   6.66    PASS    AC=1;AF=0.5;AN=2;DP=70;FS=0;MQ=21.21;MQRankSum=0.533;QD=0.49;ReadPosRankSum=-0.905;SOR=0.743;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:28,42:0.6:70:12,22:16,20:6:40,0,11:6.6596,1.2235,15.27:0,34.77,37.77:15,13,23,19:15,13,17,25

Then a couple of chunks from random bits of the file:

offset 100000 lines:

    1   82279586    .   T   C   62.12   PASS    AC=2;AF=1;AN=2;DP=25;FS=0;MQ=250;MQRankSum=0;QD=4.89;ReadPosRankSum=0;SOR=0.941;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,25:1:25:0,10:0,15:52:100,55,0:62.124,52.173,2.8992e-05:0,34.77,37.77:0,0,11,14:0,0,10,15
    1   82280153    .   G   A   62.73   PASS    AC=2;AF=1;AN=2;DP=31;FS=0;MQ=250;MQRankSum=0;QD=4.53;ReadPosRankSum=0;SOR=1.044;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,31:1:31:0,15:0,16:57:101,61,0:62.731,57.663,9.319e-06:0,34.77,37.77:0,0,18,13:0,0,14,17
    1   82280683    .   C   T   54.83   PASS    AC=2;AF=1;AN=2;DP=17;FS=0;MQ=250;MQRankSum=0;QD=5.77;ReadPosRankSum=0;SOR=1.061;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,17:1:17:0,8:0,9:45:93,49,0:54.832,45.914,0.00012581:0,34.77,37.77:0,0,10,7:0,0,11,6
    1   82280870    .   T   C   41.24   PASS    AC=1;AF=0.5;AN=2;DP=19;FS=1.808;MQ=250;MQRankSum=3.592;QD=2.58;ReadPosRankSum=1.858;SOR=1.14;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:11,8:0.421:19:8,2:3,6:40:76,0,45:41.236,0.00039866,47.817:0,34.77,37.77:4,7,2,6:5,6,6,2
    1   82280975    .   A   T   40.52   PASS    AC=1;AF=0.5;AN=2;DP=24;FS=1.569;MQ=250;MQRankSum=4.07;QD=2.02;ReadPosRankSum=1.025;SOR=1.085;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:14,10:0.417:24:5,8:9,2:40:75,0,50:40.524,0.00040643,53.142:0,34.77,37.77:7,7,4,6:8,6,4,6
    1   82281492    .   T   C   44.14   PASS    AC=1;AF=0.5;AN=2;DP=22;FS=6.923;MQ=250;MQRankSum=3.873;QD=2.2;ReadPosRankSum=-0.334;SOR=2.064;FractionInformativeReads=1    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:13,9:0.409:22:9,4:4,5:44:79,0,51:44.143,0.00018328,54.3:0,34.77,37.77:6,7,7,2:6,7,6,3
    1   82282062    .   T   C   8.34    PASS    AC=1;AF=0.5;AN=2;DP=26;FS=1.974;MQ=247.84;MQRankSum=3.383;QD=1.04;ReadPosRankSum=-1.561;SOR=0.863;FractionInformativeReads=1    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:21,5:0.192:26:9,3:12,2:8:42,0,29:8.3401,0.69098,32.665:0,34.77,37.77:9,12,3,2:9,12,1,4
    1   82282378    .   C   T   42.61   PASS    AC=1;AF=0.5;AN=2;DP=24;FS=1.685;MQ=250;MQRankSum=4.07;QD=2.02;ReadPosRankSum=2.898;SOR=1.23;FractionInformativeReads=1  GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:14,10:0.417:24:8,6:6,4:43:77,0,57:42.606,0.0002423,60.191:0,34.77,37.77:8,6,7,3:10,4,7,3
    1   82282788    .   C   T   48.99   PASS    AC=1;AF=0.5;AN=2;DP=30;FS=0;MQ=245.71;MQRankSum=2.889;QD=1.67;ReadPosRankSum=1.975;SOR=0.914;FractionInformativeReads=1 GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:14,16:0.533:30:6,9:8,7:46:84,0,46:48.985,0.00010665,49.217:0,34.77,37.77:8,6,10,6:6,8,7,9
    1   82282947    .   G   T   47.31   PASS    AC=1;AF=0.5;AN=2;DP=32;FS=10.606;MQ=250;MQRankSum=4.796;QD=1.56;ReadPosRankSum=0.68;SOR=0.725;FractionInformativeReads=1    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    0/1:15,17:0.531:32:6,8:9,9:46:82,0,48:47.312,0.00011286,51.287:0,34.77,37.77:11,4,7,10:8,7,8,9

offset 2000000 lines:

    8   142260995   .   T   C   65.76   PASS    AC=2;AF=1;AN=2;DP=28;FS=0;MQ=248.93;MQRankSum=0;QD=4.69;ReadPosRankSum=0;SOR=1.179;FractionInformativeReads=1   GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    1/1:0,28:1:28:0,14:0,14:53:104,56,0:65.759,52.946,2.3297e-05:0,34.77,37.77:0,0,11,17:0,0,13,15
    8   142261389   .   T   C   66.03   PASS    AC=2;AF=1;AN=2;DP=28;FS=0;MQ=250;MQRankSum=0;QD=4.67;ReadPosRankSum=0;SOR=0.997;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,28:1:28:0,12:0,16:54:104,57,0:66.026,54.004,1.812e-05:0,34.77,37.77:0,0,16,12:0,0,14,14
    8   142261541   .   G   A   65.34   PASS    AC=2;AF=1;AN=2;DP=30;FS=0;MQ=250;MQRankSum=0;QD=4.51;ReadPosRankSum=0;SOR=1.765;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,30:1:30:0,14:0,16:52:103,55,0:65.339,52.429,2.5886e-05:0,34.77,37.77:0,0,9,21:0,0,17,13
    8   142261618   .   A   G   64.46   PASS    AC=2;AF=1;AN=2;DP=27;FS=0;MQ=250;MQRankSum=0;QD=4.57;ReadPosRankSum=0;SOR=1.302;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB:PS   1|1:0,27:1:27:0,14:0,13:53:102,56,0:64.465,53.434,2.1227e-05:0,34.77,37.78:0,0,10,17:0,0,15,12:142261618
    8   142261621   .   C   T   60.33   PASS    AC=2;AF=1;AN=2;DP=25;FS=0;MQ=250;MQRankSum=0;QD=4.7;ReadPosRankSum=0;SOR=1.051;FractionInformativeReads=0.96    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB:PS 1|1:0,24:1:24:0,13:0,11:48:98,51,0:60.326,48.044,7.2481e-05:0,34.77,37.78:0,0,10,14:0,0,13,11:142261618
    8   142261967   .   C   T   68.75   PASS    AC=2;AF=1;AN=2;DP=37;FS=0;MQ=250;MQRankSum=0;QD=4.26;ReadPosRankSum=0;SOR=0.859;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,37:1:37:0,17:0,20:56:107,60,0:68.745,56.755,9.8367e-06:0,34.77,37.77:0,0,20,17:0,0,18,19
    8   142262327   .   A   G   68.48   PASS    AC=2;AF=1;AN=2;DP=42;FS=0;MQ=250;MQRankSum=0;QD=4.13;ReadPosRankSum=0;SOR=0.788;FractionInformativeReads=1GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB  1/1:0,42:1:42:0,21:0,21:56:106,60,0:68.476,56.786,9.8367e-06:0,34.77,37.77:0,0,22,20:0,0,11,31
    8   142262413   .   G   A   62.88   PASS    AC=2;AF=1;AN=2;DP=27;FS=0;MQ=231.86;MQRankSum=0;QD=4.75;ReadPosRankSum=0;SOR=1.302;FractionInformativeReads=1   GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    1/1:0,27:1:27:0,16:0,11:51:101,54,0:62.881,51.183,3.5205e-05:0,34.77,37.77:0,0,10,17:0,0,9,18
    8   142263074   .   G   A   51.24   PASS    AC=2;AF=1;AN=2;DP=20;FS=0;MQ=143.27;MQRankSum=0;QD=5.31;ReadPosRankSum=0;SOR=0.892;FractionInformativeReads=1   GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    1/1:0,20:1:20:0,9:0,11:36:89,39,0:51.237,35.702,0.0012015:0,34.77,37.77:0,0,11,9:0,0,10,10
    8   142263171   .   C   T   52.9    PASS    AC=2;AF=1;AN=2;DP=22;FS=0;MQ=122.2;MQRankSum=0;QD=4.81;ReadPosRankSum=0;SOR=0.693;FractionInformativeReads=0.909    GT:AD:AF:DP:F1R2:F2R1:GQ:PL:GP:PRI:SB:MB    1/1:0,20:1:20:0,11:0,9:35:91,38,0:52.899,34.894,0.0014297:0,34.77,37.77:0,0,10,10:0,0,9,11
$\endgroup$

3 Answers 3

3
$\begingroup$

Download the dbSNP VCF file (latest version) and map to your input VCF using bcftools annotate. You'll want to match reference sequence version (hg19/hg38) and use the latest version of dbSNP to get the most up to date information.

dbSNP VCF FTP folder: https://ftp.ncbi.nih.gov/snp/latest_release/VCF/

$\endgroup$
5
  • $\begingroup$ Thanks, I eventually got this working although the brittleness with which it treats input file compression standards is rather annoying and took me ages to work out, hence my delay. I've finally annotated the vcf with rsid with bcftools, so thanks very much. Since my data is refed to GRCh37, and SNPedia is made with 38, although it works I wonder if it's worth 'getting' 38 somehow. Looks like generating it from the BAM or fastq would be the cleanest way to do this. Could I ask you to comment on that, and recommend any software or guides? $\endgroup$
    – Pete
    Commented Mar 14, 2023 at 20:33
  • $\begingroup$ Yes, re-aligning to hg38 is the way to go. liftOver is an option too, but it has its shortcomings. $\endgroup$
    – Ram RS
    Commented Mar 15, 2023 at 21:22
  • $\begingroup$ So regarding the realigning, is that a massive undertaking or relatively simple? I imagine it's simple - just testing the original files (bam or fastq etc) against the new reference genome hg38. Or does it require detailed info about the original sequencing operation that I am not likely to be able to get out of Dante (as in gringers reply re illumina snpchip versions being matched with specific dbsnp builds etc)? $\endgroup$
    – Pete
    Commented Mar 16, 2023 at 22:12
  • $\begingroup$ Since it's WGS, my guess is that realigning should work just fine. Targeted sequencing might have created problems. You can ask @gringer for their take. $\endgroup$
    – Ram RS
    Commented Mar 17, 2023 at 14:14
  • $\begingroup$ RamRS - thanks. Gringer, wonder if you come by these parts again if you could comment? As an extra ask, hoping its not too close to the wind to ask on SE about software recommendations, any suggestions for a tool to do this (I'm on Linux so Windows only software is out for me)? $\endgroup$
    – Pete
    Commented Mar 20, 2023 at 12:20
2
$\begingroup$

As a partial answer to the specific questions about definitions:

  • rsIDs define a specific location and type of variation*, but the location of that variation can change within different reference genomes. This may include a change in chromosome, as well as a change in position.

  • Different dbSNP builds can add additional variation, but the same rsID will always refer to the same variant in all dbSNP builds that refer to that variation.

  • rsIDs can be invalidated in new dbSNP builds for various reasons. It's a good idea to always use the latest dbSNP build when annotating variants.

  • rsIDs are never recycled. There is a prescribed process for merging and splitting variants that ensures backwards compatibility. If a variation is merged, I think the approach is to always define the earliest-numbered variant as the primary variation. If a variation is split, I think the approach is to use a new rsID for each of the separate variations.

* by "variation", I'm referring to a type of modification in the DNA sequence, rather than a specific instance of that modification. This is explained in the dbSNP documentation (an "ss" is a submitted variant, a submission made to dbSNP):

Each ss that maps to the same position and is the same variant type is assigned to an existing RefSNP (rs) or assigned a new one.

Although that isn't entirely clear, this suggests that if a new variant is submitted that falls at the same position of an exiting one and is of the same type as the existing one, then the new one may simply be added to the same rsID as the existing one. So, for example, if we have a single nucleotide G => C variant at position N, and a new G => T variant is reported at the same position, then the new one could have the same rsID as the old one.


In response to an additional comment / question:

a post on the ANNOVAR site seems to suggest that no one really knows what rsid's are and that to different people they are different things, further muddying the waters.... Is this guy an authority on this or not?

I can't comment on authority for that guy; the authorship is unclear to me. What matters is how the dbSNP people understand dbSNP (i.e. the documentation link above), because they are the authority on it.

My understanding of it was helped because many years ago I attended a close-knit workshop that included a talk by Mike Feolo, and I had been working on a human dataset with many oddities, so I was able to ask lots of questions about corner cases. Here's the text from a section on rsIDs that I added into the appendix of my thesis; things have probably changed since then, but I expect the broad strokes are still correct:

Mike Feolo gave a talk at the [Working With The HapMap] Cambridge course on the dbSNP database. He talked about how new Single Nucleotide Polymorphisms (SNPs) were submitted, compared against currently existing SNPs, and merged if they were already in the database. The database RefSeq(rs) numbers are updated about twice a year, and whenever there is a new build of the human genome sequence. As a result of these builds, some rs numbers are merged to account for increased knowledge of mutations, a process that I was not aware of. This has an impact on my research, because I have received data from an Illumina 317k SNPchip, which I compare to the HapMap dataset. The rs numbers from the Illumina chip are for a specific dbSNP build, which may not be fully consistent with the most recent build (from which the HapMap data is retrieved). Mike has recommended that both the build number and the rs number are used when referring to a specific mutation. Sharma Buch pointed out that it was possible to ask Illumina for information about the merged rs numbers for their chips, which would allow datasets from different builds to be used together.

$\endgroup$
3
  • 1
    $\begingroup$ I think rsIDs are chr-pos specific, which is why multiple ALT alleles can exist in the same entry. $\endgroup$
    – Ram RS
    Commented Mar 14, 2023 at 19:07
  • 1
    $\begingroup$ Thanks for all of this from all 3 of you. It's gradually helping me to form a picture. If anyone feels like reading and commenting on annovar.openbioinformatics.org/en/latest/articles/dbSNP it would be really useful. Is this guy an authority on this or not? His point that rsid means different thing to different people is revealing and if correct could prevent much head scratching and wasted hours reading and pondering this. It implies that one just has to get used to how rsid is used without trying to understand it too much (a sort of quantum mechanics like detachment...) $\endgroup$
    – Pete
    Commented Mar 14, 2023 at 20:27
  • $\begingroup$ @Gringer your update is very interesting thanks. $\endgroup$
    – Pete
    Commented Mar 15, 2023 at 9:31
1
$\begingroup$

Generally, consumer WGS companies do not deliver annotated VCF's. Annotations include rsIDs, gene names or similar area designations, and others. But this can be done using the NIH dbSNP VCF files and the bcftools annotate command. The annotation files at UCSC are modified for their GenomeBrowser and are often filtered and cleaned up. While better and more consistent, they are usually not the latest release available from NIH.

Just make sure to match the build model and sequence (chromosome) naming convention between your annotation file and the original VCF. Dante usually delivers using the hs37d5 reference model. So use the .25 dbSNP VCFs files. But dbSNP VCFs use accession names whereas Dante's hs37d5 reference model uses GRCh numeric ones. So you will have to rename sequence names in one file or the other to make it work. Also, hopefully, the vendor has left aligned the InDels as that is what the dbSNP file has them defined as and what bcftools expects as well. Dante traditionally uses the DRAGEN pipeline with their Illumina sequencers; which does this. But they have recently switched to MGI DNB sequencers and the megaBOLT pipeline for customers outside the USA. We are not yet sure how they are processing InDels. Finally, Dante has historically delivered a separate SNP and InDel file. You will often want to merge them or annotate them both separately then merge. bcftools can do this.

The dbSNP VCF file has well over 1 billion rsIDs now defined. And many fields of additional information that can be annotated in. Often, if you just want to add a simple subset of the data, you can process the dbSNP VCF file to a simpler tab-separated (TSV, tab) file of entries. It can cut down the annotation file from 25 GB to 5 GB or less. See the bcftools annotation page for hints on using that tool to do this as well.

Personally, I find it easier to use online tools or services that have been tweaked. Especially if they have nice browsers. Sequencing dot com, usegalaxy.org and iobio.io are a few that come to mind. Not all allow the export of the annotated VCF.

This and similar topics of utilizing your results are heavily discussed in the consumerWGS Facebook group; if you also wish to participate there.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.