I have a plamid DNA sequenced using MiSeq PE250 with read count at 1.3 million reads. I am using different softwares to align these data. My final aim is to analyse different resequencing methodologies and validate them for future use.

I first trimmed the raw data with Trim-Galore (One level trimming) and then with Sickle (Two level trimming). When I align one level trimmed data with the reference and generate a coverage file using samtools mpileup, the coverage at position 1 is 0.

When I use two level trimmed data to align and generate coverage the coverage for position 1 in this case is 32X.

I don't understand how can the coverage for position 1 be zero with less trimming (only trim galore) and 32X with more trimmming (Trim Galore + Sickle). Following are the command which I used.

Trim Galore

trim_galore -q 30 --length 50 --no_report_file --paired read1 read2

Sickle

sickle pe --gzip-output -q 30 -l 50 -t sanger -f read1.trimmed -r read2.trimmed -o read1.sickle -p read2.sickle -s single.sickle

Alignment using Bowtie2:

bowtie2 --local -x SampleName -1 read1.sickle -2 read2.sickle -I 0 -X 1000 --fr -p 16  -S SampleName.sam

Then i sorted and indexed SAM file using Samtools sort and index command.

Then I performed BAQ recalibration using calmd . Indexed the resultant BAM file again. Used mpileup command to generate coverage for the entire sequence.

samtools mpileup -d 10000000 -aa -f ref.fasta SampleName.sorted.bam > SampleName.mpileup

cut -f1-4 SampleName.mpileup > SampleName.depth

Can anyone please explain why I received 0 coverage at position 1 with less trimming?