（一）介绍

本文档描述了如何利用 Sentieon®基因组学工具的分片能力将 DNAseq®流程分布到多台服务器上；将其他流程（如 TNseq®）进行分布遵循相同原则，因为所有 Sentieon®基因组学工具都具有相同的内置分布式处理能力。这种分布的目标是为了减少流程的总运行时间，以更快地生成结果；然而，这种分布也会带来一些额外的开销，使计算成本增加。

利用分布能力，流程的每个阶段被分成小任务；每个任务处理基因组的一部分，并可以在不同的服务器上并行运行。每个任务生成一个部分结果，需要按顺序合并为最终的单一输出；这种合并需要仔细进行，以确保考虑到边界并生成与没有分片运行的流程相同的结果。

分布的执行框架不在本文档的范围内，用户需要在保持正确的数据依赖关系的同时，分发数据/文件并启动正确的进程。

（二）分片和分片化

我们将基因组分成许多连续且不重叠的部分，每个部分称为一个分片（shard）。每个分片可以定义为单个染色体 contig 的名称，或者是遵循contig:shard_start-shard_end约定的染色体的一部分。特殊的分片 NO_COOR 用于所有没有坐标的未映射读取。

Sentieon®二进制文件支持将分片分布到多个服务器，并且可以通过添加一个或多个带参数的分片选项在单个命令中处理多个分片。在单个命令行中使用多个选项时，这些分片需要按照参考染色体列表是连续的；例如，一个命令可以包含一个覆盖 chr2 结尾的分片和一个覆盖 chr3 开头的分片，但不能同时包含一个覆盖 chr2 结尾和 chr22 开头的分片。

--shard SHARD --shard

您可以参考附录中的脚本，根据与参考相关联的字典或输入 bam 文件中的 bam 头，为基因组创建分片的示例文件。我们建议使用 100M 碱基作为分片大小。

generate_shards.sh

（三）分布式处理的数据流和数据/文件依赖关系

按照推荐的工作流程，DNAseq®的流程将一对 fastq 文件通过以下阶段进行处理：BWA 对齐生成 sorted.bam，去重生成 dedup.bam，BQSR 生成 recal.table，Haplotyper 生成 output.vcf.gz 文件。如图 2 所示，说明了这样一个流程的数据流。

 图 1：DNAseq®流程的典型数据流程

为了将上述流程分布到多个服务器上，每个阶段被划分为处理分片数据的命令，这些命令需要来自特定分片以及相邻分片的输入；但是，有两个例外情况：Dedup 阶段需要来自所有分片上的 LocusCollector 命令的所有 score 文件，而 Haplotyper 阶段需要一个完整的合并校正表。

图 2：以 4 个分片进行分布式处理的 DNAseq®流程的数据流程，不处理未映射的读取。

以图 3 为例，说明了一个以 4 个分片进行分布的流水线的数据流程；这并不是一个典型的使用案例，因为使用推荐的分片大小为 100M 碱基会导致需要超过 30 个分片。在图 3 的示例中，各个阶段需要以下输入并生成以下输出：

· 分片的 LocusCollector 阶段（去重 1）需要 sorted.bam 作为输入。该阶段生成一个文件。

i-th part_deduped.bam$shard_i.score.gz

· 分片的 Dedup 阶段（去重 2）需要 sorted.bam 输入，以及所有 LocusCollector 阶段的所有结果。该阶段生成一个文件。

i-th part_deduped$shard_i.bam

· 分片的 QualCal 阶段需要文件，以及可用的文件和文件。该阶段生成一个文件。

i-th part_deduped$shard_i.bam part_deduped$shard_i+1.bam part_deduped$shard_i-1.bam part_recal_data$shard_i.table

· 在 QualCal 之后，所有文件需要合并为一个单一的校准表文件，用于变异调用。驱动程序和 QualCal 的选项可以用于执行边界感知合并。

part_recal_data$shard_i.table --passthru --merge

· 分片的 Haplotyper 阶段需要文件，以及可用的文件和文件；此外，完整合并的校准表需要作为输入。该阶段生成一个文件。

i-th part_deduped$shard_i.bam part_deduped$shard_i+1.bam part_deduped$shard_i-1.bam part_output$shard_i.vcf.gz

· 在 Haplotyper 之后，所有文件需要合并为一个单一的输出 VCF 文件。驱动程序和 Haplotyper 的选项用于执行边界感知合并。

part_output$shard_i.vcf.gz --passthru --merge

· 如果在流程的任何阶段需要合并的输出 bam 文件，可以使用 util 二进制文件执行边界感知合并；需要在命令中添加选项，以便 util merge 不处理读取，而只是按块复制它们。

--mergemode=10

重要提示：

在合并结果时，务必按照正确的顺序依次输入部分结果，否则结果将不正确。

（四）分发 BWA

上述说明中未包含有关使用 BWA 进行分发对齐的任何信息。为了分发 BWA 对齐，您可以使用 Sentieon®工具中提供的工具为输入的 FASTQ 文件创建索引文件；然后，您可以使用 fqidx 命令的结果作为 BWA mem 的输入，在不同的服务器上处理 FASTQ 文件的特定部分；您需要确保在 BWA 命令中包含该选项，因为 fqidx 的输出包含交错的 reads 在单个输出中。这个流程的结果与在单次运行上执行的结果是相同的。

fqidx -p

BWA_K_size=100000000num_splits=10SAMPLE="sample_name" #Sample name SM tag in bamGROUP="read_group_name" #Read Group name RGID tag in bamPLATFORM="ILLUMINA"NT=$(nproc) #number of threads to use in computation, set to all threads availableFASTA="/home/regression/references/b37/human_g1k_v37_decoy.fasta"FASTQ_1="WGS-30x_1.fastq.gz"FASTQ_2="WGS-30x_2.fastq.gz"FASTQ_INDEX="WGS-30x.fastq.gz.index"#create FASTQ indicessentieon fqidx index -K $BWA_K_size -o $FASTQ_INDEX $FASTQ_1 $FASTQ_2#get the number of runs that the inputs will be split into given the sizenum_K=$(sentieon fqidx query -i $FASTQ_INDEX | cut -d' ' -f1)BWA_K_size=$(sentieon fqidx query -i $FASTQ_INDEX | cut -d' ' -f2)num_K_splits=$(expr $num_K / $num_splits + 1)#run multiple BWA on multiple serversfile_list=""for run in $(seq 0 $((num_splits-1))); do        region="$((run*num_K_splits))-$((run*num_K_splits+num_K_splits))"        #send each of these to a different server        sentieon bwa mem -R "@RG\tID:$GROUP\tSM:$SAMPLE\tPL:$PLATFORM" \          -K $BWA_K_size -t $NT -p $FASTA \          "<sentieon fqidx extract -i $FASTQ_INDEX -r $region $FASTQ_1 $FASTQ_2" | \          sentieon util sort -r $FASTA -t $NT --sam2bam -o sorted_run$run.bam -i -        file_list="$file_list sorted_run$run.bam"done#merge the resultssentieon util merge -o sorted.bam $file_list#or perform deduplication on all the intermediate BAM filesfile_list_bam=""for file in $file_list; do file_list_bam="$file_list_bam -i $file"; donesentieon driver -r $FASTA -t $NT $file_list --algo LocusCollector --fun score_info score.txt.gzsentieon driver -r $FASTA -t $NT $file_list --algo Dedup --score_info score.txt.gz deduped.bam

1.当无法创建 FASTQ 索引文件（版本 201808.02+）时进行 BWA 分发

如果无法创建 FASTQ 索引文件，则可以使用带有分数选项的 util fqidx 工具，并同时使用该选项。此方法将输入的 FASTQ 文件分割成多个读取块的片段，并从每个片段中提取每个元素以在不同的服务器上进行处理，从起始位置到结束位置。-F m/n -K n m-th m 0 n-1

请注意，此功能仅建议在文件存储速度足够快以支持每个 fqidx 进程同时读取输入 FASTQ 文件的 IT 环境中使用。这在云环境中很常见，或者在具有高带宽 NFS 系统的本地集群环境中使用。

BWA_K_size=100000000num_splits=10SAMPLE="sample_name" #Sample name SM tag in bamGROUP="read_group_name" #Read Group name RGID tag in bamPLATFORM="ILLUMINA"NT=$(nproc) #number of threads to use in computation, set to all threads availableFASTA="/home/regression/references/b37/human_g1k_v37_decoy.fasta"FASTQ_1="WGS-30x_1.fastq.gz"FASTQ_2="WGS-30x_2.fastq.gz"#run multiple BWA on multiple serversfile_list=""for run in $(seq 0 $((num_splits-1))); do        #send each of these to a different server        sentieon bwa mem -R "@RG\tID:$GROUP\tSM:$SAMPLE\tPL:$PLATFORM" \          -K $BWA_K_size -t $NT -p $FASTA \          "<sentieon fqidx extract -F $((run))/$num_splits -K $BWA_K_size $FASTQ_1 $FASTQ_2" | \          sentieon util sort -r $FASTA -t $NT --sam2bam -o sorted_run$run.bam -i -        file_list="$file_list sorted_run$run.bam"done#merge the resultssentieon util merge -o sorted.bam $file_list#or perform deduplication on all the intermediate BAM filesfile_list_bam=""for file in $file_list; do file_list_bam="$file_list_bam -i $file"; donesentieon driver -r $FASTA -t $NT $file_list --algo LocusCollector --fun score_info score.txt.gzsentieon driver -r $FASTA -t $NT $file_list --algo Dedup --score_info score.txt.gz deduped.bam

（五）为了进行大规模队列联合调用，需要对 GVCFtyper 进行分布式处理。

针对包含超过 1000 个样本的大规模联合调用，我们推荐在 GVCFtyper 中进行设置。虽然默认的基因分型模式在较小样本集中理论上更准确，但多项式模式在大样本集中同样准确，并且在大量样本中的扩展性更好。

--genotype_mode multinomial

在运行大量样本（>3000）的联合队列调用时，建议使用类似上述描述的分发技术进行分发；然而，需要考虑几个挑战：

· 联合调用的总运行时间。

· 分布式 GVCFtyper 命令将在哪些机器上运行时的资源需求。

· 存储和访问大量 GVCF 输入文件的物流问题。

· 输出 VCF 文件的文件大小。

一般推荐使用 Sentieon®的 GVCFtyper 的分片功能，将不同的基因组片段分别处理在不同的机器上。以下命令访问完整的 GVCF 输入列表，并假设这些输入文件位于共享位置，如 NFS 存储或可通过 s3 http/https 协议进行访问的远程对象存储位置。

#individual shard calling in machine 1sentieon driver -r $FASTA -t $NT --shard $shard1 --shard $shard2 \  --algo GVCFtyper $gvcf_argument GVCFtyper-shard_1.vcf.gz#individual shard calling in machine 2sentieon driver -r $FASTA -t $NT --shard $shard3 \  --algo GVCFtyper $gvcf_argument GVCFtyper-shard_2.vcf.gz…

重要提示：

在处理每个片段时，输入 GVCF 的列表需要保持一致的顺序，因为最终输出中的样本顺序将取决于输入顺序，并且合并需要所有部分文件具有相同的样本顺序。

使用--shard选项时，GVCFtyper 生成的输出 VCF 文件不是有效的 VCF 文件，因此在下面描述的合并之前不应使用它们。

在处理完所有片段后，您需要运行 GVCFtyper 的合并命令来合并结果，确保中间 VCF 输入的顺序与参考基因组一致。输入文件需要在共享位置（例如 NFS 存储或可通过 s3 http/https 协议访问的远程对象存储位置）中可用。

#merge stepsentieon driver --passthru -t $nt --algo GVCFtyper --merge joint_final.vcf.gz \  GVCFtyper-shard_1.vcf.gz GVCFtyper-shard_2.vcf.gz …

1.总体运行时间和资源要求

为了减少整体运行时间，建议使用足够的分片，允许将运行时细粒度分布为多个服务器。分片可以按分片大小或预期创建复杂性。但是，单个分片的最终合并步骤结果无法分发，需要在单个服务器中运行;此事实设置了用于 分布，因为合并可以主导整个运行时。

GVCFtyper 算法是一种非常有效的算法，可能会 I/O 受 GVCF 存储位置性能的限制 输入。因此，建议将 GVCF 输入存储在文件中 提供 600 MBps 传输速率的系统。

对于极大量样本（>10000） 联合队列调用，内存 可能会成为一个问题，某些操作系统限制也可能发挥作用。为 在这种情况下，建议执行以下操作：

· 将操作系统打开文件限制设置为足够大的数量，以允许软件以打开足够的文件句柄。这是通过命令完成的。ulimit -n NUM

· 将操作系统堆栈限制设置为足够大的数量，以允许软件为操作分配足够的内存，因为内存与样本数，并且可能太大而无法容纳在默认情况下分配的堆栈中。这是通过命令完成的。

ulimit -s NUM

·使用包含输入 GVCF 列表的文件，以防止命令过长 参数列出操作系统长度。您可以使用以下命令执行此操作：

sentieon driver ... --algo GVCFtyper output.vcf.gz - < input_files.txt

· 按照 jemalloc 应用说明中所述使用 jemalloc 内存分配，https://support.sentieon.com/appnotes/jemalloc/并通过以下方式更新 vcf 缓存大小 将以下代码添加到脚本中：

export VCFCACHE_BLOCKSIZE=4096

· 将分片大小设置为较小的数字，即 50MBases。此外，使用 GVCFtyper 命令中的驱动程序选项，以确保 所有线程都得到充分利用。

该命令将变为：--traverse_param

sentieon driver -r $FASTA -t $NT --shard $shard1 --shard $shard2 \ --traverse_param 10000/200 \ --algo GVCFtyper GVCFtyper-shard_1.vcf.gz - < input_files.txt

2.大型输出 VCF 文件的挑战

当运行非常大的队列时，输出 VCF 文件将包含大量列，其中包含每个样本的基因型信息。这么多列可能使输出文件变得难以处理，例如在完整文件上运行 VQSR 可能不切实际，因此您可能需要考虑替代 VCF 存储输出的方法。

根据您用于存储输出的方法，将输出按样本组或基因组坐标分割可能会改善大型输出文件的问题；请随时与 Sentieon 支持团队联系，告知您如何存储联合调用的输出的具体要求。

（1）按基因组区域分割输出

为了使输出 VCF 文件变小，您可以在特定的基因组子区域（例如单个染色体）上执行片段合并。您可以通过仅合并子集的中间 VCF 文件来实现此目的。例如，您可以通过仅合并涵盖每个染色体的片段来创建每个染色体的 VCF 文件：如果片段 1-4 涵盖 chr1，而片段 5 同时涵盖 chr1 和 chr2，则以下代码将创建一个仅包含 chr1 变异体的 VCF 文件：

#merge the necessary intermediate sharded VCFssentieon driver --passthru -t $nt --algo GVCFtyper --merge \  GVCFtyper_chr1_tmp.vcf.gz \  GVCFtyper-shard_1.vcf.gz GVCFtyper-shard_2.vcf.gz GVCFtyper-shard_3.vcf.gz \  GVCFtyper-shard_4.vcf.gz GVCFtyper-shard_5.vcf.gz#remove variants not in the proper contig and index the VCFbcftools view GVCFtyper_chr1_tmp.vcf.gz  -r chr1 \  -o - | sentieon util vcfconvert - GVCFtyper_chr1.vcf.gz

使用上述方法，您可以创建有效的 VCF 文件，其大小仅为完整 VCF 文件的一小部分。

为了在上述输出上运行 VQSR，您需要提供一个包含完整基因组范围内所有变异记录的 VCF 文件，因为这些记录需要用于对 VQSR 模型进行适当的校准。然而，VQSR 仅需要 VCF 数据的前 8 列，因此您无需将每个特定基因组子区域的所有 VCF 文件连接起来，可以通过提取和连接每个文件的前 8 列来创建一个包含必要信息的较小的 VCF 文件。使用下面的代码将创建所需的文件：

VarCal VarCal

vcf_list=(GVCFtyper_chr1.vcf.gz GVCFtyper_chr2.vcf.gz) # include more VCFs if applicable#extract the first 8 columns from each region to a new compressed VCFmkfifo tmp.fifosentieon util vcfconvert tmp.fifo GVCFtyper_annotations.vcf.gz &convert_pid=$!exec 3>tmp.fifo #a file descriptor to hold the fifo openbcftools view -h ${vcf_list[0]} | grep "^##" > tmp.fifobcftools view -h ${vcf_list[0]} | grep "^#CHROM" | cut -f 1-8 > tmp.fifofor vcf in ${vcf_list[@]}; do  bcftools view -H $vcf | cut -f 1-8 > tmp.fifodoneexec 3>&-wait $convert_pidrm tmp.fifo

接下来，可以在包含 VCF 的前八列的合并 VCF 上运行算法，生成用于 SNP 和 INDEL 的分段和重新校准文件。然后，可以直接将 VQSR 应用于按基因组区域分割的 VCF，使用算法进行操作：

VarCal ApplyVarCal

vcf_list=(GVCFtyper_chr1.vcf.gz GVCFtyper_chr2.vcf.gz) # include more VCFs if applicable#apply variant quality score recalibrationfor vcf in ${vcf_list[@]}; do  out_vcf=${vcf%%.vcf.gz}_snp-indel-recal.vcf.gz  sentieon driver -r $FASTA -t $nt --algo ApplyVarCal -v $vcf \    --vqsr_model var_type=SNP,tranches_file=${snp_tranches},sensitivity=99.0,recal=${snp_recal} \    --vqsr_model var_type=INDEL,tranches_file=${indel_tranches},sensitivity=99.0,recal=${indel_recal} \    $out_vcfdone

（2）按样本组分割输出

另外，Sentieon® GVCFtyper 合并工具提供了算法选项作为解决这一挑战的潜在方案。在合并步骤中使用算法选项可以生成有效的部分 VCF 文件，每个文件包含一部分样本，从而将完整的 VCF 文件分割成较小、更易处理的输出文件。使用算法选项的用法为：

--split_by_sample --split_by_sample --split_by_sample

sentieon driver --passthru -t $nt --algo GVCFtyper --merge \  --split_by_sample split.conf GVCFtyper_main.vcf.gz \  GVCFtyper-shard_1.vcf.gz GVCFtyper-shard_2.vcf.gz …

算法选项中使用的输入文件是一个以制表符分隔的文件，每行以特定样本组的输出文件开头，后跟以制表符分隔的相应样本组样本列表。您可以使用多行具有相同输出文件的方式，将多个行中的所有样本分组。下面的示例显示了两个样本组，其中样本 1 到 5 将输出到 group1 输出文件，而样本 6 到 8 将输出到 group2 输出文件：

split.conf --split_by_sample

GVCFtyper_file_group1.vcf.gz Sample1 Sample2 Sample3GVCFtyper_file_group1.vcf.gz Sample4 Sample5GVCFtyper_file_group2.vcf.gz Sample6 Sample7 Sample8

上述 GVCFtyper 合并命令将生成以下输出：

· GVCFtyper_main.vcf.gz：一个包含所有 VCF 信息，包括 INFO 列在内的部分 VCF 文件。不包含样本信息。此文件对于运行 VQSR 非常有用，因为它包含了变异校准所需的所有必要信息。

· GVCFtyper_file_group1.vcf.gz：一个包含至 INFO 列、FORMAT 列以及 group1 中所有样本的所有列的部分 VCF 文件。

· GVCFtyper_file_group2.vcf.gz：一个包含至 INFO 列、FORMAT 列以及 group2 中所有样本的所有列的部分 VCF 文件。

然后，您可以使用 bcftools 合并部分 VCF 并选择您感兴趣的样本。您可以使用以下代码来实现：

bash extract.sh GVCFtyper_main.vcf.gz \  split.conf Samle1,Sample4,Sample7

在该脚本中，第三个参数是一个逗号分隔的感兴趣样本列表，而 extract.sh 脚本如下所示：

#!/bin/bashMAIN=$1GRPS_CONF=$2SAMPLES=$3REGIONS=$4BCF=/home/release/other_tools/bcftools-1.3/bcftoolsif [ $# -eq 4 ] || [ $# -eq 3 ]; then  if [ $REGIONS == "" ]; then    BED_ARG=""  else    BED_ARG="-R $REGIONS"  fi  #parse input files from group config file  GRPS=$(grep -v "^#" $GRPS_CONF| cut -f1 | sort | uniq)  $BCF view -h $MAIN | grep -v '^#CHROM' | grep -v '^##bcftools' > out.vcf  hdr=$($BCF view -h $MAIN | grep '^#CHROM')  hdr="$hdr\tFORMAT"  arg="<($BCF view $BED_ARG -H -I $MAIN | cut -f -8)"  col=9  for g in $GRPS; do    s=$($BCF view --force-samples -s $SAMPLES -h $g 2>/dev/null | grep '^#CHROM' | \       cut -f 10-)    [ -z "$s" ] && continue    hdr="$hdr\t$s"    c="<($BCF view --force-samples -s $SAMPLES $BED_ARG -H -I $g 2>/dev/null | \      cut -f $col-)"    arg="$arg $c"    col=10  doneecho -e "$hdr" >> out.vcfeval paste "$arg" >> out.vcfelse  echo "usage $0 main_vcf_file group_config_file csv_sample_list [bed_file]"fi

这里提供的解决方案可以轻松地在 GVCFtyper_main.vcf.gz 上运行 VQSR，以创建经过校准后的文件，然后您可以在 extract.sh 脚本中使用该文件来获取 VQSR 之后的结果。

3.云环境的特殊考虑事项

在云环境中，通常没有能够容纳大量 GVCF 输入的 NFS 存储。对于联合基因分型，Sentieon® GVCFtyper 允许托管在对象存储位置（如 AWS S3）或通过 HTTP 访问的 GVCF 输入文件。然而，对于大型队列（100+），使用对象存储输入的 GVCFtyper 命令的内存需求可能会过高，因此不建议使用此方法访问输入文件。

在云环境中推荐的方法是根据节点处理的分片，将部分 GVCF 输入文件下载到计算节点上。可以使用 bcftools 进行 GVCF 输入的部分下载，但是需要在 bcftools 命令中添加--no-version 选项，以确保不同分片的标头不会因差异而导致 GVCFtyper 拒绝合并它们。

#use bcftools to download shard1 GVCF partial inputs and create the indexbcftools view --no-version -r $shard1_csv_intervals -o - $URL_sample.g.vcf.gz | \  sentieon util vcfconvert - ${sample}_s1.g.vcf.gz#or download multiple files in parallel using xargs using a listcat list_inputs.txt | xargs -P $NUM_PARALLEL_DOWNLOAD -I @ | \  sh -c "bcftools view -r $shard1_csv_intervals -o - $S3BUCKET_INPUTS/@ | \  sentieon util vcfconvert - @"

为了获得额外的效率，可以使用“瀑布式”方法，如图 4 和图 5 所示。在这种方法中，计算节点将按顺序处理多个分片，同时运行下一个分片的部分 GVCF 输入文件的下载，并与当前分片的 GVCFtyper 处理并行进行，从而更有效地共享资源，因为这是一个 CPU 密集型和 I/O 密集型的过程。流程如下：

· 下载 shard1 的部分 GVCF 输入文件。

· 启动 shard1 的 GVCFtyper，并与此同时并行下载 shard2 的部分 GVCF 输入文件。

· 在上一步完成之后，启动 shard2 的 GVCFtyper，并与此同时使用 bcftools 下载 shard3。

图 3 在云环境中将 GVCFtyper 分发到多个服务器

图 5 在云环境中将 GVCFtyper 分发到多个服务器，下载下一个分片的输入并与上一个分片的处理并行进行的瀑布式方法的详细说明

（六）Shell 示例

以下是用于分发命令的 Shell 示例，使用 1G 碱基的分片大小，仅选用于演示目的。推荐的分片大小是 100M 碱基。

# Sample file for distributing DNAseq pipeline onto 4 1GBase shards# Each stage command can be distributed to a different server for faster processing,# but the user needs to make sure that the necessary files are present in each machineFASTA="/home/b37/human_g1k_v37_decoy.fasta"FASTQ_1="WGS-30x_1.fastq.gz"FASTQ_2="WGS-30x_2.fastq.gz"FASTQ_INDEX="WGS-30x.fastq.gz.index"KNOWN1="/home/b37/1000G_phase1.indels.b37.vcf.gz"KNOWN2="/home/b37/Mills_and_1000G_gold_standard.indels.b37.vcf.gz"DBSNP="/home/b37/dbsnp_138.b37.vcf.gz"######################################## BWA mapping, distributed on 4 servers#######################################BWA_K_size=100000000num_srvr=4#get the number of runs that the inputs will be split into given the sizenum_K=$(sentieon fqidx query -i $FASTQ_INDEX | cut -d' ' -f1)BWA_K_size=$(sentieon fqidx query -i $FASTQ_INDEX | cut -d' ' -f2)num_K_srvr=$(expr $num_K / $num_srvr + 1)#run multiple BWA on multiple serversfile_list=""for run in $(seq 0 $((num_srvr-1))); do        region="$((run*num_K_srvr))-$((run*num_K_srvr+num_K_srvr))"        #send each of these to a different server        sentieon bwa mem -R "@RG\tID:$GROUP\tSM:$SAMPLE\tPL:$PLATFORM" \          -K $BWA_K_size -t $NT -p $FASTA \          "< sentieon fqidx extract -i $FASTQ_INDEX -r $region $FASTQ_1 $FASTQ_2" | \          sentieon util sort -r $FASTA -t $NT --sam2bam -o sorted_run$run.bam -i -        file_list="$file_list sorted_run$run.bam"done#merge the resultssentieon util merge -o sorted.bam $file_list######################################## define 4 shards#######################################SHARD_0="--shard 1:1-249250621 --shard 2:1-243199373 --shard 3:1-198022430 \  --shard 4:1-191154276 --shard 5:1-118373300"SHARD_1="--shard 5:118373301-180915260 --shard 6:1-171115067 --shard 7:1-159138663 \  --shard 8:1-146364022 --shard 9:1-141213431 --shard 10:1-135534747 \  --shard 11:1-135006516 --shard 12:1-49085594"SHARD_2="--shard 12:49085595-133851895 --shard 13:1-115169878 --shard 14:1-107349540 \   --shard 15:1-102531392 --shard 16:1-90354753 --shard 17:1-81195210 \   --shard 18:1-78077248 --shard 19:1-59128983 --shard 20:1-63025520 \   --shard 21:1-48129895 --shard 22:1-51304566 --shard X:1-118966714"SHARD_3="--shard X:118966715-155270560 --shard Y:1-59373566 --shard MT:1-16569 \   --shard GL000207.1:1-4262 --shard GL000226.1:1-15008 --shard GL000229.1:1-19913 \   --shard GL000231.1:1-27386 --shard GL000210.1:1-27682 --shard GL000239.1:1-33824 \   --shard GL000235.1:1-34474 --shard GL000201.1:1-36148 --shard GL000247.1:1-36422 \   --shard GL000245.1:1-36651 --shard GL000197.1:1-37175 --shard GL000203.1:1-37498 \   --shard GL000246.1:1-38154 --shard GL000249.1:1-38502 --shard GL000196.1:1-38914 \   --shard GL000248.1:1-39786 --shard GL000244.1:1-39929 --shard GL000238.1:1-39939 \   --shard GL000202.1:1-40103 --shard GL000234.1:1-40531 --shard GL000232.1:1-40652 \   --shard GL000206.1:1-41001 --shard GL000240.1:1-41933 --shard GL000236.1:1-41934 \   --shard GL000241.1:1-42152 --shard GL000243.1:1-43341 --shard GL000242.1:1-43523 \   --shard GL000230.1:1-43691 --shard GL000237.1:1-45867 --shard GL000233.1:1-45941 \   --shard GL000204.1:1-81310 --shard GL000198.1:1-90085 --shard GL000208.1:1-92689 \   --shard GL000191.1:1-106433 --shard GL000227.1:1-128374 \   --shard GL000228.1:1-129120 --shard GL000214.1:1-137718 \   --shard GL000221.1:1-155397 --shard GL000209.1:1-159169 \   --shard GL000218.1:1-161147 --shard GL000220.1:1-161802 \   --shard GL000213.1:1-164239 --shard GL000211.1:1-166566 \   --shard GL000199.1:1-169874 --shard GL000217.1:1-172149 \   --shard GL000216.1:1-172294 --shard GL000215.1:1-172545 \   --shard GL000205.1:1-174588 --shard GL000219.1:1-179198 \   --shard GL000224.1:1-179693 --shard GL000223.1:1-180455 \   --shard GL000195.1:1-182896 --shard GL000212.1:1-186858 \   --shard GL000222.1:1-186861 --shard GL000200.1:1-187035 \   --shard GL000193.1:1-189789 --shard GL000194.1:1-191469 \   --shard GL000225.1:1-211173 --shard GL000192.1:1-547496 \   --shard NC_007605:1-171823 --shard hs37d5:1-35477943"######################################## Locus collector########################################Locus Collector for shard 000$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_0 \   --algo LocusCollector --fun score_info .part_deduped.bam000.score.gz \   2> collect000.log#Locus Collector for shard 001$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_1 \   --algo LocusCollector --fun score_info .part_deduped.bam001.score.gz \   2> collect001.log#Locus Collector for shard 002$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_2 \   --algo LocusCollector --fun score_info .part_deduped.bam002.score.gz \   2> collect002.log#Locus Collector for shard 003$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_3 \   --algo LocusCollector --fun score_info .part_deduped.bam003.score.gz \   2> collect003.log#Locus Collector for shard with unmapped reads (NO_COOR)$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam --shard NO_COOR \   --algo LocusCollector --fun score_info .part_deduped.bam004.score.gz \   2> collect004.log######################################## Dedup using all score.gz files########################################Dedup for shard 000$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_0 \   --algo Dedup --score_info .part_deduped.bam000.score.gz \   --score_info .part_deduped.bam001.score.gz \   --score_info .part_deduped.bam002.score.gz \   --score_info .part_deduped.bam003.score.gz \   --score_info .part_deduped.bam004.score.gz \   --rmdup .part_deduped000.bam 2> dedup000.log#Dedup for shard 001$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_1 \   --algo Dedup --score_info .part_deduped.bam000.score.gz \   --score_info .part_deduped.bam001.score.gz \   --score_info .part_deduped.bam002.score.gz \   --score_info .part_deduped.bam003.score.gz \   --score_info .part_deduped.bam004.score.gz \   --rmdup .part_deduped001.bam 2> dedup001.log#Dedup for shard 002$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_2 \   --algo Dedup --score_info .part_deduped.bam000.score.gz \   --score_info .part_deduped.bam001.score.gz \   --score_info .part_deduped.bam002.score.gz \   --score_info .part_deduped.bam003.score.gz \   --score_info .part_deduped.bam004.score.gz \   --rmdup .part_deduped002.bam 2> dedup002.log#Dedup for shard 003$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam $SHARD_3 \   --algo Dedup --score_info .part_deduped.bam000.score.gz \   --score_info .part_deduped.bam001.score.gz \   --score_info .part_deduped.bam002.score.gz \   --score_info .part_deduped.bam003.score.gz \   --score_info .part_deduped.bam004.score.gz \   --rmdup .part_deduped003.bam 2> dedup003.log#Dedup for shard with unmapped reads (NO_COOR)$SENTIEON_FOLDER/bin/sentieon driver -t 16 -i sorted.bam --shard NO_COOR \   --algo Dedup --score_info .part_deduped.bam000.score.gz \   --score_info .part_deduped.bam001.score.gz \   --score_info .part_deduped.bam002.score.gz \   --score_info .part_deduped.bam003.score.gz \   --score_info .part_deduped.bam004.score.gz \   --rmdup .part_deduped004.bam 2> dedup004.log#Merge bam files from all shards into final output$SENTIEON_FOLDER/bin/sentieon util merge -i .part_deduped000.bam \   -i .part_deduped001.bam -i .part_deduped002.bam \   -i .part_deduped003.bam -i .part_deduped004.bam \   -o deduped.bam --mergemode=10 2> dedup_merge.log######################################## BQSR########################################QualCal for shard 000$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped000.bam \   -i .part_deduped001.bam $SHARD_0 --algo QualCal -k $KNOWN1 -k $KNOWN2 \   .part_recal_data000.table 2> bqsr000.log#QualCal for shard 001$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped000.bam \   -i .part_deduped001.bam -i .part_deduped002.bam $SHARD_1 --algo QualCal \   -k $KNOWN1 -k $KNOWN2 .part_recal_data001.table 2> bqsr001.log#QualCal for shard 002$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped001.bam \   -i .part_deduped002.bam -i .part_deduped003.bam $SHARD_2 --algo QualCal \   -k $KNOWN1 -k $KNOWN2 .part_recal_data002.table 2> bqsr002.log#QualCal for shard 003$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped002.bam \   -i .part_deduped003.bam $SHARD_3 --algo QualCal -k $KNOWN1 -k $KNOWN2 \   .part_recal_data003.table 2> bqsr003.log#QualCal for shard with unmapped reads (NO_COOR)$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped004.bam \   --shard NO_COOR --algo QualCal -k $KNOWN1 -k $KNOWN2 \   .part_recal_data004.table 2> bqsr004.log#Merge table files into complete calibration table$SENTIEON_FOLDER/bin/sentieon driver --passthru --algo QualCal \   --merge recal_data.table .part_recal_data000.table .part_recal_data001.table \   .part_recal_data002.table .part_recal_data003.table .part_recal_data004.table \   2> bqsr_merge.log######################################## Variant Calling########################################Haplotyper for shard 000$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped000.bam \   -i .part_deduped001.bam -q recal_data.table $SHARD_0 --algo Haplotyper \   -d $DBSNP .part_output000.vcf.gz 2> hc000.log#Haplotyper for shard 001$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped000.bam \   -i .part_deduped001.bam -i .part_deduped002.bam  -q recal_data.table \   $SHARD_1 --algo Haplotyper -d $DBSNP .part_output001.vcf.gz 2> hc001.log#Haplotyper for shard 002$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped001.bam \   -i .part_deduped002.bam -i .part_deduped003.bam -q recal_data.table \   $SHARD_2 --algo Haplotyper -d $DBSNP .part_output002.vcf.gz 2> hc002.log#Haplotyper for shard 003$SENTIEON_FOLDER/bin/sentieon driver -t 16 -r $FASTA -i .part_deduped002.bam \   -i .part_deduped003.bam -q recal_data.table $SHARD_3 --algo Haplotyper \   -d $DBSNP .part_output003.vcf.gz 2> hc003.log#There is no need to do variant calling on the NO_COOR shard, as there will# be no variants on unmapped reads#Merge output files into output VCF$SENTIEON_FOLDER/bin/sentieon driver --passthru --algo Haplotyper \   --merge output.vcf.gz .part_output000.vcf.gz .part_output001.vcf.gz \   .part_output002.vcf.gz .part_output003.vcf.gz 2> hc_merge.log

以下是一个用于创建分片的Shell示例。generate_shards.sh

determine_shards_from_bam(){    local bam step tag chr len pos end    bam="$1"    step="$2"    pos=1    samtools view -H $bam |\    while read tag chr len; do        [ $tag == '@SQ' ] || continue        chr=$(expr "$chr" : 'SN:\(.*\)')        len=$(expr "$len" : 'LN:\(.*\)')        while [ $pos -le $len ]; do            end=$(($pos + $step - 1))            if [ $pos -lt 0 ]; then                start=1            else                start=$pos            fi            if [ $end -gt $len ]; then                echo -n "$chr:$start-$len "                pos=$(($pos-$len))                break            else                echo "$chr:$start-$end"                pos=$(($end + 1))            fi        done    done    echo "NO_COOR"}determine_shards_from_dict(){    local bam step tag chr len pos end    dict="$1"    step="$2"    pos=1    cat $dict |\    while read tag chr len UR; do        [ $tag == '@SQ' ] || continue        chr=$(expr "$chr" : 'SN:\(.*\)')        len=$(expr "$len" : 'LN:\(.*\)')        while [ $pos -le $len ]; do            end=$(($pos + $step - 1))            if [ $pos -lt 0 ]; then                start=1            else                start=$pos            fi            if [ $end -gt $len ];then                echo -n "$chr:$start-$len "                pos=$(($pos-$len))                break            else                echo "$chr:$start-$end"                pos=$(($end + 1))            fi        done    done    echo "NO_COOR"}determine_shards_from_fai(){    local bam step tag chr len pos end    fai="$1"    step="$2"    pos=1    cat $fai |\    while read chr len other; do        while [ $pos -le $len ]; do            end=$(($pos + $step - 1))            if [ $pos -lt 0 ]; then                start=1            else                start=$pos            fi            if [ $end -gt $len ]; then                echo -n "$chr:$start-$len "                pos=$(($pos-$len))                break            else                echo "$chr:$start-$end"                pos=$(($end + 1))            fi        done    done    echo "NO_COOR"}if [ $# -eq 2 ]; then    filename=$(basename "$1")    extension="${filename##*.}"    if [ "$extension" = "fai" ]; then        determine_shards_from_fai $1 $2    elif [ "$extension" = "bam" ]; then        determine_shards_from_bam $1 $2    elif [ "$extension" = "dict" ]; then        determine_shards_from_dict $1 $2    fielse    echo "usage $0 file shard_size"fi