एक स्ट्रिंग को विभाजित करने और आर का उपयोग कर वर्णों को गिनने का तेज़ तरीका?

Question

मैं एक फास्टा फ़ाइल से पढ़ने वाले डीएनए तारों के लिए जीसी सामग्री की गणना करने के लिए एक तेज़ तरीका ढूंढ रहा हूं। यह एक स्ट्रिंग लेने और 'जी' या 'सी' अक्षर प्रकट होने की संख्या की गिनती करने के लिए उबलता है। मैं विचार करने के लिए पात्रों की सीमा भी निर्दिष्ट करना चाहता हूं।

मेरे पास एक कामकाजी फ़ंक्शन है जो काफी धीमा है, और यह मेरे कोड में बाधा उत्पन्न कर रहा है। यह इस तरह दिखता है:

> a = "GCCCAAAATTTTCCGGatttaagcagacataaattcgagg" 
> Rprof(filename="Rprof.out") 
> for(i in 1:500000){gcCount(a,1,40)}; 
> Rprof(NULL) 
> summaryRprof(filename="Rprof.out") 

        self.time self.pct total.time total.pct 
"gcCount"   77.36  76.8  100.74  100.0 
"=="    18.30  18.2  18.30  18.2 
"strsplit"   3.58  3.6  3.64  3.6 
"+"     1.14  1.1  1.14  1.1 
":"     0.30  0.3  0.30  0.3 
"as.logical"  0.04  0.0  0.04  0.0 
"as.character"  0.02  0.0  0.02  0.0 

$by.total 
       total.time total.pct self.time self.pct 
"gcCount"   100.74  100.0  77.36  76.8 
"=="    18.30  18.2  18.30  18.2 
"strsplit"   3.64  3.6  3.58  3.6 
"+"     1.14  1.1  1.14  1.1 
":"     0.30  0.3  0.30  0.3 
"as.logical"   0.04  0.0  0.04  0.0 
"as.character"  0.02  0.0  0.02  0.0 

$sampling.time 
[1] 100.74 

तेजी से इस कोड को बनाने के लिए किसी भी सलाह:

## 
## count the number of GCs in the characters between start and stop 
## 
gcCount <- function(line, st, sp){ 
    chars = strsplit(as.character(line),"")[[1]] 
    numGC = 0 
    for(j in st:sp){ 
    ##nested ifs faster than an OR (|) construction 
    if(chars[[j]] == "g"){ 
     numGC <- numGC + 1 
    }else if(chars[[j]] == "G"){ 
     numGC <- numGC + 1 
    }else if(chars[[j]] == "c"){ 
     numGC <- numGC + 1 
    }else if(chars[[j]] == "C"){ 
     numGC <- numGC + 1 
    } 
    } 
    return(numGC) 
} 

Rprof चल रहा है मेरा पीछा उत्पादन देता है?

Answer 1

बेहतर सब पर विभाजित नहीं करने के लिए, बस मैचों की संख्या:

gcCount2 <- function(line, st, sp){ 
    sum(gregexpr('[GCgc]', substr(line, st, sp))[[1]] > 0) 
} 

तेजी से परिमाण के एक आदेश है कि ।

एक छोटा सी फ़ंक्शन जो कि वर्णों पर बसता है, अभी तक तीव्रता का एक और क्रम होगा।

Answer 2

यहां एक लूप का उपयोग करने की आवश्यकता नहीं है।

इस प्रयास करें:

gcCount <- function(line, st, sp){ 
    chars = strsplit(as.character(line),"")[[1]][st:sp] 
    length(which(tolower(chars) == "g" | tolower(chars) == "c")) 
} 

Answer 3

एक एक लाइनर:

table(strsplit(toupper(a), '')[[1]]) 

Answer 4

मुझे नहीं पता कि यह कोई तेज़ है, लेकिन आप आर पैकेज seqinR - http://pbil.univ-lyon1.fr/software/seqinr/home.php?lang=eng पर देखना चाहेंगे। यह अनुक्रम विश्लेषण के लिए कई विधियों के साथ एक उत्कृष्ट, सामान्य जैव सूचना विज्ञान पैकेज है। यह सीआरएएन में है (जैसा कि मैं इसे लिखता हूं)।

जीसी सामग्री होगा:

mysequence <- s2c("agtctggggggccccttttaagtagatagatagctagtcgta") 
    GC(mysequence) # 0.4761905 

एक स्ट्रिंग से है यही कारण है कि, आप भी उपयोग कर एक fasta फ़ाइल में पढ़ सकते हैं "read.fasta()"।

Answer 5

stringi पैकेज

> stri_count_fixed("GCCCAAAATTTTCCGG",c("G","C")) 
[1] 3 5 

से इस समारोह का प्रयास करें या आप जी गिनती करने के लिए regex संस्करण का उपयोग कर सकते हैं और जी

> stri_count_regex("GCCCAAAATTTTCCGGggcc",c("G|g|C|c")) 
[1] 12 

या आप पहली बार में कार्य tolower उपयोग कर सकते हैं और फिर stri_count

> stri_trans_tolower("GCCCAAAATTTTCCGGggcc") 
[1] "gcccaaaattttccggggcc" 

समय प्रदर्शन

> microbenchmark(gcCount(x,1,40),gcCount2(x,1,40), stri_count_regex(x,c("[GgCc]"))) 
Unit: microseconds 
          expr  min  lq median  uq  max neval 
       gcCount(x, 1, 40) 109.568 112.42 113.771 116.473 146.492 100 
       gcCount2(x, 1, 40) 15.010 16.51 18.312 19.213 40.826 100 
stri_count_regex(x, c("[GgCc]")) 15.610 16.51 18.912 20.112 61.239 100 

लंबी स्ट्रिंग के लिए एक और उदाहरण।stri_dup स्ट्रिंग एन-बार

> stri_dup("abc",3) 
[1] "abcabcabc" 

आप, देखने के लिए लंबे समय तक अनुक्रम stri_count है इस पोस्ट के लिए सभी के लिए तेजी से :)

> y <- stri_dup("GCCCAAAATTTTCCGGatttaagcagacataaattcgagg",100) 
    > microbenchmark(gcCount(y,1,40*100),gcCount2(y,1,40*100), stri_count_regex(y,c("[GgCc]"))) 
    Unit: microseconds 
           expr  min   lq  median  uq  max neval 
       gcCount(y, 1, 40 * 100) 10367.880 10597.5235 10744.4655 11655.685 12523.828 100 
      gcCount2(y, 1, 40 * 100) 360.225 369.5315 383.6400 399.100 438.274 100 
    stri_count_regex(y, c("[GgCc]")) 131.483 137.9370 151.8955 176.511 221.839 100 

Answer 6

धन्यवाद सकते हैं,

replicates एक स्क्रिप्ट का अनुकूलन करने के लिए है, जिसमें मैं 200 बीपी के 100 एम अनुक्रमों की जीसी सामग्री की गणना करना चाहता हूं, मैंने यहां प्रस्तावित विभिन्न तरीकों का परीक्षण समाप्त कर दिया। केन विलियम्स की विधि ने सर्वोत्तम (2.5 घंटे) का प्रदर्शन किया, जो सेक्विन (3.6 घंटे) से बेहतर था। स्ट्रिंग str_count का उपयोग 1.5 घंटे तक कम हो गया।

अंत में मैंने इसे सी ++ में कोड किया और इसे आरसीपीपी का उपयोग करके बुलाया, जो गणना समय को 10 मिनट तक घटा देता है!

#include <Rcpp.h> 
using namespace Rcpp; 
// [[Rcpp::export]] 
float pGC_cpp(std::string s) { 
    int count = 0; 

    for (int i = 0; i < s.size(); i++) 
    if (s[i] == 'G') count++; 
    else if (s[i] == 'C') count++; 

    float pGC = (float)count/s.size(); 
    pGC = pGC * 100; 
    return pGC; 
} 

कौन सा मैं आर टाइपिंग से फोन:

sourceCpp("pGC_cpp.cpp") 
pGC_cpp("ATGCCC")