上游分析的最后一部分祸穷，這里是佳奧性穿，讓我們繼續(xù)吧！
1 獲取實驗數(shù)據(jù)（雙端測序）
raw_fq：下載sra文件雷滚，轉fq文件需曾。
clean_fq：trim_galory
ls *_1.fastqc.gz >1
ls *_2.fastqc.gz >2
paste 1 2 > config
dir='../clean_fq'
cat $config_file | while read id
do
arr=($id)
fq1=${arr[0]}
fq2=${arr[1]}
echo $dir $fq1 $fq2
trim_galore -q 25 --phred33 --length 35 -e 0.1 --stringency 3 --paired -o $dir $fq1 $fq2 &
done
qc：質量控制，并查看.html報告祈远。
2 hisat2比對
2.1 建立索引
##網上下載或者自行構建
##構建索引
hisat2-build -p 16 Arabidopsis_thaliana.TAIR10.28.dna.genome.fa genome
align：比對后的.bam/.sam文件呆万。
2.2 批量比對
index=../hisat2index
fq1=..
fq2=..
sample=${arr[0]}
if((i%$number1==$number2))
then
if[! -f $sample.bam]; then
start=$(date +%s.%N)
echo hisat2 'date'
hisat2 -P 4 -x $index -1 $fq1 -2 $fq2 | samtools sort -@ 4 -o $sample.bam
echo hisat2 'date'
dur=$(echo "$(date +%s.%N) - $start" | bc)
printf "Execution time for hisat2 : %.6f seconds" $dur
fi##endforffiles
fi##end for number1
i=$((i+1))
done ##end for $config- file
3 bamCoverage轉bw文件
##方法一
analysis_dir=$1
config_file=$2
number1=$3
number2=$4
cat $config_file | while read id
do
echo $id
echo $id
file=$(basename $id )
sample=${file%%.*}
echo $sample
if((i%$number1==$number2))
then
if [ ! -f $sample.bw ]; then
start=$(date +%s.%N)
echo bamCoverage 'date'
bamCoverage -b $id -o $sample.bw --normalizeUsing RPKM -p 4
echo bamCoverage 'date'
dur=$(echo "$(date +%s.%N) - $start" | bc)
printf "Execution time for bamCoverage : %. 6f seconds" $dur
fi##end for if files
fi##end for number1
i=$((i+1))
done##end for $config_file
##方法二、一次性全部提交车份，如果是公共服務器等著管理員來聯(lián)系你（所以會有方法一的方法）
ls *bam| while read id; do(bamCoverage -b $id -o ${id/.merge.bam/.bw} --normalizeUsing RPKM -p 4 &);done
4 featureCounts差異分析
##方法一谋减、批量bam featureCounts
gtf='/home/kaoku/rnaseq/biotree_plant/refer/Arabidopsis_thaliana.TAIR10.28.gtf.gz'
featureCounts -T 5 -p \
-a $gtf -o all.counts.txt \
/home/kaoku/rnaseq/biotree_plant/data/sam_bam_bai/*.bam
##方法二
featureCounts -T 4 -p -t exon -g gene_name -a $gtf -o all.counts.id.txt ../bams/*.bam 1>counts.id.log 2>&1
最后會獲得一個矩陣文件，all.id.counts.txt扫沼，可以再multiqc ./看一下結果出爹。
至此上游分析結束庄吼。
后續(xù)的內容便是把表觀調控的圖表重現(xiàn)出來。
我們下一篇再見严就！
