Analyzing RNA-Seq data for differential gene expression

Materials and Software

  • Data files
    • RNA-Seq data in FASTQ format (Ex: dataset1.fastq, dataset2.fastq)
  • Genome sequence files
    • Genome sequence in FASTA format (Ex: REL606.fna)
    • Genome sequence gene annotations in GFF3 format (Ex: REL606.gff3)
  • Adaptor filtering software *FAR Manual | Download
  • Read mapper software
    • Bowtie2 (first choice)
      • Download Bowtie
      • To build, just type "make" in the source code directory.
      • Add this directory to your $PATH or move bowtie, and bowtie-* star executable to your $PATH
    • BWA (alternate choice)
      • Download BWA
      • To build, just type "make" in the source code directory.
      • To install, move the executable "bwa" to somewhere in your $PATH, like $HOME/local/bin.
      • For usage see the BWA manual.
  • breseq
  • R statistics package
  • Bioconductor R modules
    • library(edgeR)
    • library(DESEQ)


Remove adaptor sequences from reads

You will need a FASTA file of adaptor sequences.

For each input file you will need to run this command (single-end data):
$far --source datasetX.fastq --target datasetX.noadaptor --adaptive-overlap --trim-end any --adapters gsaf_illumina_adapters.fasta --format fastq-sanger

There is an option to process paired-end data like this:
$far --source datasetX_R1.fastq --source2 datasetX_R2.fastq --target datasetX.noadaptor --adaptive-overlap --trim-end any --adapters gsaf_illumina_adapters.fasta --format fastq-sanger

Compile and install FAR on MacOSX

Unfortunately, FAR comes only with Windows and Linux binaries. To build FAR(2.0) for MacOSX:

  1. Install Apple Developer Tools
  2. Install cmake: $sudo port install cmake
  3. Check out code:
    $ svn co theflexibleada
  4. Compile code:
    $ cd theflexibleada
    $ cmake CMakeLists.txt
  5. Copy executable and library:
    $ cp lib/libtbb.dylib ~/local/lib
    $ cp build/* ~/local/bin
  6. Add these locations to your path with lines in ~/.profile:
    export PATH="$PATH:$HOME/local/bin"

Align reads to reference genome

Using bowtie2

First, index your genome so bowtie2 can map read to it:
$bowtie2-build REL606.fna REL606

Then, align each data set:
$bowtie2 -x REL606 -U datasetX.fastq --phred33 -S REL606.sam

Optionally, add the --local flag if your reads do not map end-to-end.

Using BWA

First, index your genome so BWA can map read to it:
$bwa index REL606.fna

Then, align each data set:
$bwa aln REL606.fna dataset1.fastq > datasetX.sai

And convert to SAM format (assumes single-end data): $bwa samse REL606.fna datasetX.sai datasetX.fastq > datasetX.sam

Count reads mapping to genes

breseq RNASEQ -f REL606.fna -r REL606.gbk -o datasetX.sam

Convert alignments to BAM

And convert to BAM format (assumes single-end data):
$samtools faidx REL606.fna
$samtools import REL606.fna datasetX.sam datasetX.unsorted.bam
$samtools sort datasetX.unsorted.bam datasetX
$samtools index datasetX.bam

Now you can use IGV to view them.

Analyze differential gene expression

Using DESeq

Topic attachments
I Attachment Action Size Date Who Comment
elsefasta gsaf_illumina_adapters.fasta manage 0.2 K 30 Jan 2012 - 16:31 Main.JeffreyBarrick  

 Barrick Lab  >  ProtocolList  >  ProtocolsRNASeqDifferentialExpression

Contributors to this topic edittopic JeffreyBarrick
Topic revision: r6 - 11 Feb 2012 - 17:41:54 - Main.JeffreyBarrick
This site is powered by the TWiki collaboration platformCopyright ©2020 Barrick Lab contributing authors. Ideas, requests, problems? Send feedback