如何将所有染色体组合在一个文件中答案

【问题标题】：How to combine all chromosomes in a single file如何将所有染色体组合在一个文件中
【发布时间】：2018-09-13 18:23:16
【问题描述】：

我下载了 VCF 格式的 1000 个基因组数据（染色体 1 -22）。如何将所有染色体组合在一个文件中？我是否应该先将所有染色体转换为 plink 二进制文件，然后再执行--bmerge mmerge-list？或者有没有其他方法可以将它们结合起来？请问有什么建议吗？

【问题讨论】：

有什么好的理由把它们结合起来吗？

标签： bioinformatics genetics vcf-variant-call-format

【解决方案1】：

您可以按照您的建议使用 PLINK。您也可以使用 BCFtools：

https://samtools.github.io/bcftools/bcftools.html

具体来说，concat 命令：

bcftools concat ALL.chr{1..22}.phase3_shapeit2_mvncall_integrated_v5a.20130502.genotypes.vcf.gz -Oz -o  ALL.autosomes.phase3_shapeit2_mvncall_integrated_v5a.20130502.genotypes.vcf.gz

如果您使用 PLINK，您可能会遇到 1000 Genomes 的问题，因为它包含与 PLINK 不兼容的多等位基因 SNP。另外，有些SNP具有相同的RS标识符，这也与PLINK不兼容。

您需要通过将多等位基因 SNP 拆分为多个记录并删除具有重复 RS 标识符的记录（或创建一个新的唯一标识符）来解决这些问题，以使 PLINK 正常工作。

此外，PLINK 二进制 PED 不支持基因型概率。我不记得 1000 Genomes 是否包含此类信息。如果确实如此并且您想保留它，则无法将其转换为二进制 PED，因为基因型概率将被硬调用，请参阅：

https://www.cog-genomics.org/plink2/input

具体来说：

由于 PLINK 1 二进制格式不能表示基因型概率，不确定性大于 0.1 的调用通常是被视为失踪，其余的被视为硬呼叫。

因此，如果您打算为输出保留 VCF 格式，我建议您不要使用 PLINK。

编辑

以下是将VCF转换为PLINK的方法：

要从 VCF 文件构建 PLINK 兼容文件，需要合并或删除重复的位置和 SNP id。在这里，我选择删除所有重复的条目。目录重复 SNP id：

grep -v '^#' <(zcat ALL.chr${chrom}.phase3_shapeit2_mvncall_integrated_v5a.20130502.genotypes.vcf.gz) | cut -f 3 | sort | uniq -d > ${chrom}.dups

使用 BCFTools，拆分多等位基因 SNP，并使用 plink 删除在上一步中找到的重复 SNP id：

bcftools norm -d both -m +any -Ob ALL.chr${chrom}.phase3_shapeit2_mvncall_integrated_v5a.20130502.genotypes.vcf.gz | plink --bcf /dev/stdin --make-bed --out ${chrom} --allow-extra-chr 0 --memory 6000 --exclude ${chrom}.dups

重要的是，这不是解决将 VCF 转换为 PLINK 问题的唯一方法。例如，您可以为重复的 RS id 唯一分配标识符。

【讨论】：

如果您决定使用 PLINK，我很高兴使用 BCFtools 和 PLINK 命令编辑此答案，以解决上述限制。
到目前为止，我尝试了 PLINK，并得到了“多等位基因位置错误”。你是对的文斯。我的目标是用以下条件修剪数据：--maf 0.01 --snps-only --indep-pairwise 50 10 0.2。所需文件应该是单个文件（最好在 PLINK 中）。我只能使用 bcftools 执行此操作吗？
我编辑了答案以包括我用来将 VCF 转换为 PLINK 的方法。在此之后，您可以使用 --merge-list 指定将染色体合并到一个文件中。
我不确定 PLINK 1.9 是否可以选择接受多等位基因位点。如果我记得它只是错误退出。 PLINK 2.0 似乎删除了它们，因此您可以使用它来读取 VCF：cog-genomics.org/plink/2.0/input#vcf。
非常感谢您的详细解答。看起来 PLINK 1.9 可以做大部分事情。我使用 bcftools concant 合并所有 22 个 VCF 文件，并使用 PLINK 1.9 使用所需的 maf 和 LD 进行修剪，因为 PLINK 使用 VCF 文件。好吧，我期待输出文件中有一些多位置变体，但我想知道，我找不到任何一个。在我看来，当我用所需的 LD 和 maf 修剪 VCF 文件时，也修剪了多位置变体。你认为它的正确方法还是我必须先使用上面的两个脚本行？

【解决方案2】：

picard GatherVcfs https://broadinstitute.github.io/picard/command-line-overview.html

将分散操作中的多个 VCF 文件收集到单个 VCF 文件中。输入文件必须按基因组顺序提供，并且不能在重叠位置包含事件。

【讨论】：