Hello everyone, I have human Rna-Seq data reading by MiniOn.I aligned with the human reference genome (GRCh38.p13) from that website= https://www.gencodegenes.org/human/ using minimap2.
I download annotation from the same website. I used SubRead for quantification and feature count. On the SubRead results, there is a problem with one of the bam files. I aligned again that fastq files with the same reference genome of the first one and same gtf file and also I checked fastq file of second one but I didn't see any problem. Only differences are these bam file has a different barcode
Do you have any suggestions for that problem?
featureCounts -T 8 -a gencode.v38.chr_patch_hapl_scaff.annotation.gtf -g 'transcript_id' -o readcouts.txt bam/*.bam
|| Total alignments : 11214480 ||
|| Successfully assigned alignments : 4051945 (36.1%) ||
|| Running time : 2.67 minutes
|| Total alignments : 0 ||
|| Successfully assigned alignments : 0 ||
|| Running time : 2.89 minutes
Salmon Quant Alignment Based results
salmon quant --ont -t reference.fa -l A -a first.bam -o salmon_quant1
Total # of mapped reads : 5465357
of uniquely mapped reads : 328808350000000
ambiguously mapped reads : 2177274
salmon quant --ont -t reference.fa -l A -a second.bam -o salmon_quant2
Completed first pass through the alignment file.
Total # of mapped reads : 3843632
of uniquely mapped reads : 2552463
ambiguously mapped reads : 1291169
Not answering your question but have a suggestion. Take a look at https://github.com/nanoporetech/pipeline-transcriptome-de if you are analyzing RNAseq data on MinION.
Have you looked at your alignments using IGV or some other browser? Are reads aligning to right areas of the genome?
Actually, I am following that pipeline but I am also trying other tools. I looked at alignments using IGV, it is also true.
Are the reads here overlapping multiple genes? Or are the reads short and are multi-mapping and thus not counted? Does QC look similar for two samples?
It's similar QC reports and both bam files have the same headers. I tried to extract the region of my interest from the second one and tried again Subread feature count now its results 5%. But I need to count whole for differential expression.
Does anyone have a different idea?
Show code please, it is completely unclear which command run ran.
featureCounts -T 8 -a gencode.v38.chr_patch_hapl_scaff.annotation.gtf -g 'transcript_id' -o readcouts.txt bam/*.bam