Chip-seq上游分析流程学习(二)

# 由于笔者这里存储的有点乱，所以会比较复杂
# 使用者在自行使用的时候，如果觉得设置路径很麻烦就把文件全放一个文件夹下面，先学会再说
# 定义变量，指定工作目录路径
# path路径需要通过pwd确定

# 可以选择写一个shell脚本
# 把下面代码输入进去
#!/bin/bash

# 也可以直接复制这些代码进行运行
path="/home/lm/Z_Projects/chipseq"
# 可以cd路径或者手动cd
# cd ${path}

mkdir clean

cat SRR_Acc_List.txt | while read id;
do
 nohup fastp \
    -i ${path}/fastq/${id}.fastq.gz \
    -o ${path}/clean/${id}.fq.gz \
    -j ${path}/clean/${id}.fastp.json \
    -h ${path}/clean/${id}.fastp.html &
done

# -i 输入的文件名及其文件路径
# -o 输出的fq.gz文件名及其文件路径
# -j 输出的json文件名及其文件路径
# -h 输出的html文件名及其文件路径

# 需要先创建clean文件夹并进入，这样有助于不同文件的归类

path="/home/lm/Z_Projects/chipseq"
# 可以cd路径或者手动cd
# cd ${path}

mkdir clean

ls ./fastq/*.gz |while read fq_res
do
nohup trim_galore -q 25 --phared33 --length 35 -e 0.1 --stringency 5 -o
$path/clean $fq_res &
done

# 下载基因组文件，这里顺便把另外几个常用的基因组文件及基因组注释文件的下载方式一并展示
cd
mkdir reference
mkdir reference/hg19 reference/hg38 reference/mm10 reference/mm39

# 下载hg19（GRCh37）
cd reference/hg19
## hg19的基因组文件
wget https://hgdownload.soe.ucsc.edu/goldenPath/hg19/bigZips/hg19.fa.gz
gzip -d hg19.fa.gz
## hg19的基因组注释文件
wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_42/GRCh37_mapping/gencode.v42lift37.annotation.gtf.gz
gzip -d gencode.v42lift37.annotation.gtf.gz

# 下载hg38（GRCh38）
cd
cd reference/hg38
## hg38的基因组文件
wget http://hgdownload.cse.ucsc.edu/goldenPath/hg38/bigZips/hg38.fa.gz
gzip -d hg38.fa.gz
## hg38的基因组注释文件
wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_42/gencode.v42.annotation.gtf.gz
gzip -d gencode.v42.annotation.gtf.gz

# 下载mm10（GRCm38）
cd
cd reference/mm10
## mm10的基因组文件
wget ftp://hgdownload.soe.ucsc.edu/goldenPath/mm10/bigZips/mm10.fa.gz
gzip -d mm10.fa.gz
## mm10的基因组注释文件
wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_mouse/release_M25/gencode.vM25.annotation.gtf.gz
gzip -d gencode.vM25.annotation.gtf.gz

# 下载mm39（GRCm39）
cd
cd reference/mm39
## mm39的基因组文件
wget ftp://hgdownload.soe.ucsc.edu/goldenPath/mm39/bigZips/mm39.fa.gz
gzip -d mm39.fa.gz
## mm39的基因组注释文件
wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_mouse/release_M32/gencode.vM32.annotation.gtf.gz
gzip -d gencode.vM32.annotation.gtf.gz

# 笔者这边稍作了修改
path="/home/lm/Z_Projects/chipseq"
mkdir -p ${path}/reference/hg38  
nohup wget -O ${path}/reference/hg38/hg38.fa.gz http://hgdownload.cse.ucsc.edu/goldenPath/hg38/bigZips/hg38.fa.gz &

# 解压缩
gzip -d hg38.fa.gz

# 笔者这边稍作了修改
path="/home/lm/Z_Projects/chipseq"
mkdir -p ${path}/reference/hg38  

nohup wget -O ${path}/reference/hg38/gencode.v47.annotation.gtf.gz https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_47/gencode.v47.annotation.gtf.gz
gzip -d gencode.v47.annotation.gtf.gz

# gff3文件
# nohup wget -O ${path}/reference/hg38/gencode.v47.annotation.gff3.gz https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_47/gencode.v47.annotation.gff3.gz

# 解压缩
gzip -d gencode.v47.annotation.gtf.gz

# 构建hg38的bowtie2索引文件
mkdir -p index/bowtie2/hg38_res

# --threads设置线程数
bowtie2-build ${path}/reference/hg38/hg38.fa  ${path}/index/bowtie2/hg38_res/hg38 --threads 8

# check一下
ls ./index/bowtie2/hg38_res/ -lh

export path="/home/lm/Z_Projects/chipseq"
export bowtie2_index="/home/lm/Z_Projects/chipseq/index/bowtie2/hg38_res/hg38"
mkdir -p ${path}/align

nohup sh -c 'cat SRR_Acc_List.txt | while read id; do
  bowtie2 -p 8 -x ${bowtie2_index} -U ${path}/clean/${id}.fq.gz | samtools sort -O bam -@ 8 -o - > ${path}/align/${id}.bam
  samtools flagstat ${path}/align/${id}.bam > ${path}/align/${id}_flagstat.txt
done' > ${path}/align/bowtie2_run.log 2>&1 &