text: user2998: After reading the GIAB paper in https://www.biorxiv.org/content/early/2018/05/25/281006 and its Figure 1, I am still having trouble understanding the data inside the GIAB VCF file for HG001 (HG001_GRCh38_GIAB_highconf_CG-IllFB-IllGATKHC-Ion-10X-SOLID_CHROM1-X_v.3.3.2_all.vcf). In particular, I need to understand the information that led to a FILTER value of PASS for some variants.
  The first case that I would like to discuss is the following variant:
  chr1    110176467   .   A   G   50  PASS    platforms=3;platformnames=Illumina,CG,10X;datasets=3;datasetnames=HiSeqPE300x,CGnormal,10XChromium;callsets=4;callsetnames=HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo;datasetsmissingcall=IonExome,SolidPE50x50bp,SolidSE75bp;callable=CS_HiSeqPE300xGATK_callable,CS_CGnormal_callable,CS_HiSeqPE300xfreebayes_callable    GT:PS:DP:ADALL:AD:GQ    0/1:.:512:114,106:145,148:760
  
  Translating the information in the line to English, this is what I get: 
  
  We considered data from the following sequencing technologies: Illumina,CG,10X,IonExome,SolidPE50x50bp, and SolidSE75bp. 
  Of those, some of them (IonExome,SolidPE50x50bp,SolidSE75bp) were missing the call, so we discarded them. 
  Therefore, the data for the variant is going to come from 3 platforms (Illumina,CG,10X), 3 datasets (HiSeqPE300x,CGnormal,10XChromium), and 4 callsets (HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo), because we analyzed the data from HiSeqPE300x with GATK and also with Freebayes.
  However, the callsets from CS_HiSeqPE300xGATK,CS_CGnormal, and CS_HiSeqPE300xfreebayes had this call in a region with low coverage of high MQ reads.
  Across all the platforms, DP=512
  ADALL (Values 114,106) is coming from “all" the datasets. Are those “all” the set (HiSeqPE300x + Cgnormal + 10XChromium), or the set (Illumina + CG + 10X + IonExome + SolidPE50x50bp + SolidSE75bp).
  AD (Values 145,148) is coming from the same “all” datasets as ADALL, but they are unfiltered, so 145>114 and 148>106.
  
  All this make sense, but I would like to know what case of Figure 1c was applied to arrive at the PASS value.
  The second variant that I would like to mention is this one:
  chr1    5705293 .   T   C   50  PASS    platforms=3;platformnames=Illumina,CG,10X;datasets=3;datasetnames=HiSeqPE300x,CGnormal,10XChromium;callsets=4;callsetnames=HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo;datasetsmissingcall=IonExome,SolidPE50x50bp,SolidSE75bp;callable=CS_HiSeqPE300xGATK_callable,CS_CGnormal_callable,CS_HiSeqPE300xfreebayes_callable;filt=CS_HiSeqPE300xGATK_filt   GT:PS:DP:ADALL:AD:GQ    0/1:.:599:121,139:49,28:627
  
  Again, my translation would be:
  
  We considered data from the following sequencing technologies: Illumina,CG,10X,IonExome,SolidPE50x50bp, and SolidSE75bp. 
  Of those, some of them (IonExome,SolidPE50x50bp,SolidSE75bp) were missing the call, so we discarded them. 
  Therefore, the data for the variant is going to come from 3 platforms (Illumina,CG,10X), 3 datasets (HiSeqPE300x,CGnormal,10XChromium), and 4 callsets (HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo), because we analyzed the data from HiSeqPE300x with GATK and also with Freebayes.
  However, the callsets from CS_HiSeqPE300xGATK,CS_CGnormal, and CS_HiSeqPE300xfreebayes had this call in a region with low coverage of high MQ reads.
  Also, the data from CS_HiSeqPE300xGATK_filt was filtered/discarded
  Across all the platforms, DP=599
  ADALL (Values 121,139) is coming from “all" the datasets. I assume (Cgnormal + 10XChromium), because CS_HiSeqPE300xGATK_filt was filtered. 
  AD (Values 49,28) are unfiltered.
  
  But notice that 49<121 and 28<139. How come that there are less unfiltered reads (AD)? Should not be the case that AD > ADALL for all instances? 
  I don’t know how to reconcile the two variants. Aren’t they providing opposite information in the AD and ADALL fields? Which specific datasets and callsets are involved in getting AD and ADALL for each variant?
  I would also appreciate answers to a couple of other questions: 
  
  Why are these two variants a PASS if there are 3 “callable” datasets in regions of low coverage? 
  Did the variants follow an arbitration process? I guess not, if the INFO flag “arbitrated” is missing.
  
  For reference, here are the descriptions in the header of the VCF file that I have been using:
  ##INFO=<ID=platforms,Number=1,Type=Integer,Description="Number of different platforms for which at least one callset called this genotype, whether filtered or not">
  ##INFO=<ID=platformnames,Number=.,Type=String,Description="Names of platforms for which at least one callset called this genotype, whether filtered or not">
  ##INFO=<ID=datasets,Number=1,Type=Integer,Description="Number of different datasets for which at least one callset called this genotype, whether filtered or not">
  ##INFO=<ID=datasetnames,Number=.,Type=String,Description="Names of datasets for which at least one callset called this genotype, whether filtered or not">
  ##FORMAT=<ID=ADALL,Number=R,Type=Integer,Description="Net allele depths across all datasets">
  ##FORMAT=<ID=AD,Number=R,Type=Integer,Description="Net allele depths across all unfiltered datasets with called genotype">
  ##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Total read depth summed across all datasets, excluding MQ0 reads">
  ##INFO=<ID=datasetsmissingcall,Number=.,Type=Integer,Description="Names of datasets that are missing a call or have an incorrect call at this location, and the high-confidence call is a variant">
  ##INFO=<ID=callable,Number=.,Type=String,Description="List of callsets that had this call in a region with low coverage of high MQ reads.">
  ##INFO=<ID=filt,Number=.,Type=String,Description="List of callsets that had this call filtered.">
  ##INFO=<ID=arbitrated,Number=1,Type=String,Description="TRUE if callsets had discordant calls so that arbitration was needed.">
  
  
  user964: The filtering criteria should be defined in the header of the VCF file. Can you include the header from your VCF file?
  VCF Format and example header: 
  FILTERs that have been applied to the data should be described as follows:
  
  ##FILTER=<ID=ID,Description=”description”>
  
  header:
  ##fileformat=VCFv4.0
  ##fileDate=20090805
  ##source=myImputationProgramV3.1
  ##reference=1000GenomesPilot-NCBI36
  ##phasing=partial
  ##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">
  ##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
  ##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">
  ##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
  ##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">
  ##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">
  ##FILTER=<ID=q10,Description="Quality below 10">
  ##FILTER=<ID=s50,Description="Less than 50% of samples have data">
  ##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
  ##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
  ##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">
  ##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">
  
  
original_text: 2998: <p>After reading the GIAB paper in <a href="https://www.biorxiv.org/content/early/2018/05/25/281006" rel="nofollow noreferrer">https://www.biorxiv.org/content/early/2018/05/25/281006</a> and its Figure 1, I am still having trouble understanding the data inside the GIAB VCF file for HG001 (HG001_GRCh38_GIAB_highconf_CG-IllFB-IllGATKHC-Ion-10X-SOLID_CHROM1-X_v.3.3.2_all.vcf). In particular, I need to understand the information that led to a FILTER value of PASS for some variants.</p>
  
  <p>The first case that I would like to discuss is the following variant:</p>
  
  <pre><code>chr1    110176467   .   A   G   50  PASS    platforms=3;platformnames=Illumina,CG,10X;datasets=3;datasetnames=HiSeqPE300x,CGnormal,10XChromium;callsets=4;callsetnames=HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo;datasetsmissingcall=IonExome,SolidPE50x50bp,SolidSE75bp;callable=CS_HiSeqPE300xGATK_callable,CS_CGnormal_callable,CS_HiSeqPE300xfreebayes_callable    GT:PS:DP:ADALL:AD:GQ    0/1:.:512:114,106:145,148:760
  </code></pre>
  
  <p>Translating the information in the line to English, this is what I get: </p>
  
  <ul>
  <li>We considered data from the following sequencing technologies: Illumina,CG,10X,IonExome,SolidPE50x50bp, and SolidSE75bp. </li>
  <li>Of those, some of them (IonExome,SolidPE50x50bp,SolidSE75bp) were missing the call, so we discarded them. </li>
  <li>Therefore, the data for the variant is going to come from 3 platforms (Illumina,CG,10X), 3 datasets (HiSeqPE300x,CGnormal,10XChromium), and 4 callsets (HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo), because we analyzed the data from HiSeqPE300x with GATK and also with Freebayes.</li>
  <li>However, the callsets from CS_HiSeqPE300xGATK,CS_CGnormal, and CS_HiSeqPE300xfreebayes had this call in a region with low coverage of high MQ reads.</li>
  <li>Across all the platforms, DP=512</li>
  <li>ADALL (Values 114,106) is coming from “all" the datasets. Are those “all” the set (HiSeqPE300x + Cgnormal + 10XChromium), or the set (Illumina + CG + 10X + IonExome + SolidPE50x50bp + SolidSE75bp).</li>
  <li>AD (Values 145,148) is coming from the same “all” datasets as ADALL, but they are unfiltered, so 145>114 and 148>106.</li>
  </ul>
  
  <p>All this make sense, but I would like to know what case of Figure 1c was applied to arrive at the PASS value.</p>
  
  <p>The second variant that I would like to mention is this one:</p>
  
  <pre><code>chr1    5705293 .   T   C   50  PASS    platforms=3;platformnames=Illumina,CG,10X;datasets=3;datasetnames=HiSeqPE300x,CGnormal,10XChromium;callsets=4;callsetnames=HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo;datasetsmissingcall=IonExome,SolidPE50x50bp,SolidSE75bp;callable=CS_HiSeqPE300xGATK_callable,CS_CGnormal_callable,CS_HiSeqPE300xfreebayes_callable;filt=CS_HiSeqPE300xGATK_filt   GT:PS:DP:ADALL:AD:GQ    0/1:.:599:121,139:49,28:627
  </code></pre>
  
  <p>Again, my translation would be:</p>
  
  <ul>
  <li><p>We considered data from the following sequencing technologies: Illumina,CG,10X,IonExome,SolidPE50x50bp, and SolidSE75bp. </p></li>
  <li><p>Of those, some of them (IonExome,SolidPE50x50bp,SolidSE75bp) were missing the call, so we discarded them. </p></li>
  <li><p>Therefore, the data for the variant is going to come from 3 platforms (Illumina,CG,10X), 3 datasets (HiSeqPE300x,CGnormal,10XChromium), and 4 callsets (HiSeqPE300xGATK,CGnormal,HiSeqPE300xfreebayes,10XGATKhaplo), because we analyzed the data from HiSeqPE300x with GATK and also with Freebayes.</p></li>
  <li>However, the callsets from CS_HiSeqPE300xGATK,CS_CGnormal, and CS_HiSeqPE300xfreebayes had this call in a region with low coverage of high MQ reads.</li>
  <li>Also, the data from CS_HiSeqPE300xGATK_filt was filtered/discarded</li>
  <li>Across all the platforms, DP=599</li>
  <li>ADALL (Values 121,139) is coming from “all" the datasets. I assume (Cgnormal + 10XChromium), because CS_HiSeqPE300xGATK_filt was filtered. </li>
  <li>AD (Values 49,28) are unfiltered.</li>
  </ul>
  
  <p>But notice that 49&lt;121 and 28&lt;139. How come that there are less unfiltered reads (AD)? Should not be the case that AD > ADALL for all instances? </p>
  
  <p>I don’t know how to reconcile the two variants. Aren’t they providing opposite information in the AD and ADALL fields? Which specific datasets and callsets are involved in getting AD and ADALL for each variant?</p>
  
  <p>I would also appreciate answers to a couple of other questions: </p>
  
  <ul>
  <li><p>Why are these two variants a PASS if there are 3 “callable” datasets in regions of low coverage? </p></li>
  <li><p>Did the variants follow an arbitration process? I guess not, if the INFO flag “arbitrated” is missing.</p></li>
  </ul>
  
  <p>For reference, here are the descriptions in the header of the VCF file that I have been using:</p>
  
  <pre><code>##INFO=&lt;ID=platforms,Number=1,Type=Integer,Description="Number of different platforms for which at least one callset called this genotype, whether filtered or not"&gt;
  ##INFO=&lt;ID=platformnames,Number=.,Type=String,Description="Names of platforms for which at least one callset called this genotype, whether filtered or not"&gt;
  ##INFO=&lt;ID=datasets,Number=1,Type=Integer,Description="Number of different datasets for which at least one callset called this genotype, whether filtered or not"&gt;
  ##INFO=&lt;ID=datasetnames,Number=.,Type=String,Description="Names of datasets for which at least one callset called this genotype, whether filtered or not"&gt;
  ##FORMAT=&lt;ID=ADALL,Number=R,Type=Integer,Description="Net allele depths across all datasets"&gt;
  ##FORMAT=&lt;ID=AD,Number=R,Type=Integer,Description="Net allele depths across all unfiltered datasets with called genotype"&gt;
  ##FORMAT=&lt;ID=DP,Number=1,Type=Integer,Description="Total read depth summed across all datasets, excluding MQ0 reads"&gt;
  ##INFO=&lt;ID=datasetsmissingcall,Number=.,Type=Integer,Description="Names of datasets that are missing a call or have an incorrect call at this location, and the high-confidence call is a variant"&gt;
  ##INFO=&lt;ID=callable,Number=.,Type=String,Description="List of callsets that had this call in a region with low coverage of high MQ reads."&gt;
  ##INFO=&lt;ID=filt,Number=.,Type=String,Description="List of callsets that had this call filtered."&gt;
  ##INFO=&lt;ID=arbitrated,Number=1,Type=String,Description="TRUE if callsets had discordant calls so that arbitration was needed."&gt;
  </code></pre>
   964: <p>The filtering criteria should be defined in the header of the VCF file. Can you include the header from your VCF file?</p>
  
  <p><a href="http://www.internationalgenome.org/wiki/Analysis/vcf4.0/" rel="nofollow noreferrer">VCF Format and example header:</a> </p>
  
  <pre><code>FILTERs that have been applied to the data should be described as follows:
  
  ##FILTER=&lt;ID=ID,Description=”description”&gt;
  </code></pre>
  
  <p>header:</p>
  
  <pre><code>##fileformat=VCFv4.0
  ##fileDate=20090805
  ##source=myImputationProgramV3.1
  ##reference=1000GenomesPilot-NCBI36
  ##phasing=partial
  ##INFO=&lt;ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data"&gt;
  ##INFO=&lt;ID=DP,Number=1,Type=Integer,Description="Total Depth"&gt;
  ##INFO=&lt;ID=AF,Number=.,Type=Float,Description="Allele Frequency"&gt;
  ##INFO=&lt;ID=AA,Number=1,Type=String,Description="Ancestral Allele"&gt;
  ##INFO=&lt;ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129"&gt;
  ##INFO=&lt;ID=H2,Number=0,Type=Flag,Description="HapMap2 membership"&gt;
  ##FILTER=&lt;ID=q10,Description="Quality below 10"&gt;
  ##FILTER=&lt;ID=s50,Description="Less than 50% of samples have data"&gt;
  ##FORMAT=&lt;ID=GT,Number=1,Type=String,Description="Genotype"&gt;
  ##FORMAT=&lt;ID=GQ,Number=1,Type=Integer,Description="Genotype Quality"&gt;
  ##FORMAT=&lt;ID=DP,Number=1,Type=Integer,Description="Read Depth"&gt;
  ##FORMAT=&lt;ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality"&gt;
  </code></pre>
  