字符串使用通配符匹配,试图Biostrings包

问题描述:

给出字符串patt字符串使用通配符匹配,试图Biostrings包

patt = "AGCTTCATGAAGCTGAGTNGGACGCGATGATGCG" 

我们可以更短的子str_col的集合:

str_col = substring(patt,1:(nchar(patt)-9),10:nchar(patt)) 

我们要匹配一个subject1

subject1 = "AGCTTCATGAAGCTGAGTGGGACGCGATGATGCGACTAGGGACCTTAGCAGC" 

在中处理“N”作为通配符(匹配subject1中的任何字母),因此str_col中的所有子字符串都与subject1匹配。

我想在字符串的大型数据库中做这种字符串匹配,我发现BiostringsBiostrings是非常有效的做到这一点。但是,为了高效,Biostrings要求您使用函数PDidc()将您的子串集合(此处为str_col)转换为类pdict的字典。您可以稍后在countPDict()等函数中使用此'词典'来计算与目标的匹配。

为了使用通配符,您必须将字典分为三部分:头(左),可信带(中)和尾(右)。您只能在头部或尾部使用通配符,例如“N”,但不能使用可信任的波段,并且不能使用宽度= 0的可信任波段。因此,例如,如果您使用的是str_col[15],则不会匹配最小宽度= 1的可信任频段:

> PDict(str_col[1:15],tb.start=5,tb.end=5) 
Error in .Call2("ACtree2_build", tb, pp_exclude, base_codes, nodebuf_ptr, : 
    non base DNA letter found in Trusted Band for pattern 15 

因为“N”在可信频段中是正确的。注意这里的字符串是DNA序列,所以“N”是“与A,C,G或T匹配”的代码。

> PDict(str_col[1:14],tb.start=5,tb.end=5) #is OK 
TB_PDict object of length 14 and width 10 (preprocessing algo="ACtree2"): 
    - with a head of width 4 
    - with a Trusted Band of width 1 
    - with a tail of width 5 

有没有办法规避Biostrings这个限制?我也试图用R基函数来执行这样的任务,但我无法想出任何东西。

我估计你需要匹配IUPAC ambiguity code的某些其他通配符,不是吗?

如果您需要完美匹配并且基本功能对您来说足够了,则可以使用与函数glob2rx()相同的技巧:只需使用gsub()进行转换即可构造匹配模式。举个例子:

IUPACtoRX <- function(x){ 
    p <- gsub("N","\\[ATCG\\]",x) 
    p <- gsub("Y","\\[CT\\]",p) #match any pyrimidine 
    # add the ambiguity codes you want here 
    p 
} 

显然,你需要为你要设定在每个歧义一条线,但它是非常简单的我会说。

这样做,你就可以如这样做:

> sapply(str_col, function(i) grepl(IUPACtoRX(i),subject1)) 
AGCTTCATGA GCTTCATGAA CTTCATGAAG TTCATGAAGC TCATGAAGCT CATGAAGCTG ATGAAGCTGA 
     TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE 
TGAAGCTGAG GAAGCTGAGT AAGCTGAGTN AGCTGAGTNG GCTGAGTNGG CTGAGTNGGA TGAGTNGGAC 
     TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE 
GAGTNGGACG AGTNGGACGC GTNGGACGCG TNGGACGCGA NGGACGCGAT GGACGCGATG GACGCGATGA 
     TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE 
ACGCGATGAT CGCGATGATG GCGATGATGC CGATGATGCG 
     TRUE  TRUE  TRUE  TRUE 

要找到匹配的数量,你可以使用如gregexpr()

> sapply(str_col,function(i) sum(gregexpr(IUPACtoRX(i),subject1) > 0)) 
AGCTTCATGA GCTTCATGAA CTTCATGAAG TTCATGAAGC TCATGAAGCT CATGAAGCTG ATGAAGCTGA 
     1   1   1   1   1   1   1 
TGAAGCTGAG GAAGCTGAGT AAGCTGAGTN AGCTGAGTNG GCTGAGTNGG CTGAGTNGGA TGAGTNGGAC 
     1   1   1   1   1   1   1 
GAGTNGGACG AGTNGGACGC GTNGGACGCG TNGGACGCGA NGGACGCGAT GGACGCGATG GACGCGATGA 
     1   1   1   1   1   1   1 
ACGCGATGAT CGCGATGATG GCGATGATGC CGATGATGCG 
     1   1   1   1 
+0

不错!非常感谢! – vitor 2014-09-02 14:14:12

+0

不确定如果您可以使用长度(gregexpr()),因为它也会返回1以表示不匹配。例如。如果你省略函数IUPACtoRX,它也会为每个人返回1。 – vitor 2014-09-02 14:49:16

+0

@agui不错,我忘记了gregexpr在不匹配的情况下返回-1。我已经改编 – 2014-09-02 16:16:48