育种数据分析之放飞自我分享 http://blog.sciencenet.cn/u/yijiaobai 关注:生物统计,数量遗传,混合线性模型,生物信息,R,Perl,Python,GWAS,GS相关方法,文章及代码

博文

snakemake-学习笔记2

已有 3514 次阅读 2019-4-3 09:18 |个人分类:便捷操作|系统分类:科研笔记

参考资料

一个稍微复杂的案例

这是另一个snakemake的案例, 之前介绍过通过简单的方法, 使用snakemake, 这里我们用另一个案例, 看看snakemake的用法.

过程介绍

  • 1, 安装snakemake
  • 2, 新建文件
  • 3, 新建一个简单的Snakemake参数文件
  • 4, 扩展, 去关联输出文件
  • 5, 使用全局变量, 关联文件
  • 6, 批量运行

1, 安装snakemake

这里需要时python3, 不支持python2

pip3 install --user snakemake pyaml

2, 新建几个FASTQ文件

这里, 我们新建两个配对的RNA-seq数据, 格式是FASTQ的文件, 然后经过下面两步处理:

  • 第一步: 数据质量控制
  • 第二部: 将基因表达合并为一个文件

创建文件

  • 创建genome.fa文件, 使用touch创建空文件即可
  • 创建fastq文件夹
  • 在fastq文件夹中, 创建Sample1.R1.fastq.gz Sample1.R2.fastq.gz Sample2.R1.fastq.gz Sample2.R2.fastq.gz四个空文件
touch genome.fa

# Make some fake data:
mkdir fastq
touch fastq/Sample1.R1.fastq.gz fastq/Sample1.R2.fastq.gz
touch fastq/Sample2.R1.fastq.gz fastq/Sample2.R2.fastq.gz

创建结果, 使用tree查看:

(base) [dengfei@localhost test]$ tree
.
├── fastq
│   ├── Sample1.R1.fastq.gz
│   ├── Sample1.R2.fastq.gz
│   ├── Sample2.R1.fastq.gz
│   └── Sample2.R2.fastq.gz
└── genome.fa

1 directory, 5 files

3, 创建snakemake参数文件

将下面代码命名为Snakemake

SAMPLES = ['Sample1', 'Sample2']

rule all:
    input:
        expand('{sample}.txt', sample=SAMPLES)

rule quantify_genes:
    input:
        genome = 'genome.fa',
        r1 = 'fastq/{sample}.R1.fastq.gz',
        r2 = 'fastq/{sample}.R2.fastq.gz'
    output:
        '{sample}.txt'
    shell:
        'echo {input.genome} {input.r1} {input.r2} > {output}'

4, 参数解释

我们下面进行代码的讲解:

    1. 这里, 定义了一个SAMPLE的数组:
      SAMPLES = ['Sample1', 'Sample2']
      数组, SAMPLES,里面有两个元素: Sample1和Sample2
    1. 定义一个rule, 名称为all, input使用expand函数, 能够将数组的内容解析给{sample}
      rule all:
      input:
         expand('{sample}.txt', sample=SAMPLES)
    1. 定义一个rule, 命名为 quantify_genes, 里面有input, output, shell, 其中{sample}是用的rule all里面的name
rule quantify_genes:
    input:
        genome = 'genome.fa',
        r1 = 'fastq/{sample}.R1.fastq.gz',
        r2 = 'fastq/{sample}.R2.fastq.gz'
    output:
        '{sample}.txt'
    shell:
        'echo {input.genome} {input.r1} {input.r2} > {output}'

5, 运行参数

  • 1, 预览命令
    使用命令:
snakemake -np

参数介绍
-n 或者—dryrun, 表示只生成命令, 但是不执行命令, 可以预览一下生成的命令.

  --dryrun, -n          Do not execute anything, and display what would be
                        done. If you have a very large workflow, use --dryrun
                        --quiet to just print a summary of the DAG of jobs.

-p 或者—printshellcmds, 表示将生成的shell打印出来

 --printshellcmds, -p  Print out the shell commands that will be executed.

注意:
-n 不执行, 只打印命令
-p 执行, 同时打印命令(shell)
两者执行的前提是结果文件还没有生成.

例子:

(snake_test) [dengfei@localhost ex2]$ snakemake -np
Building DAG of jobs...
Job counts:
    count    jobs
    1    all
    2    quantify_genes
    3

[Tue Apr  2 13:49:34 2019]
rule quantify_genes:
    input: genome.fa, fastq/Sample1.R1.fastq.gz, fastq/Sample1.R2.fastq.gz
    output: Sample1.txt
    jobid: 1
    wildcards: sample=Sample1

echo genome.fa fastq/Sample1.R1.fastq.gz fastq/Sample1.R2.fastq.gz > Sample1.txt

[Tue Apr  2 13:49:34 2019]
rule quantify_genes:
    input: genome.fa, fastq/Sample2.R1.fastq.gz, fastq/Sample2.R2.fastq.gz
    output: Sample2.txt
    jobid: 2
    wildcards: sample=Sample2

echo genome.fa fastq/Sample2.R1.fastq.gz fastq/Sample2.R2.fastq.gz > Sample2.txt

[Tue Apr  2 13:49:34 2019]
localrule all:
    input: Sample1.txt, Sample2.txt
    jobid: 0

Job counts:
    count    jobs
    1    all
    2    quantify_genes
    3
This was a dry-run (flag -n). The order of jobs does not reflect the order of execution.

相关阅读:

snakemake 学习笔记1

公众号.png



https://wap.sciencenet.cn/blog-2577109-1171194.html

上一篇:snakemake 学习笔记1
下一篇:snakemake 学习笔记3
收藏 IP: 47.75.72.*| 热度|

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-3-29 01:07

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部