100天生信-Day10
记录一下workflow
ChIP-seq workflow:
## QC: fastp
# filter fastq (IP)
fastp -i Ta_Input_1.fq.gz -I Ta_Input_2.fq.gz -o Ta_Input_1_filter.fq.gz \
-O Ta_Input_2_filter.fq.gz -w 8 > fastp_input.log 2>&1 &
# filter fastq (Input)
fastp -i Ta_IP_1.fq.gz -I Ta_IP_2.fq.gz -o Ta_IP_1_filter.fq.gz \
-O Ta_IP_2_filter.fq.gz -w 8 > fastp_IP.log 2>&1 &
## mapping: bowtie2:
# indevx
bowtie2-build Triticum_aestivum.IWGSC.dna.toplevel.fa Ta_index > bowtie2_Ta_index.log 2>&1 &
# mapping IP
bowtie2 -p 8 --very-sensitive -x /Users/lichuanshun/Documents/Bio_data/Ta_genome/Ensemble/Ta_bowtie_index/Ta_index \
-1 Ta_IP_1_filter.fq.gz -2 Ta_IP_2_filter.fq.gz -S Ta_IP.sam > bowtie_IP.log 2>&1 &
# mapping Input
bowtie2 -p 8 --very-sensitive -x /Users/lichuanshun/Documents/Bio_data/Ta_genome/Ensemble/Ta_bowtie_index/Ta_index \
-1 Ta_Input_1_filter.fq.gz -2 Ta_Input_2_filter.fq.gz -S Ta_Input.sam > bowtie_Input.log 2>&1 &
## filter bam: samtools & sambamba
# sam -> bam
samtools view -bS Ta_IP.sam > Ta_IP.bam
samtools view -bS Ta_Input.sam > Ta_Input.bam
# sort
samtools sort -@ 8 -o Ta_IP.sort.bam Ta_IP.bam
samtools sort -@ 8 -o Ta_Input.sort.bam Ta_Input.bam
# filter (rmdup)
sambamba view -h -t 8 -f bam -F 'not unmapped and not duplicate and mapping_quality >= 10' Ta_IP.sort.bam > Ta_IP.rmdup.bam
sambamba view -h -t 8 -f bam -F 'not unmapped and not duplicate and mapping_quality >= 10' Ta_Input.sort.bam > Ta_Input.rmdup.bam
# filter (uni-mapped)
sambamba view -h -t 8 -f bam -F '[XS] == null and not unmapped and not duplicate and mapping_quality >= 10' Ta_IP.sort.bam > Ta_IP.uni.bam
sambamba view -h -t 8 -f bam -F '[XS] == null and not unmapped and not duplicate and mapping_quality >= 10' Ta_Input.sort.bam > Ta_Input.uni.bam
## call peaks
# uni-mapped reads
macs2 callpeak -t Ta_IP.uni.bam -c Ta_Input.uni.bam -f BAM -n Ta_IP -g 17e9 \
--nomodel --shift -100 --extsize 200 -B -p 0.0001 --outdir peaks_result_00001_uni > MACS2_00001_uni.log 2>&1 &
# all reads
macs2 callpeak -t Ta_IP.rmdup.bam -c Ta_Input.rmdup.bam -f BAM -n Ta_IP -g 17e9 \
--nomodel --shift -100 --extsize 200 -B -p 0.0001 --outdir peaks_result_00001_rmdup > MACS2_00001_rmdup.log 2>&1 &