I am aligning RNA-seq data to a custom reference using
align.res <- subjunc(index="Repvec", file1=fwdname, readfile2=revname, output_file=bamname, nthreads = 6, sortReadsByCoordinates=TRUE)
With the last parameter, the BAM file is sorted and the .bai file is created.
Everthing is fine, when I load the BAM file into an genome viewer; however, when I want to do a quick check of the alignment with
samtools idxstat, the 'mapped read-segments' and 'unmapped read-segments' columns are empty:
samtools idxstats sample1.bam # mCherry_1495 1495 0 0 # GFP__1468 1468 0 0 # BFP2_1494 1494 0 0 # AMPR_2739 2739 0 0 # * 0 0 0
After recreating the index from the BAM file with
samtools index, the columns are populated:
samtools index sample1.bam samtools idxstats sample1.bam # mCherry_1495 1495 40412 0 # GFP__1468 1468 525 0 # BFP2_1494 1494 23843 0 # AMPR_2739 2739 3 0 # * 0 0 40681972
It would be helpful, if
align()) would add also the information required by
samtools idxstats to the .bai file.
The environment is as follows:
R version 4.3.0 (2023-04-21 ucrt) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 11 x64 (build 22621) attached base packages:  stats graphics grDevices utils datasets methods base other attached packages:  Rsubread_2.14.2 loaded via a namespace (and not attached):  compiler_4.3.0 Matrix_1.5-4.1 tools_4.3.0 grid_4.3.0 lattice_0.21-8