bioconductor和R，在尝试检索基因序列时使用getBM()函数时出现错误

Question

我正在尝试使用bioconductor的getBM()函数来检索某些基因的序列。下面是我的脚本：

setwd(".")
options(stringsAsFactors = FALSE)
cat("\014")


if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
listOfBiocPackages <- c("annotate", "GEOquery", "biomaRt")

bioCpackagesNotInstalled <- which( !listOfBiocPackages %in% rownames(installed.packages()) )
cat("package missing listOfBiocPackages[", bioCpackagesNotInstalled, "]: ", listOfBiocPackages[bioCpackagesNotInstalled], "\n", sep="")

# check there's still something left to install
if( length(bioCpackagesNotInstalled) ) {
    BiocManager::install(listOfBiocPackages[bioCpackagesNotInstalled])
}

library("easypackages")
libraries(listOfBiocPackages)  

mart <- useMart("ENSEMBL_MART_ENSEMBL")
mart <- useDataset("hsapiens_gene_ensembl", mart)

thisAnnotLookup <- getBM(mart=mart, attributes=c("affy_hg_u133a", "ensembl_gene_id", "gene_biotype", "external_gene_name", "gene_flank"), filter=c("affy_hg_u133a", "upstream_flank"), values=list(affyid=c("203423_at",   "204088_at",   "204511_at",   "204911_s_at", "205234_at")), uniqueRows=TRUE, checkFilters=FALSE)

print(thisAnnotLookup)

此脚本显示以下错误：

  'names' attribute [2] must be the same length as the vector [1]

我不能理解。

有人能帮我解决这个问题吗？谢谢!

Answer 1

如果您有两个过滤器，则需要提供两个元素的列表。然而，更大的问题是，根据您的连接，upstream_flank可能不可用，请参阅this problem cropping up for many users。对我来说，我最后一次使用的是亚洲镜像，它起作用了：

library(biomaRt)

mart <- useEnsembl(biomart = "ensembl", 
                   dataset = "hsapiens_gene_ensembl", 
                   mirror = "asia")

然后：

S = getSequence(id=c("203423_at","204088_at","204511_at","204911_s_at","205234_at"),
type="affy_hg_u133a",seqType="gene_flank",upstream = 20,mart = mart)

thisAnnotLookup <- getBM(mart=mart, 
attributes=c("affy_hg_u133a", "ensembl_gene_id", "gene_biotype", "external_gene_name"), 
filter="affy_hg_u133a", 
values=c("203423_at","204088_at","204511_at","204911_s_at","205234_at"), 
uniqueRows=TRUE, checkFilters=FALSE)

merge(thisAnnotLookup,S)
  affy_hg_u133a ensembl_gene_id   gene_biotype external_gene_name
1     203423_at ENSG00000114115 protein_coding               RBP1
2     204088_at ENSG00000135124 protein_coding              P2RX4
3     204511_at ENSG00000006607 protein_coding              FARP2
4   204911_s_at ENSG00000110171 protein_coding              TRIM3
5     205234_at ENSG00000168679 protein_coding            SLC16A4
            gene_flank
1 CTCCTCTTCCTTTGTAGGGG
2 GATTCCTCTCACCTCGGCCT
3 CGCGGGGGCTGCCGGCGGGC
4 CACCTTTCTACCCCTTAACT
5 GCCTGAGATGACATCAAATC

在我看来太挑剔了。我认为直接从基因组中获取序列会更好。