[GenoMax]

Ok. But the featureCounts command is the same as given before? or it should be changed?

Same command. You should use the rRNA GTF file instead of the original one. Adjust output file names if you need to keep those results around.

[bvk]

Same command. You should use the rRNA GTF file instead of the original one. Adjust output file names if you need to keep those results around.

OK. I used the same command. please have a look and tell me.

featureCounts -p -t exon -g gene_id -a hg19_rRNA -o MCF10A_BaP_pool1_counts.txt MCF10A_BaP_pool1_R1_postTrimming/accepted_hits.bam

//========================== featureCounts setting ===========================\\
|| ||
|| Input files : 1 BAM file ||
|| P MCF10A_BaP_pool1_R1_postTrimming/accepted_ ... ||
|| ||
|| Output file : MCF10A_BaP_pool1_counts.txt ||
|| Annotations : hg19_rRNA (GTF) ||
|| ||
|| Threads : 1 ||
|| Level : meta-feature level ||
|| Paired-end : yes ||
|| Strand specific : no ||
|| Multimapping reads : not counted ||
|| Multi-overlapping reads : not counted ||
|| Read orientations : fr ||
|| ||
|| Chimeric reads : counted ||
|| Both ends mapped : not required ||
|| ||
\\===================== http://subread.sourceforge.net/ ======================//

//================================= Running ==================================\\
|| ||
|| Load annotation file hg19_rRNA ... ||
|| Features : 1769 ||
|| Meta-features : 3 ||
|| Chromosomes/contigs : 38 ||
|| ||
|| Process BAM file MCF10A_BaP_pool1_R1_postTrimming/accepted_hits.bam... ||
|| Paired-end reads are included. ||
|| Assign fragments (read pairs) to features... ||
|| ||
|| WARNING: reads from the same pair were found not adjacent to each ||
|| other in the input (due to read sorting by location or ||
|| reporting of multi-mapping read pairs). ||
|| ||
|| Read re-ordering is performed. ||
|| ||
|| Total fragments : 56783755 ||
|| Successfully assigned fragments : 246578 (0.4%) ||
|| Running time : 3.98 minutes ||
|| ||
|| Read assignment finished. ||
|| ||
\\===================== http://subread.sourceforge.net/ ======================//

Thank you

[GenoMax]

Provided the input files were correct this is the relevant bit

Successfully assigned fragments : 246578 (0.4%)

That would mean you have 246578 reads (each from R1/R2 files) that are aligning to rRNA (if that is the only thing in your GTF file). So the ribo-depletion appears to have worked reasonably well. Does not help you if you actually were interested in rRNA.

There is also this warning to keep in mind

|| WARNING: reads from the same pair were found not adjacent to each ||
|| other in the input (due to read sorting by location or ||
|| reporting of multi-mapping read pairs). ||

[bvk]

Provided the input files were correct this is the relevant bit



That would mean you have 246578 reads (each from R1/R2 files) that are aligning to rRNA (if that is the only thing in your GTF file). So the ribo-depletion appears to have worked reasonably well. Does not help you if you actually were interested in rRNA.

There is also this warning to keep in mind

|| WARNING: reads from the same pair were found not adjacent to each ||
|| other in the input (due to read sorting by location or ||
|| reporting of multi-mapping read pairs). ||

what do you mean here (Does not help you if you actually were interested in rRNA.) I didn't get this. There are 246578 read pairs aligned to rRNA regions. Could you please tell clearly about the warning. And finally is this right or not? If not how can I do now?

[GenoMax]

Ribodepletion has worked well (only 0.4% reads align to rRNA).

Fragments that go into a library have a size distribution and generally it is possible to infer the size of the fragment based on the alignments to a reference. The warning indicates that in some cases the reads were found to be not at the expected distance from each other. This is to be expected since you are looking at a repeat region and the annotation is not perfect.

You have said you (or the person who did the experiment, if you are only analyzing the data) are interested in rRNA biology yet the samples are rRNA depleted. So those two observations are contradictory. Aside from this we don't have any information about what this experiment was for (e.g. were you only interested in finding expression levels of genes but not differential expression or are there groups of samples that need to be analyzed for differential expression) so only you have an idea as to what is to be done next.

[bvk]

Ribodepletion has worked well (only 0.4% reads align to rRNA).

Fragments that go into a library have a size distribution and generally it is possible to infer the size of the fragment based on the alignments to a reference. The warning indicates that in some cases the reads were found to be not at the expected distance from each other. This is to be expected since you are looking at a repeat region and the annotation is not perfect.

You have said you (or the person who did the experiment, if you are only analyzing the data) are interested in rRNA biology yet the samples are rRNA depleted. So those two observations are contradictory. Aside from this we don't have any information about what this experiment was for (e.g. were you only interested in finding expression levels of genes but not differential expression or are there groups of samples that need to be analyzed for differential expression) so only you have an idea as to what is to be done next.

Differential expression is already done. But the clients need the total information in a report. In that they also need the read counts aligned to different regions. I will ask the person again whether the samples are ribo depleted or not.

And if I need reads that are aligned to intron, exon or some other regions it should be the same arguments I need to use but with different gtf file which has exons and introns. Is it right?

[GenoMax]

And if I need reads that are aligned to intron, exon or some other regions it should be the same arguments I need to use but with different gtf file which has exons and introns. Is it right?

You have already done that (before you used the rRNA GTF file). If you need read summarization at feature level you will need to look at the featureCount options again (one of which would be removal of the -g gene_id option).

I hope you are making an effort to understand the finer points of the programs you are running (specially since you are going to hand the results off to someone else). Getting a program to produce an output does not always mean that the output is correct or logical. This is where an analyst like you has to use their expertise/judgement.

[bvk]

You have already done that (before you used the rRNA GTF file). If you need read summarization at feature level you will need to look at the featureCount options again (one of which would be removal of the -g gene_id option).

I hope you are making an effort to understand the finer points of the programs you are running (specially since you are going to hand the results off to someone else). Getting a program to produce an output does not always mean that the output is correct or logical. This is where an analyst like you has to use their expertise/judgement.

Thanks a lot. Actually it is only to check whether rRNA depletion is done well or not.