TECHNICAL WORKFLOW GUIDE

RNA-Seq Analysis Step by Step: From FASTQ to Biological Discovery

Author: Anjali Chaudhary, Director & Bioinformatics Scientist Published: October 2026 Last updated: October 2026 Reading time: 10 min read

RNA-Seq data analysis converts millions of short sequencing reads into biologically meaningful gene expression profiles. The standard reference-based pipeline proceeds systematically through five stages: raw read quality assessment with FastQC, adapter and low-quality trimming via fastp, splice-aware genomic alignment with STAR, gene-level feature quantification with featureCounts or Salmon, and negative binomial differential expression modeling with DESeq2.

Key Takeaways

  • Quality control before and after trimming prevents optical artifacts and adapter read-through from skewing alignments.
  • Splice-aware aligners like STAR or HISAT2 are mandatory for eukaryotic transcriptomes to span exon-exon junctions across introns.
  • Differential expression packages require raw un-normalized integer read counts as input; never feed FPKM, RPKM, or TPM into DESeq2 or edgeR.
  • Biological replicates (minimum n=3 per condition) are essential to estimate biological variation across negative binomial dispersions.

What is the Standard RNA-Seq Workflow Architecture?

Bulk RNA sequencing (RNA-Seq) is the premier tool for interrogating transcriptome dynamics across physiological states, disease models, and environmental perturbations. Modern high-throughput sequencers (e.g., Illumina NovaSeq X, NextSeq 2000) deliver tens of millions of paired-end reads per sample formatted as compressed FASTQ files.

However, transforming gigabytes of raw base calls and Phred quality scores into statistically validated lists of upregulated and downregulated genes requires a standardized computational architecture. Each phase of the pipeline builds upon the previous, demanding stringent quality checkpoints to ensure reproducible science.

How Do You Assess Raw Sequencing Quality with FastQC and MultiQC?

Before executing sequence alignment, every computational biologist must inspect raw read characteristics. FastQC assesses base-level metrics across the sequencing run, checking:

When managing multiple samples, running MultiQC aggregates individual FastQC reports into a single interactive HTML summary, highlighting outliers across biological replicates.

Why is Read Preprocessing and Trimming Critical?

Suboptimal reads compromise alignment efficiency and introduce spurious mismatches. Modern preprocessing tools like fastp execute high-speed adapter clipping, poly-A/poly-T tail removal, and sliding-window quality trimming.

Unlike older tools like Trimmomatic, fastp automatically detects adapter sequences without manual FASTA specification and generates simultaneous pre- and post-filtering QC reports. Reads that fall below length thresholds (typically <35 bp post-trimming) are discarded to prevent non-specific multi-mapping across the reference genome.

How Does Splice-Aware Alignment Function with STAR and HISAT2?

Unlike genomic DNA sequencing where reads map contiguously to chromosome assemblies, mature eukaryotic mRNA transcripts lack introns. Consequently, reads spanning exon-exon junctions cannot map linearly to genomic DNA.

Standard genomic mappers (BWA-MEM or Bowtie2) fail at splice junctions. Instead, RNA-Seq requires splice-aware aligners:

Both aligners produce coordinate-sorted BAM (Binary Alignment Map) files. Post-alignment indexing with samtools allows rapid visual inspection of coverage tracks in the Integrative Genomics Viewer (IGV).

How are Reads Quantified into Count Matrices?

Once reads are aligned, the next phase calculates the number of reads assigned to each genomic feature (genes, transcripts, or exons) based on reference annotations (GTF/GFF3).

featureCounts (part of Subread) is the gold standard for counting mapped reads against genomic intervals. Alternatively, pseudoaligners such as Salmon or Kallisto bypass full BAM generation, quantifying transcript abundances directly from raw FASTQs in minutes while correcting for GC and sequence-specific fragment biases.

Crucially, differential expression testing tools (DESeq2, edgeR) demand raw integer counts. Normalized values like TPM (Transcripts Per Kilobase Million) or FPKM (Fragments Per Kilobase of transcript per Million mapped reads) remove library-size variation mathematically required for negative binomial statistical models and should never be fed into count-based testing engines.

How Do You Model Differential Gene Expression with DESeq2?

RNA-Seq count data exhibit severe overdispersion: the variance in gene expression across biological replicates substantially exceeds the mean. Ordinary linear regression and Poisson models underestimate false discovery rates under these conditions.

DESeq2 models count data using the Negative Binomial distribution with dispersion parameters estimated via empirical Bayes shrinkage. The workflow in R proceeds as follows:

  1. Median-of-Ratios Normalization: Scales counts for differences in sequencing depth and compositional skew caused by highly expressed outlier genes.
  2. Dispersion Estimation: Shrinks gene-wise dispersion estimates towards a fitted mean-dispersion curve, providing robust statistical power even with modest sample sizes (n=3).
  3. Generalized Linear Model (GLM) Fitting: Fits logarithmic fold changes and tests for significance using the Wald test or Likelihood Ratio Test (LRT).
  4. False Discovery Rate (FDR) Correction: Adjusts p-values using the Benjamini-Hochberg procedure (yielding adjusted p-values / padj) to control type-I errors across 20,000+ simultaneous gene tests.

How Do You Interpret and Visualize Transcriptomic Discoveries?

The output of DESeq2 is an annotated table of Log2 Fold Changes (LFC) and adjusted p-values. Standard downstream exploratory figures include:

Need End-to-End RNA-Seq Analysis for Your Laboratory?

Xenoi Biotech provides complete, reproducible transcriptomics data analysis services. We take your raw FASTQ files and deliver MultiQC audits, indexed alignments, DESeq2 statistics, and publication-ready figures.

Explore Analysis Services Learn in our 30-Day Workshop

Further Reading & Authoritative Literature

AC

About the Author: Anjali Chaudhary

Director & Bioinformatics Scientist · Xenoi Biotech

Anjali Chaudhary leads genomic data analysis and curriculum design at Xenoi Biotech in Punjab, India. Her expertise encompasses next-generation transcriptomic pipelines, differential expression modeling, and scalable bioinformatics architectures.