2

我是 R 编程新手，并试图完成一项非常具体的任务。

我有一个包含 n 个样本的 fasta 序列，我读到ape：

library(ape)

matrix <- read.dna(myfasta, format="fasta", as.character=TRUE)

这创建了一个矩阵，如下所示：

|    | V1 | V2 | V3 | V4 |...
|------------------------|
|Seq1|  a |  t |  g |  c |...
|Seq2|  a |  t |  g |  a |...
|Seq3|  a |  t |  c |  c |...
|Seq4|  t |  t |  g |  a |...
|... |

其中 Seq(n) 是每个样本的 DNA 序列，V(n) 表示核苷酸位置。

如何选择在某个位置（例如“V1”）带有某个核苷酸（例如“a”）的序列，然后将这些序列作为连接字符串返回？

因此，对于 V1 位置，我想要“Seq1，Seq2，Seq3”之类的东西，而对于 V4 位置，对于相同的基础，我想要“Seq2，Seq4”

我试过了which()，filter(matrix, V1 == "a")但我很挣扎。

提前致谢！

4

1 回答 1

2

最简单的方法是选择V1 == 'a'具有布尔索引的行，然后提取rownames：

rownames(example[example[,"V1"] == "a", ]) # "No304" "No306"

你提到filter的，看起来像dplyr。使用 tidyverse 方法来操作对行名很重要的数据有点麻烦，因为默认情况下会删除行名。

如果要使用filter，则必须首先将行名保存为显式列：

library(dplyr)

as.data.frame(example) %>% 
  mutate(sequence = rownames(.), .before = everything()) %>% 
  filter(V1 == "a") %>% 
  select(sequence)

  sequence
1    No304
2    No306

数据（来自aperead.dna 文档）

library(ape)

cat(">No305",
    "NTTCGAAAAACACACCCACTACTAAAANTTATCAGTCACT",
    ">No304",
    "ATTCGAAAAACACACCCACTACTAAAAATTATCAACCACT",
    ">No306",
    "ATTCGAAAAACACACCCACTACTAAAAATTATCAATCACT",
    file = "exdna.fas", sep = "\n")

example <- read.dna("exdna.fas", format = "fasta", as.character = TRUE)
colnames(example) <- paste0("V", 1:ncol(example))

example
      V1  V2  V3  V4 ...
No305 "n" "t" "t" "c"
No304 "a" "t" "t" "c"
No306 "a" "t" "t" "c"
于 2020-12-30T16:35:54.980 回答