iSA：DNA甲基化新的分析方法（part 1）

19年的nature protocol的文章，非常值得有学习的地方！
文章：Resolution of the DNA methylation state ofsingle CpG dyads using in silico strand annealingand WGBS data
原文地址：https://www.nature.com/articles/s41596-018-0090-x
分析重点：甲基化状态，半甲基化状态，为甲基化状态，CpG二聚体

一：背景介绍：

  1.1） WGBS（Whole Genome bisulfite sequencing）是被广泛的用于组织和细胞中胞嘧啶（C）甲基化的定量。在经历了BS转化，它可以做到单碱基的分辨率。
  1.2 ）不同细胞类型的DNA甲基化组中有4 - 20%为半甲基化，其中一些可以在细胞分裂过程中遗传，这表明半甲基化是一种稳定的表观遗传标记。

PS：什么是DNA半甲基化？
  半甲基化是其中一条链添加一个甲基，而另一条链没有加上甲基基团。这种半甲基化的现象很可能是由于观察的时间点位于细胞分裂后和子链完全复制之前，被认为是甲基化和未甲基化的过渡状态。还有一些观点认为，半甲基化是甲基化过程中的随机错误。

1.3）这样有着几个情况：
1.3.1，我们在统计CpG二聚体甲基化状态的时候，只能是单独计算计算CpG的甲基化状态，然后把结果重合起来看。CpG二聚体甲基化状态在DNA的复制过程很重要。在DNA复制传代的过程中，父链和子链的CpG二聚体的甲基化状态是不清楚的。
1.3.2，我们在计算DNA半甲基化（hemi-methylated CpGs）的时候，没法定量。

二：iSA介绍：

  iSA（in silico strand annealing ） 这种是一个可基于WGBS数据，新的分析DNA甲基化的生信分析方法。它可以分析不同甲基化状态的CpG二聚体。并且可以运用到nasBS-seq，ChIP- bs -seq，TAB-seq, oxBS-seq，和fCAB-seq。
  iSA可以通过比对到同一区域的正链和负链，然后可以计算出单个CpG二聚体的甲基化状态。 将iSA与nasBS-seq和WGBS相结合，我们已证实半甲基化DNA和CTCF可共定位，作用于染色体的折叠。而且已经证明半甲基化是所有检测细胞类型中DNA甲基化组的重要组成部分，作为一种稳定的表观遗传的标致。
  和其他的分析方法比较(现在只有hairpin-bisulfite sequencing)，运用到了全基因组的数据分析。但是hairpin-bisulfite sequencing要求高测序深度。

ps:给一个hairpin-bisulfite sequencing 的 workflow

hairpin bisulfite sequencing workflow1.png
hairpin bisulfite sequencing workflow2.png
文章来自于：https://www.ncbi.nlm.nih.gov/pubmed/29224164

iSA可应用于现有WGBSdatasets的下游，将数据计算分解为单cpg甲基化。我们来看一下iSA的overview

overview of iSA

拆解步骤：
1）首先利用bismark进行mapping
2）去重复
3）用samtools或者是bedtools去不断的循环搜索（目的是为了找到正负链两端位置完全相同的基因组坐标的reads用于下游分析。 需要注意的是从细胞群中随机剪切大量染色质可能产生具有相同末端但具有不同dsDNA来源的随机的正负链。）
4）为了确定这种随机配对可能产生的干扰，通过bedtools比对多次比对的结果。
（在一个具有足够基因组覆盖和测序深度的WGBS数据中，我们通常观察到和随机配对相比20- 100倍的同端富集，表明随机配对的干扰非常小）

富集情况，通过上下游滑动10bp 20bp 30bp看结果

5）使用Samtools、Bismark和Bedtools从同端成对比对中提取DNA甲基化信息，在CpG的DNA上下reads中，一条链上的每个C都有另一条链上的反向的C对应。
6）对于CpG二聚体来说存在以下4种情况，未甲基化，半甲基化（正or负链全部甲基化），全甲基化。然后再是下游个性化分析。

Tip:
1）ISA不支持单端测序的结果分析
2）iSA不推荐分析加完adapter之后再转化的DNA文库，为什么呢？

因为在加完adapter之后，进行BS转化会对adapter和DNA序列造成一定的损伤，而且DNA的双链系列会发生改变，在mapping的时候基因组的位置会有改变。导致分析的时候不能对齐。

3）bismark读入的数据必须是进行数据过滤后具有最高质量的数据。可以用比较高的cutoff去卡这个标准。但是如果把低质量的碱基去掉了以后会对比对的结果造成影响

为什么呢？会导致其基因的坐标位置发生改变。而且在同样序列的情况下，把很多低质量的reads去掉了会造成后续去重的问题。

3：代码实战

3/1 首先下载需要分析的工具和数据

SRA  Toolkit  (http://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/2.9.2/sratoolkit.2.9.2-ubuntu64.tar.gz)

Trimmomatic(http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.38.zip)

FastQC(http://www.bioinformatics.babraham.ac.uk/projects/fastqc/fastqc_v0.11.7.zip)

Bismark(http://www.bioinformatics.babraham.ac.uk/projects/bismark/bismark_v0.20.0.tar.gz)

Bowtie 2 (ref.17,https://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.2.9/)

Bedtools15(https://github.com/arq5x/bedtools2/releases/download/v2.27.0/bedtools-2.27.0.tar.gz)

Samtools14(https://github.com/samtools/samtools/releases/download/1.9/samtools-1.9.tar.bz2)

deepTools(https://github.com/deeptools/deepTools/archive/2.5.4.tar.gz), for post-iSA visualizationof data

Example WGBS dataset (GSM1386021)
Mouse reference genome sequencefiles (Equipment setup)
lambdareference genome sequencefile (Equipment setup)