Setup RNA-Seq Pipeline from scratch: fastq (reads) to counts | Step-by-Step Tutorial

This is a detailed workflow tutorial of how to process bulk RNA-Seq reads (fastq) and generate counts matrix which can be used for various downstream analysis. In this video, I walk through setting up a pipeline in bash (linux) and perform steps to process RNA-Seq data including -
• Quality control (fastQC)
• Trimming (Trimmomatic)
• Alignment (HISAT2)
• Quantification (featureCounts)
In addition I also talk about run times, memory requirements and aligner accuracies of various aligners. I hope you find this video helpful! Leave your thoughts in the comment section below!
Link to code:
github.com/kpatel427/KZbinT...
Link to data:
▸ drive.google.com/file/d/1DGHj...
Linux Basics
▸ ubuntu.com/tutorials/command-...
▸ xie186.github.io/Novice2Exper...
▸ hackr.io/blog/basic-linux-com...
To Trim or to not Trim?
▸ www.ncbi.nlm.nih.gov/pmc/arti...
Strandedness
▸ bmcgenomics.biomedcentral.com...
▸ chipster.csc.fi/manual/librar...
▸ rseqc.sourceforge.net/#infer-e...
Chapters:
0:00 Intro
1:11 - Applications of RNA-Seq data
2:41 Schematic detailed workflow
3:49 What are splice-aware aligners?
5:59 Workflow for this tutorial
6:31 Comparison of run times, memory usage and aligner accuracies for various aligners
8:27 Which aligner should I choose?
9:05 Pre-requistes to build this pipeline (things that will not be covered in this video)
9:40 Set-up before building the pipeline
10:37 Some good practices while building a pipeline
11:21 Quality control: FastQC
15:27 To trim or to not trim?
15:58 Trimming reads: Trimmomatic
19:18 Align reads: HISAT2
24:20 Read quantification: featureCounts
Show your support and encouragement by buying me a coffee:
www.buymeacoffee.com/bioinfor...
To get in touch:
Website: bioinformagician.org/
Github: github.com/kpatel427
Email: [e-mail]
#bioinformagician #bioinformatics #rnaseq #omics #rnastar #ensembl #affymetrix #microarray #affy #ncbi #genomics #beginners #tutorial #howto #research #biology #GEO #ngs

[@Luvinlife411]

Thank you, thank you for these wonderful videos! Perfect for beginners and they've helped immensely for completing my thesis!

[@chienlee4695]

I really like your videos. This tutorial help me go over the bumpy learning curve! I am so excited that I can get the read count by myself. Thanks for putting your time to make these videos.

[@sivanim3830]

Great explanation for a beginner to get started with RNA seq! Cheers to you!

[@josefinagarciadiaz478]

I just want to say thank you very much for doing this, im just starting out in my PhD studies and you have been really helpful!

[@CHENYanZhen]

Wow, you really illustrate the workflow so systematically. As a beginner, I only knew a little about RNA-seq analysis, so this helps a lot! Besides, you mention and explain lots of practical details, unlike some other video providers. I gain so much from your videos! Many thanks to you!!!

[@aamirmalik7740]

Awesome, thank you for doing this video [much awaited] -- so nice of you Khushbu -- Please keep making such great stuff -- Many thanks

[@apedike]

This is so clear--love your content

[@codewithme_1988]

You r doing better than other. Keep it up

[@aspenkayceedowning6588]

Thank you so much for this. I am in my honors and I am the only one on my team that can do the bioinformatics and my classes haven't caught up yet. Thank you so much for this tutorial. I am liking, subscribing and commenting to show support so you can keep doing this

[@kingsleyibeabuchi3461]

Great job. Your vidoes have been very helpful. I am sincerely grateful

[@kitdordkhar4964]

KZbin is a library. You are one of the books that I frequently read. Thanks! :-)

[@Bioinformagician]

Flattered and grateful :)

[@wolveshowlhigh]

Thank you!!! You've been a lifesaver for me as a newbie to the lab. Would love to see the ATAC-seq version of this !

[@hozifaelgadal623]

thank you , this was such a help , I really loved your videos

[@felipevilicich980]

Brilliant! Thank you very much for the tutorial!

[@alvaroruiztabas5627]

Amazing helpful video. Congratulations!!!

[@christianahoyindamolajemiy5531]

Thank you so much. Very understandable and easy to follow.

[@kingsleyibeabuchi3461]

Thank you very much for this wonderful video. Keep it up. I was also able to perform most of these tasks in Galaxy

[@shrishteekandoi120]

You’re the best !

[@histephenson007]

Brilliantly done!

[@kiplimosimon1429]

Thanks for the nice and informative session.

[@juniorkika254]

Very helpful video.Thank you!

[@Ice84letters]

Thank you so much for this videos i am studying bioinformatics but we made little practice so i have many empty knowledge in some parts of this pipelines...your videos are helping me so much!! 💜

[@Gbemi78]

You are such a blessing.

[@edepoisdoenem8634]

Nice vídeo. I use R in Windows-based systems but your video was useful to understand the steps and the logic of this analysis.

[@ignaciomazzitelli4182]

Brilliant!

[@ashishbhaladhare]

thank you, its really informative for me.

[@tushardhyani3931]

Thank you for this video !!

[@ritikabassi1998]

It was a wonderful video.

[@arpanparichha6784]

This is really helpful

[@rashmitan6867]

You are doing such a great job ! Kudos! Would be willing to show the analysis of 10x spatial transcriptomics in the future ?

[@rushupatel]

Good job!!!

[@escastorage7427]

cannot be happier

[@akshayaupadhyay8063]

Hi, great videos! You mentioned that you would be detailing on how to install tools(HISAT etc) for the pre-processing steps. Could you elaborate on that, probably if there is a script for running on Mac terminal/Shell then it would be great! Thank you!

[@sdmasroor]

Thank you so much!

[@tarahsify]

You're awsomee❤❤

[@nk12987]

Thank you so much for this video it's very helpful. Can you please make a video of further steps like count matrix and finding deg's .

[@nayanbhowmik5010]

Thanks for your beautiful videos! Can you make a similar video for PE reads and include the step of generating genome index files using HISAT2? That would be so helpful.

[@kajalpanchal8239]

Thank You !

[@yaarobaltali6076]

Thank you so much for this helpful videos please keep making such a great content and I have one question, when are you going to upload a video covering FastQC analysis

[@Bioinformagician]

I shall make a video on FastQC soon, once I finish the videos I have currently lined up. Thanks!

[@ishakyusuf4616]

Thanks for your video. It is amazing.

[@ishakyusuf4616]

I wander if the purpose of --rna-strandness R and -S 2 are to accelerate the hisat2 and featurecounts works?

[@Bioinformagician]

--rna-strandedness provides information on orientation of sequencing reads. This video explains this concept really well: kzbin.info/www/bejne/lZ-uiIaOiahsl5I

[@ishakyusuf4616]

@@Bioinformagician I see. thanks🤞

[@renqiuguoli1258]

Good show, thx a lot

[@rahultanpure1338]

Great informational video. Thanks a ton.🙏 If possible please make a video on de novo Trinity pipeline.

[@Bioinformagician]

Sure, will definitely plan on creating a video covering de novo transcriptomic assembly.

[@rafaelamachadotugores4769]

thank you, girl

[@freezingtolerance7493]

Thank you for your video. it was very helpful to learn about rna-seq workflow... just quick questions: after hisat2 (alignment) is it essential to do variant calling stage before read counts? or I can just do readcounts using featurecount or htseq-count without process of SNP.

[@bastiwala0808]

Hello there. It is a very informative video kindly make a tutorial video of miRNA processing from fastq to readcount. Thanks

[@HopefulHopelessly]

Hi, wonderful video so far (I'm only a little over halfway through) I was wondering why you chose to trim it by 10 base pairs? Was that just an arbitrary decision or is it recommended? Thank you!!

[@acramulhaquekabir5852]

Great video. I was wondering when are you planning to show the next steps. From count matrix to DEA

[@Bioinformagician]

I have covered those steps in detail in two of my videos - 1. kzbin.info/www/bejne/ZpOVZaCmr7KSa68 2. kzbin.info/www/bejne/hauxq4JuptyqppI

[@user-uq3qh2cy9v]

Hi, this video is very concise. Thanks for sharing. I'm still struggling install HISAT2 onto M1 mac. How did you do it?

[@recordinghoney]

Hi! Great video. I was wondering if there is a similar quality check for CEL files. Do you happen to know?

[@debashruti4984]

Hi! Wonderful video..makes things so easy to understand. I had one question though. After alignment when do we use Cufflinks and when to use featureCounts? What differences are there in the two tools? First one being cufflinks-cuffmerge-cuffdiff and the latter being featureCounts-Deseq/EdgeR

[@louisepitcher6180]

Amazing easy to follow video! I was Just wondering if you could clarify how the lines of code would change for paired ends reads, specifically when using the TruSeq Stranded mRNA kit for bulk RNA-seq library prep. Thanks!!

[@Bioinformagician]

I am glad you liked my video! ☺️ Starting at the trimming step, Trimmomatic will be required to run on "PE" mode, followed by HISAT (change -U parameter for paired end reads). There is a strandedness parameter in both - HISAT (--rna-strandness parameter) and featurecounts (-S parameter). You can test strandedness of your protocol by running infer_experiment.py script from RSeQC package (rseqc.sourceforge.net/#infer-experiment-py).

[@truehighs7845]

Love it, thanks!!! I have a question, indeed not only RNA will yield a sequence that can be difficult to reconnect the exome itself because expression of a given mRNA, the problem that persists is correlating what expression for what mRNA, otherwise whatever counting operation you make, it will not distinguish between various RNA outputs, it will only be a composite signal of all the RNA present in your sample - even if you drill down to a single cell assay, you still have a composite signal hard to get any significance from?

[@sofiapuvogel3663]

Thank you very much for this video! I have a question for paired end data. After doing samtools sort, the reverse read and the forward read are not consequetive reads in the sorted bam anymore. Is that okey?

[@user-vk8bd1re8c]

❤❤❤

[@adeelabbasi7375]

please make a video on De Novo RNA seq with paired end reads - its a hot topic

[@user-ns7ng2hz9w]

Great and clear… Thanks. I got two questions though: 1. what’s the difference between the count matrix and TPM/FPKB/RPKB? Noting that TPM/FPKB/RPKB can’t be use for DEG, while count matrix can. 2. Why Does bowtie can’t be used for alignment to genome?

[@maytelopez-cascales6113]

Hi, great tutorial, i have a question. I have aligned with star, bulk RNAseq from voxels clones, and i have 52% uniquely mapped reads, 30% of multimapped, and the rest of not aligned. Can you explain different aspects of the output of star, if they have good quality or not? with FASTQC i cant find any adapter.

[@subarnathakur6973]

Please make a video of setting up a pipeline using Trinity for denovo sequence assembly

[@parnakanodia2262]

Hi thank you for your fabulous lecture, could you please post a lecture (on urgent basis if possible) on how to remove duplicates from fastq files? Using oxford academic bioinformatics' minirmd?

[@kruthiiirao]

❤

[@anishanna5125]

Thank you for the very useful tutorial. Could you please advise how to exclude the multi-mapped reads that were mapped by HISAT2? I knocked out a gene (and confirmed the KO by qPCR and western blot before sequencing) but the KO gene was not differentially expressed in my RNAseq. The only explanation that I could think of it that the KO gene sequence is very close to multiple other genes, thus multi-mapped sequences from other similar genes could masked the gene deletion. Thank you!

[@Bioinformagician]

Did you quantify your reads using featureCounts or any quantification tool that excludes multi-overlap reads? Also, can you check how many reads map uniquely to the gene in your BAM files (in both WT and KO) using sam flags?

[@ynpan910]

Hi I have watched couple of your videos and they are so informative! For this pipeline, can I ask if there are R packages that can do the same thing? bc I am not familiar with linux or bash scripts but I will be learning them. Thanks for all the efforts you put in the videos.

[@chrisdoan3210]

Hi Bioinformatician. Would you please tell me why we have many pipelines for RNA-Seq? So these steps below are not considered pipeline (workflow): 1. Check quality with FastQC 2. Trim reads with Trimmomatic 3. Align reads to the reference with STAR 4. Calculate gene hit counts with FeatureCounts 5. Compare hit counts between groups with DESeq2

[@amitrupani9898]

Thanks for this great tutorial! Question - If you have high duplicates rate (say 70-80% in case of FFPE) pre alignment, would you remove them? There are some wet lab QC metrics like DV200, RIN, etc that also determines the effectiveness of downstream analysis might be worth considering/including for your tutorial on QC. Keep it up and Cheers!

[@Bioinformagician]

A general consensus about dealing with duplicate reads is to NOT remove them. The level of duplication could be technical (FFPE for example) or could be due to higher coverage or expression levels of certain genes. My concern with duplicate reads are 1. how do we distinguish reads based on source of duplication and 2. In the process of removing those reads, we could risk losing out on true expression levels. Thank you for bringing in the wet lab aspects that determine RNA quality to my notice. Definitely very important, will surely talk about it in my RNA-Seq QC video. Thanks again!

[@amitrupani9898]

@@Bioinformagician Thanks for sharing your thoughts. Looking forward to upcoming content. Cheers!

[@francescosilvestro2092]

Thanks a lot for this tutorial, it helps me a lot to start with. What's your math background? Have you studied advanced statistic or something related to it?

[@serrasonmez7912]

hello, thank you for your videos, and can you please make video about the github, so uploading and presenting the projects on github?

[@sMr_Borgov]

Hi, quick question, the count table display gene ID's and their counts, how do I get more gene information? (such as gene name, symbol, and annotation added ). Thanks

[@gcbicca]

I run another dataset of rnaseq and the results of fastqc for this dataset is the same but with Sequence Duplication levels with a flag, this is significant in experiment of Rnaseq that analyze control vs treatment?

[@CostanzoPadovano]

Thank you for this tutorial! I hope that you will create a patreon to support your channel :)

[@Bioinformagician]

Working on creating one. Until then, you could support me by buying me a coffee - www.buymeacoffee.com/bioinformagic

[@grsbiosciences]

Hi madam, how to convert raw data.gz file to fpkm file. You have used in your other video fpkn file and proceseed later using deseq2 library. There you have taken directly tab file. Buti have fastp.no.tab file. How i should convert yhis to fpkm tab file

[@user-mn4gu7kv2d]

Hi, thanks for the video! But not all the commands are working. I am using the terminal in macOS and changed to the bash shell. What's the problem and how to fix it?

[@sutrishakundu2147]

I don't have the genome sequencing data...I need to do de novo assembly of my transcriptome data...can you please elaborate how can I do this? What additional steps I need to follow before hisat2?

[@Tekofilic]

Hi can someone guide me about the error: Failed to open the annotation file ../celegan/caenorhabditis.gtf, or it’s format is incorrect or it contains no ‘exon’ features

[@tolga1292]

When will you release the video about post alignment qc? :)

[@nourarifi2642]

thank you for your video. I have a question : Can we use exactly the same pipeline to analyse single cell RNA sequencing data or it's not the same ?? thank u in advance

[@mostafakarami7687]

I have the same question! I'm curious what would be different from bulk sequencing.

[@shreyaslabhsetwar9055]

Hey, I wish to quantify human micro RNAs from double stranded fastq files. Which reference shall I use?

[@vargas4762]

Thanks for this tutorial! Could u help me on how can i create an index based on other genome like Aedes aegypti. i want to make RNA-seq analysis of infected mosquitoes, wich file should i use as reference genome ?

[@Bioinformagician]

If a reference genome is not available, you could assemble the genome first using the sequence assembly from NCBI. Once assembled, you could use it as a reference genome.

[@anamikapandey4769]

first of all i would like to thank you for your youtube videos which are very interesting and understandable.please keep posting more videos like this. i have one issue while i am preparing my counts table for chip seq analysis. as when i my giving the last command of features counts then its gives error all the time and the error is "failed to open the annotation file , format is incorrect or it contains no exon features" . i am just brainstorming with it from past couple of months. i will be heartily obliged to you if you could help me regarding this .

[@shivanikushwaha3839]

I get the same issue. Did you find a solution?

[@linkandash8025]

@Bioinformagician can you do an Orthofinder tutorial ?

[@deepanshugarg4257]

How to use featureCounts command for paired end reads? and hisat2 as well?

[@mohdfahad6224]

if i have replicates, at which step i will have to merge the two rna seq data?

[@sonalidey1454]

Maam why raw sequence number and after adapter removal sequence number varies ... as all raw read have adapter

[@quinattasneemrafique536]

The link to the raw reads is not attached in the description. Kindly provide please

[@RaviKumar-cb7tw]

Can you please tell me how to install Hisat2 in my computer? I am using Macbook AIr M2 chip and Sonoma 14.5 OS

[@stanyang4321]

Are you using R for this tutorial ?

[@anxu4239]

Could you please make a video of how to trans counts data into FPKM or RPKM or TPM?

[@Bioinformagician]

You could use tools like RSEM for read quantification which gives you FPKM and TPM values along with expected counts.

[@zulfakhan2655]

can we convert fasta instead of fastq to count data?

[@chrisdoan3210]

Would you please guide me on how to set up the tools in this tutorial for a Mac M1 user because I got a lot of difficulties and errors for about a month? Thank you so much!

[@Bioinformagician]

I set this pipeline up in Mac M1 pro, you should be able to reproduce this without major difficulties. Can you share what errors/difficulties you faced?

[@chrisdoan3210]

@@Bioinformagician would you please tell me how to update bash on an M1 mac? I tried to do it and asked the Internet but still couldn't figure out how. I got this error: Error: Cannot install in Homebrew on ARM processor in Intel default prefix (/usr/local)! when I try: brew install bash. The default bash version on Mac is old. Thank you so much!

[@chrisdoan3210]

@@Bioinformagician Would you please tell me how to set up bash on Mac M1? The default bash is obsolete and I got error when trying to upgrade bash.

[@Bioinformagician]

@@chrisdoan3210 You are not required to update bash. It should work fine on an M1 Mac. What error do you see when you try to run this pipeline in your bash?

[@chrisdoan3210]

@@Bioinformagician I don't have an error now but I think it is better to use bash in the future than zsh.

[@JaskaranSingh-om5mv]

how to download the tools required for the pipeline on a M2 MacBook?

[@chrisdoan3210]

Hi Bioinformatician. I use STAR instead of hisat2 so what command I should use to run your data?

[@Bioinformagician]

Did you try following the steps mentioned in STAR’s manual?

[@chrisdoan3210]

@@Bioinformagician I found this: STAR --runThreadN 6 \ --runMode genomeGenerate \ --genomeDir chr1_hg38_index \ --genomeFastaFiles /reference_data_ensembl38/Homo_sapiens.GRCh38.dna.chromosome.1.fa \ --sjdbGTFfile /reference_data_ensembl38/Homo_sapiens.GRCh38.92.gtf \ --sjdbOverhang 99 But I got difficulty replacing the attribute of STAR with your corresponding input from hisat2.

[@chrisdoan3210]

Thank you for your videos! When I run this code: hisat2 -q --rna-strandness R -x HISAT2/grch38/genome -U data/demo_trimmed.fastq. I got error message: hisat2: command not found. Would you tell me how to fix this? Does the way we install hitsat2 matter? What is location of hisat2 should be?

[@Bioinformagician]

Error is telling you it cannot find hisat2 executable. You need to provide the complete path to hisat2 (path to where hisat2 is downloaded).

[@chrisdoan3210]

Hi@@Bioinformagician, I tried to follow your suggestion: /Users/Chris/Bash_script/HISAT2/hisat2 -q --rna-strandness R -x HISAT2/grch38/genome -U data/demo_trimmed.fastq | samtools sort -o HISAT2/demo_trimmed.bam. I got this message: dyld[1437]: Library not loaded: @rpath/libcrypto.1.0.0.dylib. I looked for the answer on the Internet but still couldn't fix that. I stucked at this code. Would you tell me how to fix it? Thank you so much!

[@ShivangiSinghMBB]

Can please post a tutorial for vcf file analysis with R

[@Bioinformagician]

Definitely have plans on covering this topic. Thanks for the suggestion :)

[@juanete69]

What happens with reads aligned multiple times?

[@quinattasneemrafique536]

Can you please make a detailed video on Metagenomics

[@Bioinformagician]

Thanks for the suggestion. I am currently working on some metagenomics content. Please stay tuned!

[@quinattasneemrafique536]

​@@Bioinformagicianthat's great news, waiting eagerly 😍

[@juliangrandvallet5359]

So... after this, how to analyze it on R?

[@anamikapandey4769]

hello ma'am , i just wanna say that when i am running your provided pipelines with code then at the time of creating features counts , it is showing error "failed to open the annotation file , format is incorrect or it contains no exon features" . please correct me , as i have tried many ways but didn't get the answer.

[@seatballmub]

I figured it out. I extracted the 'Homo_sapiens.GRCh38.106.gtf.gz' file to my working directory under the filename 'homo38.106.gtf' for simplicity. Then, I entered the following command exactly into the terminal: featureCounts -S 2 -a homo38.106.gtf -o quants/demo_featurecounts.txt HISAT2/demo_trimmed.bam I don't know why this worked/was even necessary for me. You'll notice this omits the ../ preceding the filename in her example code that I've copied below. featureCounts -S 2 -a ../hg38/Homo_sapiens.GRCh38.106.gtf -o quants/demo_featurecounts.txt HISAT2/demo_trimmed.bam

[@kathirvelm755]

Great mam. .I need help for building cnv pipeline.. will be a great help if u tell anyway for contacting you .

[@Bioinformagician]

You can find my email ID under every video.

[@mostafaismail4253]

can you make a tutorial for NEXTFLOW or snakemake?

[@Bioinformagician]

I can surely consider making a video on it. Thanks for the suggestion!

[@mostafaismail4253]

@@Bioinformagician Thanks

[@andrewzlobin7011]

Did you listen your presentation?

[@mohammedfarhanlakdawala4326]

Do you plan to have a similar pipeline in R ?

[@Bioinformagician]

You mean executing these bash commands using system() in R or using R packages for command line versions of these tools?

[@mohammedfarhanlakdawala4326]

@@Bioinformagician using R packages for command line.

[@Bioinformagician]

@@mohammedfarhanlakdawala4326 I prefer to not perform alignment and quantification in R because R is not very efficient in memory management.

[@mohammedfarhanlakdawala4326]

@@Bioinformagician thank you so much!