Hello, I am trying to use 12 fastq file of tomato chloroplast against a reference genome. I am trying to use the samtools mpileup command to compare between 12 samples and use bcftools call to get only the variants among the samples. When I am trying to load the vcf files into the IGV browser, I could only visualize the homozygous variants but not herozyguous regions. Can anyone suggest me what I could do to visualize both homozygous and heterozygous variants on IGV?
Here is my script:
Ref=NC_007898.3 bowtie2-build $Ref.fasta $Ref
for f in *.fastq; do
echo "File -> $f"
b=basename -s .fastq $f
echo "Working on -> $b"
bowtie2 -x $Ref -U $b.fastq | samtools sort - -o $b.bam
samtools mpileup -v -u -f NC_007898.3.fasta $b.bam |bcftools call -v -c --ploidy 1 - > $b.vcf