近期接触到一个新的测序类型 -- LACE-seq (linear amplification of complementary DNA ends and sequencing)。也许，我们对ChIP-seq并不陌生，知道它可以用来研究某个转录因子或者组蛋白在基因组范围的结合位点。那么，LACE-seq技术可以类比ChIP-seq，用来研究RNA结合蛋白 (RNA-binding protein, RBP) 在基因组范围内所绑定的RNA靶标。
目前，鉴定RBP靶标常用方法主要有RIP-seq和CLIP-seq，但是由于这两种方法均依赖于特异性抗体富集RBP，且需要百万数量级的细胞来制备文库。因此，限制了这些方法在稀有细胞类型及临床穿刺样本中的应用。针对这些缺点，LACE-seq实现了在微量细胞中鉴定RBP作用靶点，通过线性扩增逆转录酶在RBP结合位点处的终止信号，实现了在单碱基分辨率和单细胞层面精准鉴定RBP的结合位点。
研究表明，RBP在早期生殖、个体发育、细胞分化、增殖和凋亡等生理过程中发挥了关键的调控作用，而LACE-seq技术则为研究RBP在胚胎发育和生殖疾病中的功能机制奠定了一些基础。
从上面的建库示意图可知，在免疫沉淀后，使用含有poly-T的序列来做PCR扩增第二条链，这说明该文库富集了mRNA。另外，序列中引入了两个额外的碱基，这也是后续分析需要注意的事项。
软件准备
cutadapt去接头，bowtie比对，piranha用于peakcaling。前面两个软件比较常见，这里就不多说了。主要来说piranha，该软件安装需要依赖系统库GSL (GNU Scientific Library)，没有的话需要提前安装：
wget http://mirrors.ustc.edu.cn/gnu/gsl/gsl-2.7.tar.gz
tar -zxf gsl-2.7.tar.gz
cd gsl-2.7
./configure --prefix=/path/gsl-2.7
make
make check
make install
安装后，需要将其添加到环境变量里面，为了使用以后方便使用可以将其添加到home目录下的.bashrc文件里面：
export PATH=$PATH:/path/gsl-2.7/bin
export C_INCLUDE_PATH=$C_INCLUDE_PATH:/path/gsl-2.7/include
export CPLUS_INCLUDE_PATH=$CPLUS_INCLUDE_PATH:/path/gsl-2.7/include
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/path/gsl-2.7/lib
export LIBRARY_PATH=$LIBRARY_PATH:/path/gsl-2.7/lib
添加后，使用命令source ~/.bashrc立刻加载修改后的设置，使之立即生效。然后就可以安装软件piranha软件了：
wget http://smithlabresearch.org/downloads/piranha-1.2.1.tar.gz
tar -xzf piranha-1.2.1.tar.gz
./configure
make all
make install
上面安装的piranha输入文件格式只能是bed，如果想直接接受bam文件，需要先安装BAMTools，安装方法可以参考：https://github.com/pezmaster31/bamtools/wiki，接受bam格式的安装方法如下：
./configure --with-bam_tools_headers="/path/BAMTools/include" --with-bam_tools_library="/path/BAMTools/lib"
make all
make install
为了省事起见，也可以选择使用conda来安装piranha：
conda create -n piranha -c bioconda piranha
不过，虽然用conda安装软件很省事，但安装的piranha也不接受bam格式的输入文件。
分析流程
1、去接头
先去除序列中的接头，然后再去除poly-A尾巴：
cutadapt -f fastq -q 30,0 -a ATCTCGTATGCCGTCTTCTGCTT -m 18 --max-n 0.25 --trim-n -o sample.rmadp.fq.gz sample.fastq.gz
cutadapt -f fastq -a A{15} -m 18 -n 2 -o sample.trimed.fq.gz sample.rmadp.fq.gz
2、比对
先比对到pre-RNA (NCBI网站可以下载) 上面去除RNA，将未必对上的序列再比对到基因组上。比对设定的参数允许两个错配的作用，应该就是矫正建库时引入的两个额外碱基：
bowtie -p 5 -v 2 -m 10 --best -strata --un sample.unmap.fastq -S rnaidx sample.trimed.fq.gz sample.rna.sam
bowtie -p 5 -v 2 -k 10 --best -strata --un -S refidx sample.unmap.fastq sample.sam
samtools view -@ 5 -Sb -o sample.bam sample.sam
samtools sort -O BAM -@ 5 -o sample.sort.bam sample.bam
samtools index -@ 5 sample.sort.bam
3、peakcalling
最后，使用piranha软件来做callpeak，如果输入文件bed，可以用bedtools软件的bamtobed将bam转为bed格式：
# bed输入
bedtools bamtobed -i sample.sort.bam >sample.bed
Piranha -s -p 0.001 -b 20 -d ZeroTruncatedNegativeBinomial -o sample.peak.bed sample.bed
# bam输入
Piranha -s -p 0.001 -b 20 -d ZeroTruncatedNegativeBinomial -o sample.peak.bed sample.sort.bam
过程还是挺简单的，至此，我们就得到了富集的peak文件，基于此就可以做后续的分析了。有点奇怪的地方，按照软件的说明参数-s指定是否需要排序，我这里使用的是sort后的bam文件，还是需要这个参数，否则会报错。也许，要求的不是coordinate排序方式。
结果展示
下面展示一下，本人使用文献数据分析的部分结果，数据来自GSE137925。
meta-analysis of Ddx4 LACE-seq：
IGV track of Chr17(qB1)：
IGV track of D10Wsu102e and Pou5f1：
从上面的结果来看，基本重现了文献的结果，只是本人做的图显得有些潦草了，不过，结果可以跟文献吻合说明咱们的数据处理过程没有错误，这才是最重要的！
参考资料
http://smithlabresearch.org/software/piranha
<<Global profiling of RNA-binding protein target sites by LACE-seq>>
往期回顾
可视化：两种方法绘制桑基图
可视化：density与ridgeline
单细胞：不同亚群的基因平均表达热图
可视化：bubble
可视化：scatterplot
