双端测序，我先下载了一个样本的试试
注意，trimglore的时候要一次同时处理两个R1R2文件，不能一个一个来
我还不会写循环，写了服务器也挂掉，气人，制作config.raw 参考一下https://www.jianshu.com/p/54f3f348d488
第一步：Trim_galore质控及去接头
~~~
$trim_galore -q 25 --phred33 --length 35 -e 0.1 --stringency 4 --paired -o $输出路径 R1.fastq.gz R2.fastq.gz
~~~
参考：http://www.manongjc.com/detail/57-quhcpsrvumrjbaz.html
第二步：Bowtie2比对并用samtools转bam
组里服务器下有mm10,看一下路径方便用
index在 /storage2/anlei/reference/index/bowtie2/mm10/mm10
~~~
bowtie2 -p 5 --very-sensitive -X 2000 -x /storage2/anlei/reference/index/bowtie2/mm10/mm10 -1 test1.fq -2 test2.fq -S test.sam|samtools sort -@ 5 -O bam -o test.bam -
~~~
ps.我文件夹安装了大家可以直接用 sambamba
第三步：去除PCR重复、线粒体基因
#排序后去重
~~~
samtools sort -@ 50 -n Aging2.bam -o Aging2.sort.bam
samtools fixmate -@ 50 -m Aging2.sort.bam Aging2.sort.fixed.bam
samtools sort -@ 50 Aging2.sort.fixed.bam -o Aging2.position.fixed.bam
samtools markdup -@ 50 -r Aging2.position.fixed.bam Aging2.rmdup.bam
samtools index Aging2.rmdup.bam
samtools flagstat Aging2.rmdup.bam > Aging2.rmdup.stat
#去除低质量q30和chrM并转bed
samtools view -h -f 2 -q 30 Aging2.rmdup.bam |grep -v chrM |samtools sort -O bam -@ 5 -o - > Aging2.last.bam
samtools index Aging2.last.bam
~~~
#计算线粒体DNA比例
ps 诺唯赞的技术回的：一般在测序时会注意线粒体DNA的污染率情况。比如少量测序后，线粒体 DNA 的污染率如果没有超过 50%，则不需要进行线粒体 DNA的去除，那么推荐数据量测到 50 M reads/样（15G），以保证有足够的数据量用于后续分析；若线粒体 DNA 的污染率超过 50%，则推荐进行线粒体 DNA 的去除，去除完之后数据量达到 20-30 M reads/样本（6-9G）即可。https://www.jianshu.com/p/30ebe8ad88fa
~~~
# Calculate %mtDNA:
## mtReads=samtools idxstats Aging2.rmdup.bam | grep 'chrM' | cut -f 3
## totalReads=samtools idxstat Aging2.rmdup.bam | awk '{SUM += $3} END {print SUM}'
## echo '==> mtDNA Content:' $(bc <<< "scale=2;100*$mtReads/$totalReads")'%'
#转bed
bedtools bamtobed -i Aging2.last.bam > Aging2.bed
##查看两次过滤去除情况
samtools idxstats Aging1.bam | grep -v "_" > 1
samtools idxstats Aging1.rmdup.bam | grep -v "_" > 2
samtools idxstats Aging1.last.bam | grep -v "_" > 3
paste 1 2 3 | cut -f 1,3,7,11
~~~
第五步：macs2做call peak
~~~
macs2 callpeak -t Aging2.bed -g mm --nomodel --shift -100 --extsize 200 -n "Aging2" --outdir ../peaks/
~~~
六、bam转bw
~~~
bamCoverage -p 5 --normalizeUsing CPM -b Aging2.last.bam -o Aging2.bw
~~~
七：计算插入片段长度，FRiP，IDR值与deeptools可视化
准备基因组注释文件中提取的bed文件，共三列，染色体，起始坐标和终止坐标
如何下载参考https://mp.weixin.qq.com/s/POPN8kzMQT1jcil8ICvPxg
可视化用deeptools处理数据
computeMatrix reference-point --referencePoint TSS -p 50 -b 2000 -a 2000 -R /storage2/anlei/wangwenjing/ChIP-seq/deeptools/bw/mm10_NCBI_refseq.bed -S *.bw --skipZeros --missingDataAsZero -o 2k_TSS.gz
画profile（--perGroup画在一起加这个参数）
plotProfile -m 2k_TSS.gz -out merge.pdf
画热图
plotHeatmap -m 2k_TSS.gz -out merge.pdf
1\插入片段
samtools view Aging2.last.bam | cut -f 9 >Aging2.txt
用R画图
A2=read.table('Aging2.txt')
dim(a)
hist(abs(as.numeric(a[,1])), breaks=100)
2、FRiP值的计算
raw.bam转raw.bed
bedtools intersect -a Aging1.raw.bed -b Aging1_peaks.narrowPeak |wc -l|awk '{print $1}'
3 、IDR计算
即样本间重复性
