从 R 中的局部比对识别氨基酸替换

Question

我想确定我的序列中特定感兴趣区域的位置和氨基酸变化，并将该信息存储在 table。

是否可以通过使用 bioconductor 包在 R 中做这样的事情？我设法用 DECIPHER 包中的 AlignSeqs() 对序列进行了简单的比对，但我无法自动提取序列中的差异。我从 FASTA 文件开始。

我想以这样的方式结束：

Isolate ID    Reference_AA   Sample_AA   Pos
1             S              T           254
2             T              D           200
3             L              A           230

我有一个 74 AA 长的参考序列，我想查看查询序列（比参考长得多）与参考序列的差异，并将结果列在 table.位置列与参考序列中的位置有关，与查询序列中的位置无关。我希望参考序列中的第一个 AA 从 68 开始，而不是 1。

我发现很难为此添加示例序列，因为它们往往很长，但这里有一些更短的东西可以使用（与上面的 table 无关）：

>ref
VGRALPDVR

>query1
KSSYLDYAMSVIVGTALPDVRDGLKPVHRRVLY

>query2
ELKSSYLDYAMSVIVGRAAPDVRDGLKPV

预期输出：

ID       Reference_AA    Sample_AA   Pos
query1   R               T           70
query2   A               L           72

Answer 1

既然你想要在参考中的位置，你可以只对参考序列使用一系列成对比对。 biostrings 包含一个 mismatchTable 函数，可以为您提供包含所有所需信息的数据框。使用 dplyr:

进行一些重新格式化

library(Biostrings)
library(dplyr)

seqs<-readAAStringSet("test.fa")

mismatches <- function(query, ref) {
    pairwiseAlignment(ref, query, substitutionMatrix = "BLOSUM50",
                      gapOpening = 3, gapExtension = 1) %>%
      mismatchTable() %>%
      mutate(ID=names(query), 
             Pos=PatternStart+67, 
             Reference_AA=as.character(PatternSubstring),
             Sample_AA=as.character(SubjectSubstring)) %>% 
             select(ID, Reference_AA, Sample_AA, Pos)
}  

bind_rows(mismatches(seqs[2], seqs[1]), mismatches(seqs[2], seqs[1]))

#>      ID Reference_AA Sample_AA Pos
#>1 query1            R         T  70
#>2 query2            L         A  72

编辑

以下是使用 lapply 循环输入的方法：

bind_rows(lapply(seq_along(seqs[-1]), function(i) mismatches(seqs[i+1], seqs[1])))
#>   ID Reference_AA Sample_AA Pos
#>1 ref            R         T  70
#>2 ref            L         A  72

从 R 中的局部比对识别氨基酸替换

Identifying amino acid substitutions from local alignments in R

r

bioinformatics

bioconductor