Pipeline from FASTQ to DEG

This document gives the instructions for aligning the FASTQ sequence files from 3'-DGE mRNA sequencing experiments, counting the number of uniquely tagged and aligned sequence fragments, and identifying differentially expressed genes (DEGs) from their UMI counts.

Directory Structure

The structure of this computational pipeline including both datasets and programs:

LINCSData
├── Datasets                                         # Input and output datasets
│   ├── Alignment                                    # Sequence, alignment and counts
│   │   └── LINCS.Dataset
│   │       └── LINCS.Dataset.Gene.20150409          # Plate 1
│   │           ├── Aligns                           # SAM alignment files
│   │           │   └── hg38
│   │           ├── Counts                           # UMI counts files
│   │           │   └── hg38
│   │           └── Seqs                             # Demultiplexed FASTQ files
│   │               └── Demux
│   └── Difference                                   # Data for DEG calculation
│       └── LINCS.Dataset
│           └── LINCS.Dataset.Gene.LINCS.20150409    # Plate 1
│               ├── Configs                          # Control file for DEG calculation
│               ├── Counts                           # UMI counts files
│               ├── Params                           # Parameter files
│               ├── Results                          # Calculation results
│               └── Scripts                          # Dataset-specific R programs
├── Programs
│   ├── DEG                                          # R programs for DEG calculation
│   │   └── Library                                  # R utlity functions
│   ├── STAR                                         # Shell scripts for STAR
│   ├── Subread                                      # Shell scripts for featureCounts
│   └── UMICount                                     # C++ program for unique UMI counts
└── References                                       # Reference genome
    └── UCSC                                         # UCSC human genome library
        └── Human
            └── hg38                                 # Version hg38
                ├── Annotation                       # GTF annotation file
                ├── Sequence                         # Genome file
                └── STAR-Index-2.5.3a-RL46           # STAR index files for alignment

Program Files

This pipeline makes use of the following software and programs.

Software Prerequisites

The following software programs need to be installed before running the pipeline:

• STAR
• Subread
• R
  • edgeR package
  • matrixStats package
  • gplots package

Custom programs

The in-house programs developed for the pipeline include:

Programs/
├── DEG
│   ├── Compare-Molecule-Expression.R                # Identify DEGs
│   └── Extract-Gene-Expression-Samples.R            # Post-process UMI counts
├── STAR
│   ├── STAR-Align-Seq.sh                            # Call STAR to align FASTQ files
│   └── STAR-Align-Seqs.Single.sh                    # Align plate-1 FASTQ files
├── Subread
│   ├── FeatureCounts-Count-Aligns.sh                # Call featureCounts for counting
│   ├── Launch-FeatureCounts-Count-Aligns.STAR.sh    # Count STAR's SAM
│   └── Launch-FeatureCounts-Count-Aligns.UMI.sh     # Count SAM-Alignment-Counter' SAM
└── UMICount
    ├── Move-Report-Aligns.sh                        # Move featureCounts' reported SAM
    ├── SAM-Alignment-Counter                        # Select alignments with unique UMI
    └── SAM-Count-Aligns.sh                          # Script for SAM-Alignment-Counter

Pipeline Workflow

This computational pipeline consists of three steps: setup the pipeline directory, align sequence fragments, and compare gene expressions.

Setup Pipeline Directory

1. Create a top directory for the pipeline using the command:

   mkdir "${HOME}/DToxS-LINCS-Dataset"

2. Download the directory LINCSData with all its contents to the top pipeline directory created above.

3. Set the environmental variable TOP_DIR to the top pipeline directory :

   export TOP_DIR="${HOME}/DToxS-LINCS-Dataset"

Alignment of Sequence Fragments

This part is to align the FASTQ sequence files produced by mRNA sequencing experiment to the human reference genome library and generate a table of gene-wise UMI counts for the aligned sequences tagged with unique UMI barcodes.

Align FASTQ Sequences

Run the following command to align the demultiplexed FASTQ files:

${TOP_DIR}/LINCSData/Programs/STAR/STAR-Align-Seqs.Single.sh [Number of Threads]

where [Number of Threads]is the number of computing threads to use for alignment.

Count STAR Alignments

Run the following command to count the alignments from STAR output:

${TOP_DIR}/LINCSData/Programs/Subread/Launch-FeatureCounts-Count-Aligns.STAR.sh [Number of Threads]

where [Number of Threads] is the number of computing threads to use for alignment.

Move Reported Alignments

Run the following command to move the SAM alignment files with addtional gene assignment information reported by featureCounts from the counts directory to the alignment directory:

${TOP_DIR}/LINCSData/Programs/UMICount/Move-Report-Aligns.sh "${TOP_DIR}/LINCSData/Datasets/Alignment/LINCS.Dataset/LINCS.Dataset.Gene.20150409/Counts/hg38" "${TOP_DIR}/LINCSData/Datasets/Alignment/LINCS.Dataset/LINCS.Dataset.Gene.20150409/Aligns/hg38"

Select Unique UMIs

Run the following command to select the alignments tagged with unique UMIs from the SAM files reported by featureCounts:

${TOP_DIR}/LINCSData/Programs/UMICount/SAM-Count-Aligns.sh "${TOP_DIR}/LINCSData/Datasets/Alignment/LINCS.Dataset/LINCS.Dataset.Gene.20150409/Aligns/hg38"

Count Unique Alignments

Run the following command to count the SAM alignments tagged with unique UMIs from SAM-Alignment-Counter output:

${TOP_DIR}/LINCSData/Programs/Subread/Launch-FeatureCounts-Count-Aligns.UMI.sh [Number of Threads]

where [Number of Threads]is the number of computing threads to use for alignment.

Comparison of Gene Expressions

This part is to study the change of gene expresion level caused by the drug treatment on cell samples using the gene-wise UMI counts table obtained in the previous part, and identify the differentially expressed genes at each treatment condition.

Copy Counts File

Run the following command to copy the gene-wise unique UMI counts file from the Alignments' Counts directory to the Difference's Counts directory:

copy -a "${TOP_DIR}/LINCSData/Datasets/Alignment/LINCS.Dataset/LINCS.Dataset.Gene.20150409/Counts/hg38/RNAseq-Read-Counts.UMI.hg38.20150409.txt" "${TOP_DIR}/LINCSData/Datasets/Difference/LINCS.Dataset/LINCS.Dataset.Gene.LINCS.20150409/Counts/"

Extract LINCS Samples

Run the following command to extract the UMI counts of LINCS samples from the UMI counts table and synchronize it with the experiment design file:

Rscript "${TOP_DIR}/LINCSData/Programs/DEG/Extract-Gene-Expression-Samples.R" "${TOP_DIR}/LINCSData/Datasets/Difference/LINCS.Dataset/LINCS.Dataset.Gene.LINCS.20150409/Counts" "${TOP_DIR}/LINCSData/Programs/DEG/Library"

Identify Differentially Expressed Genes

Run the following command to identify the DEGs at each drug treatment condition of plate 1:

Rscript "${TOP_DIR}/LINCSData/Programs/DEG/Compare-Molecule-Expression.R" "${TOP_DIR}/LINCSData/Datasets/Difference/LINCS.Dataset/LINCS.Dataset.Gene.LINCS.20150409/Configs/Configs.LINCS.Dataset.Gene.LINCS.20150409.tsv" "${TOP_DIR}/LINCSData/Datasets/Difference/LINCS.Dataset/LINCS.Dataset.Gene.LINCS.20150409" "${TOP_DIR}/LINCSData/Programs/DEG/Library"

All of the result data files are generated in the Results directory.