Entering edit mode
12.5 years ago
Lee Katz
★
3.2k
Hi all, I couldn't find if this has been asked here. I realize that samtools assumes that you are working with a diploid and will output heterozygous snps. How could I modify my script to make it not output heterozygous snps? In other words, how do I trick it into thinking that it is looking at a haploid genome? My snippet of code is in Perl and uses a command() subroutine which runs a string and prints the command to the screen.
command("samtools mpileup -uf $$settings{assembly} $bam > $out1");
command("bcftools view -cg - < $out1 > $out2");
command("vcfutils.pl vcf2fq -d $minDepth -D $maxDepth < $out2 > $fastqOutNonstandard");
command("vcfutils.pl varFilter -d $minDepth -D $maxDepth < $out2 > $variantsFile");
Is there perhaps something like GATK which could do it instead?