随着学习的不断深入,分析的数据越来越多。你会发现,日常生信分析不过是调用一些相同的函数或者包分析不同的数据,换汤不换药。
那么,如何把分析过程流程化,让数据像工厂的流水线一样自动被处理?
最简单的法子,就是用shell命令写个for循环。
进阶还可以设置下多进程(参考梨酱:[Linux 1] Shell“ 多线程”,提高工作效率)。
但是,用shell写流程存在一些问题:
- shell 命令一旦断掉,如何判断中间文件的完整性?
- 下次想要修改流程,重新修改shell 脚本太麻烦!
- 想删掉/保留中间文件,但是代码好复杂?
......
有没有更简捷、少掉头发的办法,把分析的过程“组装”到一起,让电脑“自行”处理数据呢?
那就赶紧用 Snakemake 吧!!!
(Snakemake 翻成中文是“盘蛇”?Σ(っ °Д °;)っ误)
为什么要用 Snakemake ?
Snakemake的优点:
- 支持多线程
- 支持断点运行
- 支持shell命令,还可以和 Python 库结合使用
行,现在我们知道 Snakemake 比 shell 更好用了。
搓手手,那么如何使用 Snakemake 呢?
来试试栗子 ~
0. 安装 snakemake
conda install -c bioconda snakemake -y # -y 表示安装时不询问
1.简单入门,合并两个文件
# 创建两个文本文件
echo "a" > a.txt
echo "b" > b.txt# 创建一个 snakefile 文件
$ touch snakefile# 打开 snakefile 并写入以下代码,然后保存退出
rule test: # 类似Python中的def func,用rule定义一个名为“test”的流程input: # 定义"test"流程中的输入文件expand("{file}.txt", file=["a", "b"]) # snakemake流程中的替换文件的规则,类似于shell中的*.txtoutput: # 定义“test”流程中的输出文件merged.txt # 输出文件名shell: # 表示我们执行的shell命令cat {input} > {output} # 将我们定义的输入文件用cat命令合并输出到output定义的文件# 在shell中运行
$ snakemake
你会发现,Snakemake 定义的 rule 分为 input、output 和 shell 三个部分。
但是,它只支持 shell 命令吗?
不,只要你将 shell 改为 run,也可以使用 Python 或者 R 代码。如下↓
rule test1: # 用rule定义一个名为“test1”的流程input: # 定义"test1"流程中的输入文件expand("{file}.txt", file=["a", "b"]) # snakemake流程中的替换文件的规则,类似于shell中的*.txtoutput: # 定义“test”流程中的输出文件merged.txt # 输出文件名run: # 输入并执行 Python 代码data1 = open("a.txt").read() #读取 a 文件中的所有内容data2 = open("b.txt").read() #读取 b 文件中的所有内容 data = data1 + data2 # 合并到 data 中
你还可以将 run 改为 scripts,执行自定义的 Python 脚本。
2. 将上例延伸——如何合并两个gtf文件?
rule merge_novel_and_known:input: novel="samples/cuffmerge/novel_transcript_gn.gtf", known=GTFoutput: "samples/new_annotation/all_transcripts.gtf"params: gtf=GTF, fa=FASTAmessage: "--- Merging known and novel transcripts."shell: """cat {input.novel} {input.known} > {output}"""
你能说出每一行代码的意思吗?
提示:
params 指定程序运行的参数
message 运行该rule时,终端给出的提示信息
其实,每一个 rule 就像是一个代码模块积木。
我们利用 Snakemake,先搭建基础的rule积木,接下来只要按顺序把它们放在一起,就可以拼凑出完整的流程!
之后在处理其他类似数据时,只需要略微更改目录等信息,便可复用流程。
何乐而不为?
3. 进阶,如何用 Snakemake 写 Chip-Seq 分析流程?
提示:
1. [软件使用 3] 使用MACS2分析ChIP-seq数据,快速入门!
2. Writing a RNA-Seq workflow with snakemake. Denis Puthier & Aitor Gonzalez. 2 Nov 2015.
如果有任何问题或者建议,可以在评论中留言~
疫情期间,希望大家都能保护好自己!
在家也可以好好学习写代码呀~
ヾ(◍°∇°◍)ノ゙
:
梨酱:[R 01] 不要一直用ggplot2啦,尝试用ggpubr画图吧!zhuanlan.zhihu.com梨酱:[Python 3] 5分钟学会3种方法给模块添加路径!zhuanlan.zhihu.com梨酱:[论文写作 1] 如何用word批量制作三线表?zhuanlan.zhihu.com参考:
- Snakemake Tutorial.
- Writing a RNA-Seq workflow with snakemake. Denis Puthier & Aitor Gonzalez. 2 Nov 2015.
- snakemake使用笔记. 井底蛙蛙呱呱呱. 简书. 2018.
- snakemake--我最喜欢的流程管理工具. 徐洲更hoptop. 简书. 2018.
- 生物信息学100个基础问题 —— 番外5:使用Snakemake快速搭建生信分析流程. 孟浩巍. 知乎. 2019.