Making a Pipeline

Combining cells can make things easier

Shell Variables

Assign the variables in this notebook.

source bioinf_intro_config.sh
mkdir -p $TRIMMED $STAR_OUT

Combining the Steps

FASTQ="21_2019_P_M1_S21_L002_R1"
echo "---------------- TRIMMING: $FASTQ ----------------"
fastq-mcf $MYINFO/neb_e7600_adapters.fasta \
    $RAW_FASTQS/${FASTQ}_001.fastq.gz \
        -q 20 -x 0.5 \
        -o $TRIMMED/${FASTQ}_001.trim.fastq.gz

echo "---------------- MAPPING: $FASTQ ----------------"
STAR \
    --runMode alignReads \
    --twopassMode None \
    --genomeDir $GENOME_DIR \
    --readFilesIn $TRIMMED/${FASTQ}_001.trim.fastq.gz \
    --readFilesCommand gunzip -c \
    --outFileNamePrefix ${STAR_OUT}/${FASTQ}_ \
    --quantMode GeneCounts \
    --outSAMtype None \
    --runThreadN 2

---------------- TRIMMING: 21_2019_P_M1_S21_L002_R1 ----------------
Command Line: /home/jovyan/work/scratch/bioinf_intro/myinfo/neb_e7600_adapters.fasta /data/hts_2019_data/hts2019_pilot_rawdata/21_2019_P_M1_S21_L002_R1_001.fastq.gz -q 20 -x 0.5 -o /home/jovyan/work/scratch/bioinf_intro/trimmed_fastqs/21_2019_P_M1_S21_L002_R1_001.trim.fastq.gz
Scale used: 2.2
Phred: 33
Threshold used: 751 out of 300000
Adapter Adapter (AGATCGGAAGAGCACACGTCTGAACTCCAGTCA): counted 2515 at the 'end' of '/data/hts_2019_data/hts2019_pilot_rawdata/21_2019_P_M1_S21_L002_R1_001.fastq.gz', clip set to 6
Files: 1
Total reads: 2437108
Too short after clip: 1347
Clipped 'end' reads: Count: 44977, Mean: 15.55, Sd: 8.27
Trimmed 288960 reads by an average of 1.70 bases on quality < 20
---------------- MAPPING: 21_2019_P_M1_S21_L002_R1 ----------------
Jun 26 16:20:46 ..... started STAR run
Jun 26 16:20:46 ..... loading genome
Jun 26 16:20:46 ..... started mapping
Jun 26 16:22:09 ..... finished successfully

And let’s check the result

ls ${STAR_OUT}

21_2019_P_M1_S21_L001_R1_short_introns_Aligned.sortedByCoord.out.bam
21_2019_P_M1_S21_L001_R1_short_introns_Aligned.sortedByCoord.out.bam.bai
21_2019_P_M1_S21_L001_R1_short_introns_Log.final.out
21_2019_P_M1_S21_L001_R1_short_introns_Log.out
21_2019_P_M1_S21_L001_R1_short_introns_Log.progress.out
21_2019_P_M1_S21_L001_R1_short_introns_ReadsPerGene.out.tab
21_2019_P_M1_S21_L001_R1_short_introns_SJ.out.tab
21_2019_P_M1_S21_L002_R1_Aligned.out.bam
21_2019_P_M1_S21_L002_R1_Log.final.out
21_2019_P_M1_S21_L002_R1_Log.out
21_2019_P_M1_S21_L002_R1_Log.progress.out
21_2019_P_M1_S21_L002_R1_ReadsPerGene.out.tab
21_2019_P_M1_S21_L002_R1_short_introns_Aligned.sortedByCoord.out.bam
21_2019_P_M1_S21_L002_R1_short_introns_Aligned.sortedByCoord.out.bam.bai
21_2019_P_M1_S21_L002_R1_short_introns_Log.final.out
21_2019_P_M1_S21_L002_R1_short_introns_Log.out
21_2019_P_M1_S21_L002_R1_short_introns_Log.progress.out
21_2019_P_M1_S21_L002_R1_short_introns_ReadsPerGene.out.tab
21_2019_P_M1_S21_L002_R1_short_introns_SJ.out.tab
21_2019_P_M1_S21_L002_R1_SJ.out.tab
genome_Log.out
multiqc_data
multiqc_report.html

head ${STAR_OUT}/${FASTQ}_ReadsPerGene.out.tab

N_unmapped      46060   46060   46060
N_multimapping  35006   35006   35006
N_noFeature     12466   2145291 18783
N_ambiguous     203327  820     316
CNAG_04548      0       0       0
CNAG_07303      0       0       0
CNAG_07304      6       0       6
CNAG_00001      0       0       0
CNAG_07305      1       0       1
CNAG_00002      51      0       51