R的字符串处理很慢?

ryusukekenji

谢老师和yanlinlin82 前辈的程序完全看不懂，太强了！！！

<br />
>seqChar = strsplit(seqStr, "")[[1]]<br />

看了好几遍，运行了一下还是摸不着头脑，seqStr经过strsplit后是list 吗？为什么要加上[[1]] 呢？

<br />
reverseAndComplement <- function(strseq, doReverse = T,<br />
    doComplement = T) {<br />
    if (!doReverse)<br />
        rev <- c<br />
    ifelse(doComplement, paste(c(A = "T", T = "A", C = "G", G = "C")[rev(strsplit(toupper(strseq),<br />
        "")[[1]])], collapse = ""), rev(s))<br />
}<br />
system.time(reverseAndComplement(seqStr))<br />
#   user  system elapsed<br />
#  0.900   0.000   0.902<br />

rev <- c 和 rev(s) 也看不懂？！ @_@

yanlinlin82

幸好被你指出，果然是笔误，导致无法正确处理反转。

改进的代码如下（这次还加了个names的技巧）：

<br />
reverseAndComplement <- function(strseq, doReverse = T, doComplement = T) {<br />
    map <- c('A','C','G','T')<br />
    names(map) <- ifelse(rep(doComplement, length(map)), rev(map), map)<br />
    if (!doReverse) rev <- c<br />
    paste(map[rev(strsplit(toupper(strseq), "")[[1]])], collapse = "")<br />
}<br />
strseq <- "CAGAACATTTAGCGGTAAAA"<br />
reverseAndComplement(strseq)       # "TTTTACCGCTAAATGTTCTG"<br />
reverseAndComplement(strseq, F)    # "GTCTTGTAAATCGCCATTTT"<br />
reverseAndComplement(strseq, T, F) # "AAAATGGCGATTTACAAGAC"<br />
reverseAndComplement(strseq, F, F) # "CAGAACATTTAGCGGTAAAA"<br />

结果看起来似乎很好。

那个“惊世骇俗”的“if (!doReverse) rev <- c”，我在开始写时就特别想笑，这纯粹是钻了R语言的空子，用一个局部的rev变量遮盖掉系统的rev函数，从而使反转失效。此法纯属搞笑，如非万不得已，尤其是在头脑不清醒时，请勿模仿使用，否则后果自负，哈哈 [s:11]

yanlinlin82

是的，seqStr经过strsplit后是list，所以用[[1]]来提取第一项。strsplit其实能够同时切分多个字符串，因此要用列表来存放每个字符串的切分结果。

ryusukekenji

[quote]引用第9楼yanlinlin82于2010-02-20 01:19发表的回 7楼(ryusukekenji) 的帖子 :

是的，seqStr经过strsplit后是list，所以用[[1]]来提取第一项。strsplit其实能够同时切分多个字符串，因此要用列表来存放每个字符串的切分结果。

[/quote]

> test <- c('A','B','C')<br />
> strsplit(test,'')<br />
[[1]]<br />
[1] "A"<br />
<br />
[[2]]<br />
[1] "B"<br />
<br />
[[3]]<br />
[1] "C"<br />

if (!doReverse) rev <- c 哈哈哈哈，了解了。strsplit 也了解了。谢谢yanlinlin82版主。

yihui

R的随意性跟C比起来就像个村夫，C一向都是西装笔挺，R灰头土脸云游四方

yanlinlin82

其实C也有夸张的要命的写法：

http://www.cnblogs.com/exclm/archive/2008/12/18/1357586.html

http://en.wikipedia.org/wiki/International_Obfuscated_C_Code_Contest

随意性与语言无关，而与人有关。所以才有那句话：好的代码是先给人读，然后才给机器读的。

记得老大曾经不时地给出一些关于写R程序的规则，诸如避免循环、避免attach之类的，确实为金玉良言。现在，为避免误导公众，我想“rev<-c”也该被收录列为禁忌之一了。 [s:13]

yihui

发现一个简单的调换字符的方法：

> chartr("ATCG", "TAGC", "CAGAACATTTAGCGGTAAAA")<br />
[1] "GTCTTGTAAATCGCCATTTT"<br />

yanlinlin82

回复第15楼的谢益辉：又领教了！一比较，前面的所有代码都显得很傻很天真了……

ryusukekenji

回复第14楼的 yanlinlin82：

#记得老大曾经不时地给出一些关于写R程序的规则，诸如避免循环、避免attach之类的，确实为金玉良言....

哦，attach原来也会拖慢速度吗？在R-help中还常看到一些编码attach了数据才运算。

回复第15楼的谢益辉：看了一下chartr的代码，才几个字竟然就那么好用，厉害。

yihui

attach()的问题不在于速度，而是它会把命名空间搞乱，可能会给一些人带来混淆。例如：

> x = sleep<br />
> attach(x)<br />
# 对x中的变量操作<br />
> extra = extra + 1<br />
# 但x本身的变量并没有受影响<br />
> identical(x$extra, extra)<br />
[1] FALSE<br />
> detach()<br />

为了不要把环境弄乱，一般尽量用with(data, {...})，比如

with(sleep, plot(extra ~ group))

这样能拥有attach()的优点（不必写dataname$varname，只需写varname），而且不会对命名空间产生影响。

当然，萝卜白菜各有所爱，也有很多人愿意用attach()，怎么方便怎么来。

ryusukekenji

回复第18楼的谢益辉：哦，原来是这样子，看到一些洋人研究人员都爱使用attach。这样的话不就出现另一个问题了？举个下列的例子，那么我们是否就该避免使用同样的变量名了？

>#例子<br />
>> x <- c(10:15)<br />
> smp <- data.frame(x=c(1:5),y=c(6:10))<br />
> x<br />
[1] 10 11 12 13 14 15<br />
> attach(smp)</p>
<p>        The following object(s) are masked _by_ .GlobalEnv :</p>
<p>         x<br />

guo.bailing

回复第18楼的谢益辉：正被这个问题困扰，感谢谢老大

月珥

完了我觉得我的C思想已经根深蒂固了。看到这贴首先想到的是：

直接把每个ATCG直接储存成两位二进制数，比如说令A=00，T=10，C=01，G=11，然后合并起来，得到一个00100111...这样的二进制链，最后每隔一位取反就行……

不知R里头做类似的工作方便不……

suckbunny

回复第21楼的月珥：用1减偶数项就好了。

winthander

请问一下给位老大，谢老大在三楼的程序中

[seqStr = paste(sample(c("A", "T", "C", "G"), 3e+06, TRUE), collapse = "")

seqChar = strsplit(seqStr, "")[[1]]]

先paste后strsplit得到list,

如果直接使用[sample(c("A", "T", "C", "G"), 3e+06,

TRUE)]得到的向量好像也可以用于分析，为什么呢

aubonj

来拜读一下，学习了。

可是还是有些云里雾里的，各位老大能推荐一下哪里可以补补课，怎样将显性循环替换为高效的R语句。。。

我经常是写出好几个for的嵌套，或者需要进行场景的步进模拟和monte carlos，于是那个运算就非常耗时了。。。先谢谢了！！[s:19]

enthumelon

回复第21楼的月珥：你的思路我总结下：

1. ACGT的字符编码就用十进制表示: 0,1,2,3（文本替换就行了，纯粹R的操作嘛，干脆转为factor）;

2. 然后转码后用sfsmisc:::digitsBase做成0-1向量合并成矩阵;

3 用1减去偶数项再转义成ACGT。

这样做的话，在C里面做取反自然快（算法复杂度在O(N)），在R里面嘛...不是个好思路。

因为无非就是字符转换，基因文件处理时（如果太大了）可以分段将文件读入内存（一次放不进内存的情况多了去了），然后做脚标操作即可(虽然算法复杂度也是O(N),但肯定N钱的常数要大些，因为要记录和替换脚标)。