问题标签 [vcf-variant-call-format]

我尝试使用一种工具，我需要在输入中使用通配符。

这是一个例子：

我有这个错误：

所以我有一个带有案例和控制的功能。我无法提取代码。

为了加快某个蛇形步骤，我想：

• 使用此结果将我的 bamfile 每个染色体拆分
  bamtools split -in sample.bam --reference
  为名为的文件sample.REF_{chromosome}.bam
• 对每个结果执行变体调用，例如sample.REF_{chromosome}.vcf
• 使用 vcf-concat (VCFtools) 重新组合获得的 vcf 文件
  vcf-concat file1.vcf file2.vcf file3.vcf > sample.vcf

问题是我不知道哪些染色体可能在我的 bam 文件中。所以我无法准确指定bamtools split. 此外，我不确定如何输入vcf-concat以获取所有 vcf 文件。

我想过使用 samples.fofn 并做类似的事情

并使用相同fofn的方法连接获得的 vcf 文件。但这感觉就像一个非常尴尬的黑客，我很感激你的建议。

编辑 20180409

正如@jeeyem 所建议的，我尝试了这些dynamic()功能，但我无法弄清楚。

我完整的蛇文件在GitHub 上，动态部分在第 99-133 行。

我得到的错误是：（ InputFunctionException in line 44 of /home/wdecoster/DR34/SV-nanopore.smk: KeyError: 'anon___snakemake_dynamic' Wildcards: sample=anon___snakemake_dynamic 使用anon匿名 {sample} 标识符）

使用 --debug-dag 运行会给出（出错前的最后一部分）： candidate job cat_vcfs wildcards: sample=anon candidate job nanosv wildcards: sample=anon___snakemake_dynamic, chromosome=_ candidate job samtools_index wildcards: aligner=split_ngmlr, sample=anon___snakemake_dynamic.REF__ candidate job split_bam wildcards: sample=anon___snakemake_dynamic, chromosome=_ InputFunctionException in line 44 of /home/wdecoster/DR34/SV-nanopore.smk: KeyError: 'anon___snakemake_dynamic' Wildcards: sample=anon___snakemake_dynamic

这表明通配符被误解了？

干杯，沃特

我bcftools consensus用来从 vcf 文件中提取单倍型。给定输入文件：

创建以下输出文件：

我目前有两条规则可以做到这一点。它们的区别仅在于输出文件和 shell 命令中的数字 1 和 2。代码：

虽然这段代码有效，但似乎应该有一种更好、更 Pythonic 的方式来只使用一个规则来做到这一点。是否可以将其合并为一个规则，或者我目前的方法是最好的方法？

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

目前我正在使用signeR允许您创建体细胞突变特征的包。首先，我想从vignette重现他们的结果。

如您所见，输入数据文件有两个选项，.vcf并且先前已预处理 mutation counts file，这是通过它们的genCountMatrixFromVcf()功能获得的简单数据框..实际上来自.vcf...

当然，您可以在这里找到一些示例.vcf文件和这个突变计数矩阵，以供进一步分析。

在小插图中，他们使用（如我所读）这些数据：SUBSTITUTIONS_13Apr2012_snz.txt。所以我决定.vcf从这个数据框创建一个文件，并完成小插图的所有步骤。

所以我做了几个简单的步骤来构建.vcf上述数据：

结果我们有：

文件的所有必填列.vcf都存在。我以前制作过这种类型的文件，比我使用的readVcfAsVRanges()一切都好。

然而，在signeR包装中存在一些问题......因为如果我不创建这个突变计数矩阵，我会遇到一个错误：

所以我试图找出我的.vcf文件与“工作”文件相比有什么问题......

我的文件：

他们的“工作”示例文件（实际上在运行genCountMatrixFromVcf()功能之后也有错误但不同......此时并不重要）。

正如你所看到的dim那样，我在0那里的价值是不同的，所以这可能是问题所在。所以我想问你如何修复我的.vcf文件？

我想从 vcf 文件中获取每个种群中的个体数量，按照种群的读取顺序。我文件的字段如下所示

在此处查看示例文件vcf 文件

例如，在我链接到的文件中，我有两个群体，Chalifour 2003 和 Chalifour 2015。个人有一个前缀“CHALIFOUR_2003 ...”来标识这一点。

我希望能够提取如下内容： Chalifour_2003* 35 Chalifour 2015* 45

用“35”和“45”表示每个种群中的个体数量（尽管这些数字是虚构的）。我根本不关心输出的格式，我只需要数字，重要的是要按照将人口读入文件的顺序列出人口。

任何有关尝试获取此信息的途径的建议将不胜感激。

我从之前在较小数据集上多次使用过的代码中遇到了一种从未见过的奇怪行为。我正在使用 Pandas 数据帧 read_table 解析 VCF 文件。VCF 文件有一个标题，然后是 9 列 + 多列个人。在我以前for row in genomes_df.itertuples():遍历数据帧的每一行之前，我可以用row.SVLEN. 当我检查type(row)它是一个熊猫对象时。今天我在相同 VCF 格式的更大文件（350 列 vs 之前 10 列）上运行我的脚本，它给了我AttributeError: 'tuple' object has no attribute 'SVLEN'，因为现在type(row)是一个元组！

这里发生了什么？列名不同（NWD107911.mark_dupesvs NWD107911），但我检查了名称中没有空格（在另一篇文章中阅读它可能导致不同的行为）。

我有一个基于 apache-beam 的数据流作业，可以使用vcf 源从单个文本文件（存储在谷歌云存储中）读取，将文本行转换为数据存储Entities并将它们写入数据存储接收器。工作流程运行良好，但我注意到的缺点是：

• 写入数据存储的速度最多约为每秒 25-30 个实体。
• 我尝试使用--autoscalingAlgorithm=THROUGHPUT_BASED --numWorkers=10 --maxNumWorkers=100，但执行似乎更喜欢一名工作人员（见下图：目标工作人员曾经增加到 2，但“基于在当前运行步骤中并行化工作的能力”减少到 1）。

我没有使用祖先路径作为键；所有的实体都是一样的kind。

管道代码如下所示：

因为我有数百万行要写入数据存储区，所以以 30 个实体/秒的速度写入会花费很长时间。

问题：输入只是一个巨大的 gzip 文件。我是否需要将其拆分为多个小文件以触发多个工作人员？有没有其他方法可以加快导入速度？我错过了num_workers设置中的某些内容吗？谢谢！

所以我有一些 Python 代码（运行 Python 2.7.12），它使用 VEP 针对特定的转录本注释 vcf 文件。

当我通过双击它（或从命令提示符运行它）运行脚本时，它会给出以下错误：

C:\annotatetsca_KH\annotate.py:364: UnicodeWarning: Unicode 相等比较未能将两个参数转换为 Unicode - 将它们解释为不相等

如果指定的成绩单中的成绩单：

然后，这会导致某些变体未针对指定的转录本进行注释。但是，如果我将脚本加载到 IDLE 中，然后通过 IDLE 运行完全相同的代码，我不会收到此错误，并且变量已正确注释。

有谁知道为什么通过 IDLE 运行脚本会产生不同的结果？据我所知，以及与我交谈过的一些同事，IDLE 和 Python 应该以完全相同的方式对待一切。有没有人知道一种方法来绕过必须通过 IDLE 运行脚本以避免此错误，因为该脚本被许多不熟悉 Python/IDLE 的人使用，并且会发现只运行脚本要容易得多通过双击它（像任何其他 Windows 应用程序一样）。

编辑（已解决）：

我找到了问题的原因。成绩单列表中的几份成绩单末尾有一个空格。似乎 IDLE 会去除尾随空格，但直接运行 python 脚本不会。因此，当直接运行它时，尾随空格被转换为“\xa0”，因此当它试图将脚本与指定列表中的一个匹配时，它不匹配。我已经删除了所有空格，它工作正常，我还将添加一行以在脚本列表中加载时删除所有空格！

我正在尝试在不存在的行中添加一个“chr”字符串。此操作仅在没有'##' 的行中是必需的。起初我使用 grep + sed 命令，如下所示，但我想运行覆盖原始文件的命令。

因此，要在文件中运行命令，我会这样写：

这是 vcf 文件的内容。