Almost all steps of a short-read analysis pipeline can use fastq.gz (gzipped) files, rather than fastq, to save storage space. For example,
# Untar from solexa_public, gunzip, and gzip
tar -xOzf /lab/solexa_public/LAB/RUN/QualityScore/Foo.1.tar.gz | gzip -f > fastq/My_sample.1.fq.gz
tar -xOzf /lab/solexa_public/LAB/RUN/QualityScore/Foo.2.tar.gz | gzip -f > fastq/My_sample.2.fq.gz
# [Whitehead fastq files only] Untar from solexa_public, gunzip, remove WI GTC-specific header suffixes and gzip
tar -xOzf /lab/solexa_public/LAB/RUN/QualityScore/Foo.1.tar.gz | /nfs/BaRC_Public/BaRC_code/Perl/rm_WI_Illumina_suffix/rm_WI_Illumina_suffix.STDIN.pl | gzip -f > fastq/My_sample.1.fq.gz
tar -xOzf /lab/solexa_public/LAB/RUN/QualityScore/Foo.2.tar.gz | /nfs/BaRC_Public/BaRC_code/Perl/rm_WI_Illumina_suffix/rm_WI_Illumina_suffix.STDIN.pl | gzip -f > fastq/My_sample.2.fq.gz
# Run FastQC
fastqc fastq/My_sample.1.fq.gz
fastqc fastq/My_sample.2.fq.gz
# Trim adapters
gunzip -c fastq/My_sample.1.fq.gz | fastx_clipper -v -z -l 30 -a GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGATGTATCTCGTATGCCGTCTTCTGCTTG -o fastq/My_sample.noVec.1.fq.gz
gunzip -c fastq/My_sample.2.fq.gz | fastx_clipper -v -z -l 30 -a GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGATGTATCTCGTATGCCGTCTTCTGCTTG -o fastq/My_sample.noVec.2.fq.gz
fastqc fastq/My_sample.noVec.1.fq.gz
fastqc fastq/My_sample.noVec.2.fq.gz
# Drop reads of low quality
gunzip -c fastq/My_sample.noVec.1.fq.gz | fastq_quality_filter -v -q 20 -p 75 -z -o fastq/My_sample.noVec.qgt20.1.fq.gz
gunzip -c fastq/My_sample.noVec.2.fq.gz | fastq_quality_filter -v -q 20 -p 75 -z -o fastq/My_sample.noVec.qgt20.2.fq.gz
# Trim low-quality bases at end of reads
gunzip -c fastq/My_sample.noVec.qgt20.1.fq.gz | fastq_quality_trimmer -v -t 20 -l 30 -o fastq/My_sample.noVec.qgt20.t.1.fq
gunzip -c fastq/My_sample.noVec.qgt20.2.fq.gz | fastq_quality_trimmer -v -t 20 -l 30 -o fastq/My_sample.noVec.qgt20.t.2.fq
# After trimming, get reads that are still paired (and output as regular fastq for next (bowtie) step)
# We hope to modify cmpfastq.pl to accept gz files as input
/nfs/BaRC_Public/BaRC_code/Perl/cmpfastq/cmpfastq.pl fastq/My_sample.noVec.qgt20.t.1.fq fastq/My_sample.noVec.qgt20.t.2.fq
# Map paired reads with bowtie
# bowtie cannot use gzipped input files
bowtie -l 30 -n 1 -X 1000 --best -S /nfs/genomes/GENOME/bowtie/INDEX -1 fastq/My_sample.noVec.qgt20.t.1.fq-common.out -2 fastq/My_sample.noVec.qgt20.t.2.fq-common.out mapped_reads/My_sample.30.1.sam
# Sort and index mapped reads
/nfs/BaRC_Public/BaRC_code/Perl/SAM_to_BAM_sort_index/SAM_to_BAM_sort_index.pl mapped_reads/My_sample.30.1.sam
# Delete SAM file (since BAM file has all the same info)
rm -f mapped_reads/My_sample.30.1.sam
# Get summary counts
zcat fastq/My_sample.1.fq.gz | echo $((wc -l/4))
zcat fastq/My_sample.noVec.1.fq.gz | echo $((wc -l/4))
zcat fastq/My_sample.noVec.2.fq.gz | echo $((wc -l/4))
zcat fastq/My_sample.noVec.qgt20.1.fq.gz | echo $((wc -l/4))
zcat fastq/My_sample.noVec.qgt20.2.fq.gz | echo $((wc -l/4))
more fastq/My_sample.noVec.qgt20.t.1.fq | echo $((wc -l/4))
more fastq/My_sample.noVec.qgt20.t.2.fq | echo $((wc -l/4))
zmore fastq/My_sample.noVec.qgt20.t.1.fq-common.out.gz | echo $((wc -l/4))
samtools flagstat mapped_reads/My_sample.30.1.sorted.bam > mapped_reads/My_sample.30.1.flagstat.txt
grep ' mapped (' mapped_reads/My_sample.30.1.flagstat.txt
