问题标签 [samtools]

我正在尝试从 qsub 调用运行一些 samtools 命令（在集群上运行）。由于某种原因，这些命令似乎无法被识别。但是，如果我复制粘贴命令并直接从终端集群运行它，它就可以正常工作。有没有人遇到过这样的问题或知道我做错了什么？谢谢，

帕特里克

我的 qsub（这不起作用）：

如果我直接从终端运行命令，它可以工作：

samtools index /data/test.bam /data/test.bai

假设我从 shell 运行了以下命令

对于那些不熟悉 samtools 视图的人（因为这是 stackoverflow）。这本质上是在创建一个具有新标头的新 bam 文件。bam 文件通常是大型压缩文件，因此在某些情况下，即使通过文件也可能很耗时。另一种方法是执行 command2，然后使用 samtools reheader 切换标题。这会通过大文件两次。上述命令通过 bam 一次，这对于较大的 bam 文件很有用（即使压缩后它们也会大于 20GB - WGS）。

我的问题是如何使用子进程在 python 中实现这种类型的命令。

我有以下内容：

任何帮助是极大的赞赏。谢谢。

努力定义导致此错误的原因。

安装了anaconda并使用了conda install pysam，工作了一段时间，但突然出现此错误

看了一下文件里面对应的代码是

知道是什么原因造成的吗？？？

我正在尝试解析来自 samtools mpileup 的文本输出。我从一个字符串开始

每当我+后面有一个整数n时，我想选择该整数后面的 n 个字符并将整个内容替换为*. 所以对于这个测试用例，我会有

我有正则表达式 \+[0-9]+[ACGTNacgtn]+，但这会导致输出.$......*.*.*和尾随的 G 也丢失。如何选择 n 个字符，其中 n 事先不知道但在字符串本身中指定？

我正在尝试创建一个脚本，从 sam 文件的标题中删除读取组。从命令行运行的代码如下。

这是输入文件的示例

运行上述命令后，输出为

基本上，我只是删除以“@RG”开头且不是 ID e2_20 的行。

问题是，如果我在 bash 脚本中运行此命令，则该命令不起作用。

脚本如下。

该脚本在 SGE 集群中运行，这就是脚本顶部出现奇怪语法的原因。脚本的名称是 reHeadBams.bash。我通过输入以下命令在 shell 上运行脚本：

该命令的参数是输入文件，然后是输出文件，最后是我要查找的样本或读取组。

脚本的输出如下所示：

因此脚本删除了所有读取组，而不是 ID 为 e2_9 的组。

我从脚本中回显了命令，输出是

虽然我可能弄错了，但问题似乎是 awk 使用我的命令行参数而不是输入文件的列作为 $1 和 $2。有谁知道为什么会这样？

很抱歉对我的问题进行了冗长的描述。如果您需要任何澄清，请告诉我。

我很快意识到生物信息学并不是一个术语定义明确且易于访问的学科。我与我的一些结果有明显的差异。

我samtools view -b -h -f 8 fileName.bam > mateUnmapped.bam在几个 BAM 文件上使用过。我的印象是该命令仅提取其伙伴与草稿基因组不对齐的读取（还包括标题；输出为 BAM 格式）

当我samtools 'flagstat'在结果文件上使用时，我得到了一个有趣的结果：“单例”的数量与读取的总数不匹配……这对我来说似乎很奇怪。

我能找到的唯一和解是在这里：

http://seqanswers.com/forums/showthread.php?t=46711

一位回答本论坛提出的问题的人声称，单例有时被定义为根本没有伙伴读取的序列。但是，这仍然不能解释我的结果。Flagstat 说我大约 40% 的读取是单例，但我觉得根据我使用的“查看”命令，它们应该都是单例。

经验丰富的生物信息学家可以帮助我吗？

我有一个间歇性地从 htsjdk.samtools.reference.IndexedFastaSequenceFile 获取异常的问题：

htsjdk.samtools.SAMException：序列字典和索引包含不同数量的重叠群或 htsjdk.samtools.SAMException：无法获取...

问题是我一遍又一遍地调用相同的代码，只是间歇性地得到这个错误。

我认为错误消息实际上具有误导性，尤其是在序列字典响应的情况下。

我正在尝试从这里访问 htsjdk.jar 提供的方法： https ://samtools.github.io/htsjdk/

并记录在这里： https ://samtools.github.io/htsjdk/javadoc/htsjdk/index.html

使用 jython。我需要访问/查询 BAM 文件索引（BAI 文件）以获取二进制 BAM 文件中的起始位置的方法。测试 BAM 和 BAI 文件可以从以下网址获得： https ://github.com/samtools/htsjdk/tree/master/testdata/htsjdk/samtools/BAMFileIndexTest

在放入 Jython 注册表后的 jython 2.7.0 中：

我可以访问一些方法，例如 bai_foo.getNumberOfReferences() 等，但所需的方法
getBinsOverlapping(int referenceIndex, int startPos, int endPos) 在 BrowseableBAMIndex 接口中。

但是当谈到在 Jython 中对 Java 类进行子类化时，我迷失了方向。任务是获取与给定基因组位置相对应的 BAM 文件块列表。对于测试 BAM/BAI 文件，即对于 chrM 10000-15000（染色体、start_pos、end_pos），我使用现成的 samtools 独立程序而不是 htsjdk 获得 11 个映射读取：

非常感谢您的帮助

达雷克

编辑：重新成为 groovy 部分 Groovy 版本：2.4.4

上面的代码在 groovy 中工作。我的问题不是这段代码不起作用，而是我不知道从处理 BAI 文件格式的 Java 类中访问方法的正确方法。我确实在 htsjdk/src/java/htsjdk/samtools/*java 文件（来自 repo@github 的 git clone）中搜索了 AbstractBAMFileIndex，但仍然不清楚我需要做什么。

使用 Sam-tools Java API 如何获得 seq。在开始和结束之间，我正在编写以下方法，它会抛出 Exception 。此方法写入以从文件中获取序列。

方法在这里

这是发生堆栈跟踪异常

我正在尝试在 openSUSE 上安装 samtools，我这样做了：

工作得很好。

工作得很好。

对于 samtools：

产生这个输出：

我不明白为什么这行不通。我知道如何使工具工作，但我是 linuxOS 的新手。