双端测序，我先下载了一个样本的试试

注意，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可视化

bamCoverage --normalizeUsing CPM -b Aging1.rmdup.bam -o /storage2/anlei/MLZ/ATAC-seq/bw/Aging1.bw -p 40

可视化

computeMatrix reference-point --referencePoint TSS -p 50  \

-b 10000 -a 10000    \

-R /storage2/anlei/wangwenjing/ChIP-seq/deeptools/bw/mm10_NCBI_refseq.bed  \

-S /storage2/anlei/MLZ/ATAC-seq/bw/*.bw  \

--skipZeros  -o matrix1_TSS.gz  \

--outFileSortedRegions regions1_test_genes.bed

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计算

即样本间重复性