首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >「翻译」在生物信息学中使用 GNU-Parallel

「翻译」在生物信息学中使用 GNU-Parallel

作者头像
王诗翔呀
发布2020-09-25 18:16:28
1.1K0
发布2020-09-25 18:16:28
举报
文章被收录于专栏:优雅R优雅R

原文出处:https://www.danielecook.com/using-gnu-parallel-for-bioinformatics/

GNU Parallel[1] 是一个用于加速生信分析不可或缺的一个工具。它允许你非常简单地对命令并行化处理。下面我将介绍一些如何使用它以及如何将它应用于生信。

很多高性能计算平台节点已经预先安装了它。你可以从 homebrew[2] 或其他包管理器找到和安装它。

基本用法

让我们从一个简单的例子开始:

seq 1 5 | parallel -j 4 echo

这里我们 (1) 打印了数字 1 到 5,且 (2) 将该序列数据通过管道传进了 parallel 命令。我们提供了一个命令 echo ,它将通过 -j=4 的选项指定进行并行化。我们可以通过添加 --dry-run 打印将要运行的命令。

seq 1 5 | parallel --dry-run -j 4 echo
echo 3
echo 4
echo 5
echo 2
echo 1

结果是乱序的!这是并行化的本质:不是所有的任务都会花费相同的时间,所以有的结束的早,有的结束的晚,因此输出顺序并不一致。我们可以使用 -k 选项强制程序执行“先入先出”准则。让我们看一下下面的结果:

seq 1 5 | parallel -j 4 -k echo
1
2
3
4
5

像其他任何命令一样,我们可以将输出保存到文件中:

seq 1 5 | parallel -j 4 -k echo > out.txt

-j

为了让 GNU Parllel 工作,你需要一个多核 CPU。指定超过所拥有的核心数会让性能变得糟糕。因此,调节 -j 选项以便于命令更好地工作是非常重要的。

幸运地是,parallel 运行你通过 -j 指定计算占有的 CPU 比例或相对数量。

parallel -j 100% # 使用 100% 的核心数
parallel -j -1 # 使用比所有核心数少 1 个的核心数
parallel -j +1 # 使用比所有核心数多 1 个的核心数

使用 ::: 传递参数

使用 ::: 指定并行指定的命令参数(列表来源)。

parallel -j 4 -k echo ::: `seq 1 5`

「注意」,上面这种情况能够传递的参数数量是有限的,通过管道传递参数或像下面一样通过文件传递参数可能更好:

seq 1 5 | parallel -j 4 -k echo

使用 :::: 传递文件内的参数

For large argument lists you can specify a file with a list of arguments. Specify a file of arguments (one per line) using ::::.

如果参数列表很大,你可以通过文件指定,文件每一行对应要并行的一个参数:

parallel -j 4 -k echo :::: my_args.txt

使用 `

默认 parallel 假定参数放在输入命令的结尾,但是你可以通过 {} 显式地指定:

parallel --dry-run -j 4 -k echo \"{} "<-- a number"\" ::: `seq 1 5`
echo "1 <-- a number"
echo "2 <-- a number"
echo "3 <-- a number"
echo "4 <-- a number"
echo "5 <-- a number"

注意上面我们使用转移符号输出了双引号。

组合

你可以组合 :::::: 来添加额外的参数,然后它们会生成所有可能的组合。这对测试有不同参数组合的命令非常有用:

parallel --dry-run -k -j 4 Rscript run_analysis.R {1} {2} ::: `seq 1 2` ::: A B C
Rscript run_analysis.R 1 A
Rscript run_analysis.R 1 B
Rscript run_analysis.R 1 C
Rscript run_analysis.R 2 A
Rscript run_analysis.R 2 B
Rscript run_analysis.R 2 C

并行化函数

在一些情况下,你想要执行一系列的命令。例如,下面的代码计算了输入DNA序列互补配对的碱基计数:

echo "ATTA" |  tr ATCG TAGC | \
    python -c "import sys; o=sys.stdin.read().strip(); print(o, o.count('T'), o.count('G'), o.count('C'), o.count('A'))"

这个命令有 2 个操作。但我们可以将它整合为 'one-liner':创建一个 bash 函数,导出它,然后使用它作为输入:

function count_nts {
    # $1 is the first argument passed to the function
    echo $1 | tr ATCG TAGC | \
    python -c "import sys; o=sys.stdin.read().strip(); print(o, o.count('T'), o.count('G'), o.count('C'), o.count('A'))"
}

# Use the `-f` flag to export functions
export -f count_nts

parallel -j 4 count_nts ::: TAAT TTT AAAAT GCGCAT | tr ' ' '\t'

有了上面的基础,让我们看看如何使用它加速生信分析。

使用 GNU Parallel 进行 Variant Calling

当处理 BAMs 或 VCFs 时,你可以并行处理所有的染色体。大多数变异检测软件或注释工具允许你通过指定区间一次处理一个染色体。这允许我们使用「拆分-应用-组合」策略到该分析中。

从 BAM 中分割染色体

chrom_list=`samtools idxstats in.bam | cut -f 1 | grep -v '*'`

# For c. elegans you can would see the following 7
# I
# II
# III
# IV
# V
# X
# MtDNA

我们可以创建一个 bash 函数:

function bam_chromosomes {
    samtools idxstats $1 | cut -f 1 | grep -v '*'
}

应用操作到每一条染色体

下面是处理代码:

#!/bin/bash

genome=path/to/genome.fa
export genome # This is critical!

function parallel_call {
    bcftools mpileup \
        --fasta-ref ${genome} \
        --regions $2 \
        --output-type u \
        $1 | \
    bcftools call --multiallelic-caller \
                  --variants-only \
                  --output-type u - > ${1/.bam/}.$2.bcf
}

export -f parallel_call

chrom_set=`bam_chromosomes test.bam`
parallel --verbose -j 4 parallel_call sample_A.bam ::: ${chrom_set}

一些重要的注意事项:

  • 你必须导出 export 所有并行化函数中使用到的变量,例如上面的 genome
  • --output-type=u 是出于效率考虑。
  • Finally, I use a variable substitution to remove the extension from the bam and to generate a <sample>.<chromsome>.bcf filename:

组合变异检测结果

一旦我们完成工作,接着我们使用 bash 数组和组合所有结合并将其廉洁为单个文件。

# Generate an array of the resulting files
# to be concatenated.
sample_name="sample_A"
set -- `echo $chrom_set | tr "\n" " "`
set -- "${@/#/${sample_name}.}" && set -- "${@/%/.bcf}"
# This will generate a list of the output files:
# sample_A.I.bcf sample_B.II.bcf etc.

set -- "${@/#/test.}" && set -- "${@/%/.bcf}"

# Output compressed result
bcftools concat $@ --output-type b > $sample_name.bcf

# Remove intermediate files
rm $@

为了确保所有中间文件被删除,这里使用了bash trap[3]。

总结

GNU Parallel 可以极大提高简单并行场景任务处理效率。虽然需要编写额外的代码用于处理拆分和组合两步,但这可以得到极大的效率提升。

Reference

[1]

GNU Parallel: https://www.gnu.org/software/parallel/

[2]

homebrew: https://www.danielecook.com/using-gnu-parallel-for-bioinformatics/brew.sh

[3]

bash trap: http://redsymbol.net/articles/bash-exit-traps/

本文参与 腾讯云自媒体分享计划,分享自微信公众号。
原始发表:2020-09-19,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 优雅R 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体分享计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 基本用法
    • -j
      • 使用 ::: 传递参数
        • 使用 :::: 传递文件内的参数
          • 使用 `
            • 组合
              • 并行化函数
              • 使用 GNU Parallel 进行 Variant Calling
                • 从 BAM 中分割染色体
                  • 应用操作到每一条染色体
                    • 组合变异检测结果
                      • Reference
                  • 总结
                  相关产品与服务
                  高性能计算平台
                  高性能计算平台(TencentCloud High Performance Computing,THPC)是一款腾讯云自研的高性能计算资源管理服务,集成腾讯云上的计算、存储、网络等产品资源,并整合 HPC 专用作业管理调度、集群管理等软件,向用户提供弹性灵活、性能卓越、自助化的计算服务。可以帮助您高效地管理云上高性能计算资源,实现弹性使用云上高性能计算资源的需求。
                  领券
                  问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档