We've moved!
This site is now read-only. You can find our new documentation site and support forum for posting questions here.
Be sure to read our welcome blog!

# Does the CalculateGenotypePosteriors change parents' genotype?

torontoMember

For example, after the gatk recalibration step, the genotype of a parent, which was 0/0, was changed to 0/1 after running CalculateGenotypePosteriors?

@zwang
Hi,

Yes, it is possible for the genotype to change after CalculateGenotypePosteriors. Have a look at the documentation for more information.

-Sheila

• torontoMember

We analyzed a trio following the Genotype Refinement workflow: https://www.broadinstitute.org/gatk/guide/article?id=4723.
We found a weird record.

Before applying CalculateGenotypePosteriors, the mother has genotype 0/0 (2nd to the last column)
chr1 96925200 . A T 530.13 PASS AC=1;AF=0.167;AN=6;BaseQRankSum=0.599;ClippingRankSum=1.20;DP=87;FS=0.000;GQ_MEAN=486.00;MLEAC=1;MLEAF=0.167;MQ=60.00;MQ0=0;MQRankSum=-8.750e-01;NCC=0;QD=18.93;ReadPosRankSum=0.921;SOR=0.874;VQSLOD=16.04;culprit=MQ GT:AB:AD:DP:GQ:PL 0/1:0.460:13,15:28:99:561,0,486 0/0:.:25,0:25:34:0,34,969 0/0:.:31,0:31:68:0,68,1097

After applying CalculateGenotypePosteriors, the mother's genotype is 0/1,
chr1 96925200 . A T 530.13 PASS AC=2;AF=0.333;AN=6;BaseQRankSum=0.599;ClippingRankSum=1.20;DP=87;FS=0.000;GQ_MEAN=486.00;MLEAC=1;MLEAF=0.167;MQ=60.00;MQ0=0;MQRankSum=-8.750e-01;NCC=0;QD=18.93;ReadPosRankSum=0.921;SOR=0.874;VQSLOD=16.04;culprit=MQ GT:AB:AD:DP:GQ:JL:JP:PL:PP 0/1:0.460:13,15:28:99:0:24:561,0,486:527,0,545 0/1:.:25,0:25:25:0:24:0,34,969:25,0,935 0/0:.:31,0:31:34:0:24:0,68,1097:0,34,1063

But we know, for sure, that the mother's genotype is homo reference, the kid has a de novo call. Here is the command I used for running the CalculateGenotypePosteriors

java -Xmx5g -jar GenomeAnalysisTK.jar
-T CalculateGenotypePosteriors
-R ucsc_hg19.fa
-V joint_recal.vcf
--skipPopulationPriors
-ped fam.ped
-o joint_recal_gposteriors.vcf

We used samtools mpileup to look at the reads at that location on both parents, there are no supporting reads for alternative allele.

Any idea on what is going on?

@zwang
Hi,

Can you tell me how you generated the trio VCF? Also, can you please post a screenshot of the bam file in that region? Please also post the bamout file in the region.

Thanks,
Sheila

• torontoMember

For each sample, here the command to call

java -Xmx11G -Djava.io.tmpdir=/localhd/468805[9]/tmp.eWZ7LiizrR
-jar GenomeAnalysisTK.jar
-T HaplotypeCaller
-R ucsc_hg19.fa
-o child.chr8.htc.vcf.ungenotyped.VCF
-GQB 10 -GQB 30 -stand_emit_conf 10 -A AlleleBalance -A AlleleBalanceBySample -A BaseCounts -A ClippingRankSumTest -A Coverage -A DepthPerAlleleBySample -A DepthPerSampleHC -A FisherStrand -A GCContent -A GenotypeSummaries -A HardyWeinberg -A HomopolymerRun -A InbreedingCoeff -A LowMQ -A MappingQualityRankSumTest -A MappingQualityZero -A MappingQualityZeroBySample -A NBaseCount -A QualByDepth -A RMSMappingQuality -A ReadPosRankSumTest -A SpanningDeletions -A StrandBiasBySample -A StrandOddsRatio -A TandemRepeatAnnotator -A VariantType
--dbsnp dbsnp_138.ucsc.hg19.vcf --genotyping_mode DISCOVERY --emitRefConfidence GVCF --variant_index_type LINEAR --variant_index_parameter 128000 -L chr8 -nct 8 -I input.bam

Since I made calls by chromosome, I next merged the calls.

Then the joint genotyping is:
java -Xmx5g -jar GenomeAnalysisTK.jar -R ucsc_hg19.fa -T GenotypeGVCFs -o joint.vcf --dbsnp dbsnp_138.ucsc.hg19.vcf -nt 8 mother.g.vcf --variant:VCF father.g.vcf --variant:VCF child.g.vcf

I will generate the bamout file shortly.

• torontoMember

If we upload the bamout files, are they private? Basically, there was no reads from mother and father.

The uploads go to a private directory that can only be accessed by our staff. We can't 100% guarantee we'll never get hacked, but we do make every reasonable effort to prevent data leaks. And we do delete user data when we're done with it, unless we are explicitly given permission to use some of it for systematic testing purposes.

• torontoMember

Ok. How to upload the files?

@zwang
Hi,

If you post the bam file screenshot to the forum, they are viewable by the public. Geraldine is referring the uploading to the FTP when we ask you to upload a bug report. For now, I am requesting a small screenshot of the bam file and bamout file of the region in question. I want to see why the call was changed from homozygous reference to heterozygous. To do that, I first need to see what the region looks like in the bam file. The GQ from HaplotypeCaller is 34 which is not very high, so I wanted to see what is going on in the region.

How do you know for sure the mother is homozygous reference?

Thanks,
Sheila

@zwang
Hi again,

To upload the files, you can click on "Attach a file" below the comment box. Then select the screenshots you would like to upload. We really prefer this way because it is easy for us to view, and the region is not easily identifiable. However, if you prefer to keep things private, you can submit the files to our FTP. Instructions are here. You don't have to submit all the files required in the instructions, just the relevant bam files.

-Sheila

• torontoMember

Here is the screen shot, as I mentioned before there is no reads from mother. And one thing, when I used the -bamout option, I did not use -nct 8; the program says it does not support that option with bamout.

• torontoMember

Here is the problematic site: before applying CalculateGenotypePosteriors the mother's genotype (the bold) is 0/0, after CalculateGenotypePosteriors it became 0/1

joint_recal.vcf:chr1 96925200 . A T 530.13 PASS AC=1;AF=0.167;AN=6;BaseQRankSum=0.599;ClippingRankSum=1.20;DP=87;FS=0.000;GQ_MEAN=486.00;MLEAC=1;MLEAF=0.167;MQ=60.00;MQ0=0;MQRankSum=-8.750e-01;NCC=0;QD=18.93;ReadPosRankSum=0.921;SOR=0.874;VQSLOD=15.95;culprit=MQ GT:AB:AD:DP:GQ:PL 0/1:0.460:13,15:28:99:561,0,486 0/0:.:25,0:25:34:0,34,969 0/0:.:31,0:31:68:0,68,1097

joint_recal_gposteriors.vcf:chr1 96925200 . A T 530.13 PASS AC=2;AF=0.333;AN=6;BaseQRankSum=0.599;ClippingRankSum=1.20;DP=87;FS=0.000;GQ_MEAN=486.00;MLEAC=1;MLEAF=0.167;MQ=60.00;MQ0=0;MQRankSum=-8.750e-01;NCC=0;QD=18.93;ReadPosRankSum=0.921;SOR=0.874;VQSLOD=15.95;culprit=MQ GT:AB:AD:DP:GQ:JL:JP:PL:PP 0/1:0.460:13,15:28:99:0:24:561,0,486:527,0,545 0/1:.:25,0:25:25:0:24:0,34,969:25,0,935 0/0:.:31,0:31:34:0:24:0,68,1097:0,34,1063

BTW, the last three columns are in the order of child, mother, and father.

#### Issue · Github March 2016 by Sheila

Issue Number
722
State
open
Last Updated
Assignee
Array
Milestone
Array
• torontoMember

Do I need to upload more data?

@zwang
Hi,

Sorry for the late response. I need to think about this/check with the team. I will get back to you soon.

-Sheila

• torontoMember

Any update on this issue?