一、deepmod
1、三代测序得到的fast5文件是muti fast5,一个fast5文件里面有4000条fast5序列,deepmod不支持muti fast5,需要拆分成singal fast5. 使用ont官方的程序ont_fast5_api进行转换
https://github.com/nanoporetech/ont_fast5_api
$ multi_to_single_fast5 --input_path fast5 --save_path single_fast5 --recursive
2、需要使用Albacore对fast5文件进行basecalling,写入events文件
1)Albacore可以用docker来运行,https://hub.docker.com/r/robegan21/albacore
软件提供的训练集是用albacore2.3.1完成的,我们也用这个版本进行basecalling
2)read_fast5_basecaller.py -l 来查看不同的试剂盒和芯片对应的config文件,我的芯片是FLO-MIN106,试剂盒是SQK-109,对应的config文件是r94_450bps_linear.cfg
$ docker run -it --rm -v $PWD:/data robegan21/albacore:2.3.1 bash
$ read_fast5_basecaller.py -i single_fast5 -f FLO-MIN106 -k SQK-109 -r -t 10 -c /opt/albacore/r94_450bps_linear.cfg -o fast5,fastq -s albacore
3、Deepmod进行m6A检测
1)Deepmod安装参照官网,
https://github.com/WGLab/DeepMod/blob/master/docs/Install.md
tensorflow、matplotlib直接使用conda安装可能会出问题,建议安装时指定版本,按照安装说明中的版本,tensorflow安装1.7.0,matplotlib安装2.1.2
2)github上下载Deepmod,这个软件是免安装的,下载下来,只要上面的依赖都安装好了,就可以直接用
$ python Deepmod/bin/Deepmod.py detect --wrkBase albacore/workspace/pass --Ref ref/wt.fasta --Base A --FileID wt --modfile Deepmod/train_mod/rnn_conmodA_P100wd21_f7ne1u0_4/mod_train_conmodA_P100wd21_f3ne1u0 --threads 10 --outFolder wt
如果报错的话,加上 --basecall_1d Basecall_1D_001 再试一下
二、mCaller
1、安装参照说明即可,没有遇到什么问题
https://github.com/al-mcintyre/mCaller
2、multi fast5 转成single fast5
https://github.com/nanoporetech/ont_fast5_api
$ multi_to_single_fast5 --input_path fast5 --save_path single_fast5 –recursive
3、输入文件需要带event信息的fast5文件和对应的fastq文件,我们原始的fast5是不带event信息的,需要重新进行basecalling。用Albacore或者guppy均可,albacore用上面deepmod里面说的docker程序运行,guppy的话可以下载安装,用guppy进行basecalling
$ guppy_basecaller -i single_fast5 -s single_fast5 -c dna_r9.4.1_450bps_fast.cfg --fast5_out on --num_callers 4 --cpu_threads_per_caller 3
### --fast5_out 默认是off,这里要选成on,这样就可以输出带event信息的fast5文件了
4、将basecalling得到的fastq文件合并,并建立索引
$ cat *.fastq >> wt.fastq
$ nanopolish index -d workspace -s sequencing_summary.txt wt.fastq
## workspace 是guppy basecalling输出的带event信息的fast5文件的文件夹
5、后面就按部就班来就好了
$ bwa index wt.fasta
$ bwa mem -x ont2d -t 12 wt.fasta wt.fastq | samtools view -Sb - | samtools sort -T /tmp/wt.sorted -o wt.sorted.bam
$ samtools index wt.sorted.bam
$ nanopolish eventalign -t 12 --scale-events -n -r wt.fastq -b wt.sorted.bam -g wt.fasta > wt.eventalign.tsv (不要丢掉wt.fasta后面的>)
$ mCaller.py -m A -r wt.fasta -d r94_model_NN_6_m6A.pkl -e wt.eventalign.tsv -f wt.fastq -b A ( -m 参数可以指定motif,比对GATC,如果想要检测所有的A的话,直接-m A就好了)
然后会得到一个wt.eventalign.diffs.6 的文件,里面有每个A是否发生甲基化修饰的信息了。后面还可以通过调整阈值线来判断哪些是真的甲基化修饰
三、tombo
1、安装
https://github.com/nanoporetech/tombo
使用conda安装,没有问题
2、tombo不支持multi fast5,同样使用ont_fast5_api进行转换
https://github.com/nanoporetech/ont_fast5_api
$ multi_to_single_fast5 --input_path fast5 --save_path single_fast5 –recursive
3、resquiggle 必须使用single_fast5文件夹为输入文件夹
$ tombo resquiggle single_fast5 wt.fasta –processes 4 –num-most-common-errors 5
4、检测6mA。tombo检测碱基修饰有3种方式,推荐使用alternative model
$ tombo detect_modifications alternative_model --fast5-basedirs single_fast5 --alternate-bases 6mA --statistics-file-basename wt_6ma
1)--statistics-file 是指定输出文件的名字
2)--alternate-bases还支持CpG dam dcm检测
得到wt_6ma.6mA.tombo.stats,后续的motif分析以及作图都是基于这个文件。
motif检测
$ tombo text_output signif_sequence_context --statistics-filename wt_6ma.6mA.tombo.stats --genome-fasta wt.fasta --num-regions 1000 --num-bases 20
得到检测到甲基化的A的周围20bp的序列,使用MEME进行motif检测