Hi, I am using Rsubread to map Illumina RNA-seq reads, to investigate alternative splicing events. Using: 1000 Synthetic reads (size of reads are 210 bp) and a reference chromosome concatenated from all the reads and 7000 more.

My aim is to use subjunc and learn about the exon-exon junctions in BED files. I have noticed the percent of mapped reads in the align function is much higher than in the subjunc function (align is 70%, subjunc is 40%).

Why does subjunc map less reads than the align function? As I read on the manual, the align is much more flexible, but I still can't understand what is there in the subjunc function that doesn't map as much as the align function, and how can increase the percentage?

I tried:

1. Playing around with the parameters (indels, mismatch, and nbestlocations)
2. Merging paired end reads (using PEAR)

Unfortunately the percentage of mapped reads is still the same.

Any suggestions?

align.stat <- align(index="/Users/files_intron_rent_chrs/INDEX/intron_retention_reference_chromosomes",
                    readfile1=mergedreads,
                    output_format = "SAM",
                    output_file = paste("alignedmerged", "SAM", sep="."),
                    phredOffset = 33,
                    indels = 1,
                    maxMismatches = 1,
                    detectSV = TRUE,
                    keepReadOrder = FALSE,
                    minFragLength = 20,
                    maxFragLength = 500,
                    unique = FALSE,
                    nBestLocations = 10)

## Align function Results ##
# Total reads : 932 
# Mapped : 686 (73.6%)  
# Uniquely mapped : 676 
# Multi-mapping : 10 
# Unmapped : 246 
# Fusions : 128
# Indels : 166

----------


ex_junctions <- subjunc(index ="/Users/files_intron_rent_chrs/INDEX/intron_retention_reference_chromosomes",
                     readfile1=mergedreads,
                     output_format = "SAM",
                     output_file = paste("alignedmerged", "SAM", sep="."),
                     phredOffset = 33,
                     indels = 1,
                     maxMismatches = 1,
                     reportAllJunctions = TRUE,
                     minFragLength = 20,
                     maxFragLength = 700,
                     nBestLocations = 16)


## subjunc function Results ##
# Total reads : 932  
# Mapped : 381 (40.9%)  
# Uniquely mapped : 374  
# Multi-mapping : 7 
# Unmapped : 551 
# Junctions : 131
# Fusions : 0  
# Indels : 89  

sessionInfo( )

Thanks in advance.