Discussion about RNA-seq ref-based tutorial

[lldelisle]

Update from GTN We've converted this issue to a "Discussion" to support sharing multiple approaches and discussing/voting/commenting on each.

---

Dear authors of the Reference-based RNA-Seq data analysis tutorial,
Namely @bebatut @malloryfreeberg @MoHeydarian @erxleben @pavanvidem @blankclemens @mblue9 @nsoranzo @pvanheus.

With @hrhotz we will organize a galaxy training on RNA-seq (https://www.sib.swiss/training/course/20220427_GALXY) end of April.

I noticed that in the tutorial you estimate the strandness with infer experiment and then count the reads per genes with featureCounts.
I personally usually use STAR to count the reads per gene and then use multiQC to determine the strandness.

I was wondering what would be the best way to talk about this alternative.

• Should I add a paragraph after "Counting the number of reads per annotated gene"?
• Or could we imagine a way to do a sort of 'alternative' so users can choose the 'STAR count option' or 'featureCount option' at the beginning of the tutorial and then the html changes at the places where the 2 options use different tools/different options. I don't know if it is feasible but something like html tag which would say if a section is 'STAR count' specific or 'featureCount' specific... Maybe @hexylena or

@shiltemann can tell me if this is technically impossible... I imagine that as the number of contributor increases this situation may happen in other tutorials...
Thanks for your input...

PS: If you're not interested in this topic but you received notification, I am sorry. You can unsuscribe by clicking at the end of the email notification or on the github on the right menu.

[hexylena]

Fancy Javascript + New GTN Markup

I've had mockups for Choose Your Own Adventure type tutorials in the past! #1120 The code was kinda ugly to write as a tutorial author, but the result was quite nice for students. There you had to force content into blocks (either manually using divs or blockquotes/tip boxes/etc) and then tag those with the path name, e.g. {: .openstack} and it would be shown or hidden if the User selected Path matched that

[shiltemann]

Multiple Hands On Boxes

Hi @lldelisle!

I think it would probably work to just have two hands-on boxes in that section with the two different approaches? If you label them "Option 1: Hands-on: ..." and option 2, it should be clear? you could even put them in collapsible boxes to make it even clearer if you wanted

More complex solutions like you propose are possible, but maybe not necessary? That said with the 16S tutorial we do somthing similar with long or short version, some shared bits, and some different sections depending on which version you are viewing. A bit more work to maintain however, but a possibility.

[malloryfreeberg]

Collapsible Boxes

👍 to the idea of having options! I think this is a great way to encourage people to explore other approaches and figure out which is the best for their data. Even better if there is some information in the tutorial about why one might choose one option over the other.

I like the idea of putting alternatives in collapsible boxes. These boxes could also mention that if this alternative option is chosen, some of the downstream results might be slightly different but that it's ok! Hopefully the overall conclusions are similar :)

[hexylena]

Micro Tutorials

We could also explore micro tutorials where you use the follow_up_trainings to lead to option 1/2/3 at a given choice point, before having them go back to using the same end-half tutorial.

[pavanvidem]

Hi @lldelisle,
I wonder why do you have to determine the strandness if you already have counts from STAR.

SureSTAR count is a good alternative to featureCounts. IMO, there are pros and cons with both these approaches. STAR count simplifies 3 steps of mapping, estimating strandness and counting to a single STAR run step. This is a huge improvement. But there are the following problems with this approach:

• There will be no flexibility of parameter turning for counting. It can be either using a different MAPQ or counting only reads instead of fragments for a bad paired-end sample or most importantly, STAR does not allow choosing different feature types and gene identifiers for non-standard GTF files.
• We also need a feature lengths file (we use featureCounts to generate this), which is later needed for goseq. We have to write an additional step to generate a length file.
• I normally do not use STAR count feature for bulk RNA-Seq data. Do you get more or less the same number of counts as in featureCounts?

I am not against having an alternative, but we have to tackle the above-mentioned problems and write a proper description of when to choose STAR count alternative. I like the idea of having different tool based workflows in a single tutorial 👍

[lldelisle]

• For the MAPQ etc this is totally true, you loose flexibility.
• For the length, when I use goseq (locally) I use TxDb.Mmusculus.UCSC.mm10.ensGene (which is probably less precise) so I don't know how this could be implemented in the goseq tool... I need to think about it.
• For the comparison, I will check and tell you.

[lldelisle]

Hi,
Few updates here.

• I did in parallel STAR counting and featureCounts. As you can see in this plot, there is no difference for SR but there is for PE:
  
  This is because featureCounts is clever than STAR (see --quantMode GeneCounts for paired-end alexdobin/STAR#1505). I think the difference is bigger with unstranded than stranded as it is linked to mates mapping to 2 genes at a time.
• For the length, there is a tool in galaxy called "get_length_and_gc_content". The current version failed because I realized that the gtf is an ensembl GTF with "chr" in front of chromosome names. While it works for the main chromosomes, it is a big mess for others. Anyway, I updated the tool so it will work with the current GTF. The new version should be on usegalaxy.eu next week so I will test it.

[lldelisle]

Hi there,
I just noticed that I did not have the same results when using the new freshly computed counts following the tutorial with featureCounts and the counts provided on zenodo:

The one where it deviates are the PE sequenced samples where the counts I optained are roughly half what is on zenodo. My suspicion is that the one on zenodo were 'read counts' while on the tutorial you propose to use 'fragment counts'.

• Should I update the zenodo to match the tutorial?
  • Then the number of DE change and the GO categories also.
  • There is no more under represented KEGG pathway and for the over-represented there are 4 if you use featureCounts and 2 if you use STAR.

Another point, I think for RNAseq analysis, when you have multiple samples, the use of collection is for me the best use case (For me the collections are the best recent improvement of galaxy). Do you have any objection to a switch of a 'collection-based' tutorial for the first part (reads to counts)?

[hrhotz]

WRT the first point (which I can confirm): as alternatives to discuss:

• change the first part of the tutorial, i.e. change the option "Count fragments instead of reads" in the two featureCounts steps to "Disabled". This way, we don't have to up-date all the downstream stuff.
• just add a a comment to the tutorial, like: "The featureCounts table you download from zenodo (or from the Data Library) were generated with the option "Count fragments instead of reads" set to "Disabled". This way you roughly get twice the number of counts as each each read (from a pair) is counted individually.

WRT the other point: I support this suggestion.

[nsoranzo]

I'd vote for HRH's option 2 (i.e. adding a a comment to the tutorial).

[lldelisle]

The PR is available at #3330 if anyone wants to review...