आर

Question

में दो बड़े लॉजिकल वैक्टरों को पार करने के लिए सबसे तेज़ तरीका दो लॉजिकल वैक्टर, x और y, लंबाई> 1E8 के लिए, 2x2 क्रॉस सारणी की गणना करने का सबसे तेज़ तरीका क्या है?

मुझे संदेह है कि इसका जवाब सी/सी ++ में लिखना है, लेकिन मुझे आश्चर्य है कि आर में कुछ ऐसा है जो इस समस्या के बारे में पहले से ही काफी स्मार्ट है, क्योंकि यह असामान्य नहीं है।

उदाहरण कोड, 300 एम प्रविष्टियों के लिए (अगर 3E8 बहुत बड़ा है तो एन = 1 ई 8 देने के लिए स्वतंत्र महसूस करें; मैंने कुल आकार केवल 2.5 जीबी (2.4 जीबी) के तहत चुना है। मैंने 0.02 की घनत्व को लक्षित किया है, इसे और अधिक बनाने के लिए दिलचस्प (एक, एक विरल वेक्टर इस्तेमाल कर सकते हैं कि अगर मदद करता है, लेकिन जैसे रूपांतरण समय लग सकता है)

set.seed(0) 
N = 3E8 
p = 0.02 
x = sample(c(TRUE, FALSE), N, prob = c(p, 1-p), replace = TRUE) 
y = sample(c(TRUE, FALSE), N, prob = c(p, 1-p), replace = TRUE) 

कुछ स्पष्ट तरीके:।

1. table
2. bigtabulate
3. सरल लॉजिकल ऑपरेशंस (उदा। sum(x & y))
4. वेक्टर गुणा (बू)
5. data.table
6. ऊपर से कुछ multicore पैकेज से parallel (या नए parallel पैकेज के साथ)

मैं पहली बार में एक चाकू लिया है तीन विकल्प (मेरा जवाब देखें), लेकिन मुझे लगता है कि कुछ बेहतर और तेज होना चाहिए।

मुझे लगता है कि table बहुत धीरे-धीरे काम करता है। bigtabulate तार्किक वैक्टर की एक जोड़ी के लिए ओवरकिल की तरह लगता है। आखिरकार, वेनिला लॉजिकल ऑपरेशंस करना एक क्लज की तरह दिखता है, और यह प्रत्येक वेक्टर को कई बार (3 एक्स? 7 एक्स?) देखता है, यह उल्लेख नहीं करता कि यह प्रसंस्करण के दौरान बहुत सारी अतिरिक्त मेमोरी भरता है, जो एक भारी समय बर्बाद है।

वेक्टर गुणा आमतौर पर एक बुरा विचार है, लेकिन जब वेक्टर स्पैस होता है, तो उसे इस तरह के भंडारण से बाहर लाभ हो सकता है, और फिर वेक्टर गुणा का उपयोग कर लाभ हो सकता है।

N और p भिन्न होने के लिए स्वतंत्र महसूस करें, यदि यह सारणीकरण कार्यों के किसी भी दिलचस्प व्यवहार का प्रदर्शन करेगा। :)

---

अद्यतन 1. मेरा पहला जवाब तीन अनुभवहीन तरीकों पर समय देता है, जो विश्वास table धीमी है का आधार है। फिर भी, एहसास करने की मुख्य बात यह है कि "तार्किक" विधि पूरी तरह से अक्षम है। यह क्या कर रहा है को देखो:

• 4 तार्किक वेक्टर संचालन
• 4 प्रकार रूपांतरण (पूर्णांक या तार्किक एफपी - sum के लिए) तार्किक ऑपरेशन के लिए
• 4 वेक्टर summations
• 8 कार्य (1, सारांश के लिए 1)

न केवल यह, लेकिन यह संकलित या समांतर भी नहीं है। फिर भी, यह अभी भी table के पैंट को धड़कता है।ध्यान दें कि bigtabulate, के साथ एक अतिरिक्त प्रकार का रूपांतरण (1 * cbind...) अभी भी table धड़कता है।

अद्यतन 2. किसी को भी ऐसा न हो कि का कहना है कि आर समर्थन NA में तार्किक वैक्टर, और कहा कि कि इन पार tabulations के लिए प्रणाली में एक रिंच (जो ज्यादातर मामलों में सच है) हो जाएगा, मैं कहना चाहिए कि मेरी वैक्टर आ is.na() या is.finite() से। :) मैं NA और अन्य गैर-सीमित मूल्यों को डिबग कर रहा हूं - they've been a headache for me recently। यदि आप नहीं जानते कि क्या है या अपने प्रविष्टियों के सभी नहीं कर रहे हैं NA, आप any(is.na(yourVector)) साथ परीक्षण कर सकते हैं - यह बुद्धिमान होगा इससे पहले कि आप & ए

---

अद्यतन 3 इस क्यू में उत्पन्न होने वाली कुछ विचारों की अपनाने। ब्रैंडन बर्टेलसन ने टिप्पणियों में एक बहुत ही उचित सवाल पूछा: एक उप-नमूना (प्रारंभिक सेट, आखिरकार, नमूना है ;-)) क्रॉस-टैबलेशन बनाने के प्रयोजनों के लिए पर्याप्त हो सकता है, तो इतना डेटा क्यों उपयोग करें? आंकड़ों में बहुत दूर नहीं जाना है, लेकिन डेटा उन मामलों से उत्पन्न होता है जहां TRUE अवलोकन दोनों चर के लिए बहुत दुर्लभ हैं। एक डेटा विसंगति का परिणाम है, दूसरा कोड में संभावित बग के कारण (संभावित बग क्योंकि हम केवल कम्प्यूटेशनल परिणाम देखते हैं - परिवर्तनीय x को "कचरा इन" और y के रूप में "कचरा आउट" के रूप में सोचें। परिणामस्वरूप सवाल यह है कि क्या कोड के कारण उत्पादन में समस्याएं पूरी तरह से उन मामलों में हैं जहां डेटा असंगत है, या ऐसे कुछ अन्य उदाहरण हैं जहां अच्छा डेटा खराब हो जाता है? (यही कारण है कि मैंने stopping when a NaN, NA, or Inf is encountered के बारे में एक प्रश्न पूछा।)

भी बताते हैं यही कारण है कि मेरी उदाहरण TRUE मूल्यों के लिए एक कम संभावना है, इन वास्तव में काफी समय 0.1% से कम होते हैं

इस एक अलग समाधान पथ दर्शाती है हाँ: यह पता चलता है कि हम दो सूचकांकों का उपयोग कर सकते हैं।?(अर्थात। प्रत्येक सेट में TRUE के स्थान) और गिनती सेट चौराहे। मैंने सेट चौराहे से परहेज किया क्योंकि मुझे थोड़ी देर पहले मैटलैब द्वारा जला दिया गया था (हाँ, यह आर है, लेकिन मेरे साथ भालू), जो पहले एक सेट के तत्वों को छेड़छाड़ करने से पहले सॉर्ट करेगा। ,

Answer 1

यहाँ तार्किक विधि, table के लिए परिणाम हैं, और bigtabulate एन = 3E8 के लिए: (O(n^2) बजाय O(n log n) जैसे मैं थोड़ा याद जटिलता और भी शर्मनाक था।):

  test replications elapsed relative user.self sys.self 
2  logical   1 23.861 1.000000  15.36  8.50 
3 bigtabulate   1 36.477 1.528729  28.04  8.43 
1  table   1 184.652 7.738653 150.61 33.99 

इस मामले में , table एक आपदा है। , इस बिंदु पर

  test replications elapsed relative user.self sys.self 
2  logical   1 0.220 1.000000  0.14  0.08 
3 bigtabulate   1 0.534 2.427273  0.45  0.08 
1  table   1 1.956 8.890909  1.87  0.09 

ऐसा लगता है कि एक ही तार्किक कार्यों लेखन सबसे अच्छा है भले ही उस हनन sum, और प्रत्येक तार्किक वेक्टर कई बार जांच करता है:

तुलना के लिए, यहां एन = 3E6 है। मैंने अभी तक कार्यों को संकलित करने की कोशिश नहीं की है, लेकिन इससे बेहतर परिणाम मिलेंगे।

अद्यतन 1 अगर हम bigtabulate मूल्यों कि पहले से ही पूर्णांक हैं, अर्थात, अगर हम जैसे रूपांतरण bigtabulate की 1 * cbind(v1,v2) बाहर कर देते हैं, तो एन = 3E6 कई 1.80, 2.4 के बजाय है। "तार्किक" विधि से संबंधित एन = 3 ई 8 एकाधिक 1.53 के बजाय केवल 1.21 है।

---

अद्यतन 2

के रूप में यहोशू उलरिच बताया गया है, थोड़ा वैक्टर में परिवर्तित करने में एक महत्वपूर्ण सुधार है - हम आवंटन और एक बहुत कम डेटा चारों ओर जा रहे हैं: आर के तार्किक वैक्टर प्रति 4 बाइट का उपभोग प्रवेश ("क्यों?", आप पूछ सकते हैं ... Well, I don't know, but an answer may turn up here.), जबकि थोड़ा सा वेक्टर खपत करता है, ठीक है, एक बिट, प्रति प्रविष्टि - यानि 1/32 जितना अधिक डेटा। तो, x 1.2e9 बाइट्स का उपभोग करता है, जबकि xb (नीचे दिए गए कोड में थोड़ा संस्करण) केवल 3.75e7 बाइट्स का उपभोग करता है।

मैंने table और bigtabulate अद्यतन बेंचमार्क (एन = 3e8) से भिन्नताएं गिरा दी हैं। ध्यान दें कि logicalB1 मानता है कि डेटा पहले से ही थोड़ा वेक्टर है, जबकि logicalB2 टाइप रूपांतरण के लिए दंड के साथ एक ही ऑपरेशन है। चूंकि मेरे तार्किक वैक्टर अन्य डेटा पर संचालन के परिणाम हैं, इसलिए मुझे थोड़ा वेक्टर से शुरू करने का लाभ नहीं है। फिर भी, भुगतान करने का जुर्माना अपेक्षाकृत छोटा है। ["तार्किक 3" श्रृंखला केवल 3 तार्किक परिचालन करती है, और फिर एक घटाव करता है। चूंकि यह पार सारणीकरण है, हम कुल पता है, के रूप में DWin टिप्पणी की है।]

 test replications elapsed relative user.self sys.self 
4 logical3B1   1 1.276 1.000000  1.11  0.17 
2 logicalB1   1 1.768 1.385580  1.56  0.21 
5 logical3B2   1 2.297 1.800157  2.15  0.14 
3 logicalB2   1 2.782 2.180251  2.53  0.26 
1 logical   1 22.953 17.988245  15.14  7.82 

हम अब यह और भी कई सकल अक्षमताओं के साथ ही 1.8-2.8 सेकंड ले रही है, को तेज कर दिया है। कोई संदेह नहीं है इसे 1 सेकंड से कम में करने के लिए व्यवहार्य होना चाहिए, इसमें से एक या अधिक: सी कोड, संकलन, और मल्टीकोर प्रोसेसिंग शामिल हैं। सभी 3 (या 4) विभिन्न लॉजिकल ऑपरेशंस स्वतंत्र रूप से किए जा सकते हैं, भले ही यह अभी भी गणना चक्रों का अपशिष्ट है।

सबसे अच्छे चुनौतीपूर्ण, logical3B2 का सबसे समान, table से लगभग 80X तेज है। यह निष्क्रिय तार्किक ऑपरेशन से लगभग 10 एक्स तेज है। और इसमें अभी भी सुधार के लिए बहुत सी जगह है।

---

उपर्युक्त उत्पादन करने के लिए यहां कोड है। नोट मैं कुछ ऑपरेशन या वैक्टरों पर टिप्पणी करने की सलाह देता हूं, जब तक आपके पास बहुत सी रैम न हो - x, x1, और xb का निर्माण y ऑब्जेक्ट्स के साथ, मेमोरी का एक उचित हिस्सा उठाएगा।

इसके अलावा, ध्यान दें: 1 के बजाय bigtabulate के लिए पूर्णांक गुणक के रूप में मुझे 1L का उपयोग करना चाहिए था। कुछ बिंदु पर मैं इस परिवर्तन के साथ फिर से चलाऊंगा, और bigtabulate दृष्टिकोण का उपयोग करने वाले किसी भी व्यक्ति को उस परिवर्तन में अनुशंसा करता हूं।

library(rbenchmark) 
library(bigtabulate) 
library(bit) 

set.seed(0) 
N <- 3E8 
p <- 0.02 

x <- sample(c(TRUE, FALSE), N, prob = c(p, 1-p), replace = TRUE) 
y <- sample(c(TRUE, FALSE), N, prob = c(p, 1-p), replace = TRUE) 

x1 <- 1*x 
y1 <- 1*y 

xb <- as.bit(x) 
yb <- as.bit(y) 

func_table <- function(v1,v2){ 
    return(table(v1,v2)) 
} 

func_logical <- function(v1,v2){ 
    return(c(sum(v1 & v2), sum(v1 & !v2), sum(!v1 & v2), sum(!v1 & !v2))) 
} 

func_logicalB <- function(v1,v2){ 
    v1B <- as.bit(v1) 
    v2B <- as.bit(v2) 
    return(c(sum(v1B & v2B), sum(v1B & !v2B), sum(!v1B & v2B), sum(!v1B & !v2B))) 
} 

func_bigtabulate <- function(v1,v2){ 
    return(bigtabulate(1*cbind(v1,v2), ccols = c(1,2))) 
} 

func_bigtabulate2 <- function(v1,v2){ 
    return(bigtabulate(cbind(v1,v2), ccols = c(1,2))) 
} 

func_logical3 <- function(v1,v2){ 
    r1 <- sum(v1 & v2) 
    r2 <- sum(v1 & !v2) 
    r3 <- sum(!v1 & v2) 
    r4 <- length(v1) - sum(c(r1, r2, r3)) 
    return(c(r1, r2, r3, r4)) 
} 

func_logical3B <- function(v1,v2){ 
    v1B <- as.bit(v1) 
    v2B <- as.bit(v2) 
    r1 <- sum(v1B & v2B) 
    r2 <- sum(v1B & !v2B) 
    r3 <- sum(!v1B & v2B) 
    r4 <- length(v1) - sum(c(r1, r2, r3)) 
    return(c(r1, r2, r3, r4)) 
} 

benchmark(replications = 1, order = "elapsed", 
    #table = {res <- func_table(x,y)}, 
    logical = {res <- func_logical(x,y)}, 
    logicalB1 = {res <- func_logical(xb,yb)}, 
    logicalB2 = {res <- func_logicalB(x,y)}, 

    logical3B1 = {res <- func_logical3(xb,yb)}, 
    logical3B2 = {res <- func_logical3B(x,y)} 

    #bigtabulate = {res <- func_bigtabulate(x,y)}, 
    #bigtabulate2 = {res <- func_bigtabulate2(x1,y1)} 
) 

Answer 2

आप विशाल तार्किक वैक्टर पर कार्रवाई का एक बहुत कुछ कर रहे हैं, bit पैकेज पर एक नज़र डालें। यह बूलियन को सही 1-बिट बूलियन के रूप में संग्रहीत करके स्मृति का एक टन बचाता है।

यह table के साथ मदद नहीं करता है; यह वास्तव में इसे और भी खराब बनाता है क्योंकि बिट वेक्टर में इसके अद्वितीय निर्माण के कारण अधिक अद्वितीय मूल्य हैं। लेकिन यह वास्तव में तार्किक तुलना के साथ मदद करता है।

# N <- 3e7 
require(bit) 
xb <- as.bit(x) 
yb <- as.bit(y) 
benchmark(replications = 1, order = "elapsed", 
    bit = {res <- func_logical(xb,yb)}, 
    logical = {res <- func_logical(x,y)} 
) 
#  test replications elapsed relative user.self sys.self user.child sys.child 
# 1  bit   1 0.129 1.00000  0.132 0.000   0   0 
# 2 logical   1 3.677 28.50388  2.684 0.928   0   0 

Answer 3

एक अलग रणनीति लाभ यह है कि नमूने बहुत पक्षपाती (यानी ज्यादातर FALSE) ले, TRUE मूल्यों का सूचकांक का उपयोग कर चौराहों सिर्फ सेट पर विचार करना है।

उस अंत में, मैं func_find01 पेश करता हूं और एक अनुवाद जो bit पैकेज (func_find01B) का उपयोग करता है; उपरोक्त उत्तर में दिखाई देने वाले सभी कोड नीचे चिपकाए गए हैं।

मैंने func_find01B का उपयोग करने के लिए भूल गए, पूर्ण एन = 3e8 मूल्यांकन फिर से चलाया; मैं दूसरे पास में, इसके खिलाफ तेज विधियों को फिर से चलाता हूं।

  test replications elapsed relative user.self sys.self 
6 logical3B1   1 1.298 1.000000  1.13  0.17 
4 logicalB1   1 1.805 1.390601  1.57  0.23 
7 logical3B2   1 2.317 1.785054  2.12  0.20 
5 logicalB2   1 2.820 2.172573  2.53  0.29 
2  find01   1 6.125 4.718798  4.24  1.88 
9 bigtabulate2   1 22.823 17.583205  21.00  1.81 
3  logical   1 23.800 18.335901  15.51  8.28 
8 bigtabulate   1 27.674 21.320493  24.27  3.40 
1  table   1 183.467 141.345917 149.01 34.41 

बस "तेज" तरीके:

 test replications elapsed relative user.self sys.self 
3  find02   1 1.078 1.000000  1.03  0.04 
6 logical3B1   1 1.312 1.217069  1.18  0.13 
4 logicalB1   1 1.797 1.666976  1.58  0.22 
2 find01B   1 2.104 1.951763  2.03  0.08 
7 logical3B2   1 2.319 2.151206  2.13  0.19 
5 logicalB2   1 2.817 2.613173  2.50  0.31 
1  find01   1 6.143 5.698516  4.21  1.93 

तो, find01B तरीकों कि पूर्व परिवर्तित बिट वैक्टर का उपयोग नहीं करते, एक स्लिम मार्जिन (2.327 सेकंड बनाम 2.099 सेकंड) द्वारा बीच सबसे तेज है। find02 कहां से आया? मैंने बाद में एक संस्करण लिखा जो पूर्व-गणना वाले बिट वैक्टर का उपयोग करता है। यह अब सबसे तेज़ है।

सामान्यतः, "सूचकांक विधि" दृष्टिकोण का चलने का समय सीमांत & संयुक्त संभावनाओं से प्रभावित हो सकता है। मुझे संदेह है कि यह विशेष रूप से प्रतिस्पर्धी होगा जब संभावनाएं कम होंगी, लेकिन किसी को प्राथमिकता, या उप-नमूना के माध्यम से जानना होगा।

---

अद्यतन 1. मैं भी जोश ओ ब्रायन के सुझाव का समय समाप्त हो गया है, table() के बजाय tabulate() का उपयोग कर। परिणाम, 12 सेकंड बीतने पर, लगभग 2 एक्स find01 और लगभग आधे bigtabulate2 हैं। अब जब कि सर्वोत्तम तरीकों 1 सेकंड तक पहुंच जाते हैं, तो यह भी अपेक्षाकृत धीमी है:

user system elapsed 
7.670 5.140 12.815 

---

कोड:

func_find01 <- function(v1, v2){ 
    ix1 <- which(v1 == TRUE) 
    ix2 <- which(v2 == TRUE) 

    len_ixJ <- sum(ix1 %in% ix2) 
    len1 <- length(ix1) 
    len2 <- length(ix2) 
    return(c(len_ixJ, len1 - len_ixJ, len2 - len_ixJ, 
      length(v1) - len1 - len2 + len_ixJ)) 
} 

func_find01B <- function(v1, v2){ 
    v1b = as.bit(v1) 
    v2b = as.bit(v2) 

    len_ixJ <- sum(v1b & v2b) 
    len1 <- sum(v1b) 
    len2 <- sum(v2b) 

    return(c(len_ixJ, len1 - len_ixJ, len2 - len_ixJ, 
      length(v1) - len1 - len2 + len_ixJ)) 
} 

func_find02 <- function(v1b, v2b){ 
    len_ixJ <- sum(v1b & v2b) 
    len1 <- sum(v1b) 
    len2 <- sum(v2b) 

    return(c(len_ixJ, len1 - len_ixJ, len2 - len_ixJ, 
      length(v1b) - len1 - len2 + len_ixJ)) 
} 

func_bigtabulate2 <- function(v1,v2){ 
    return(bigtabulate(cbind(v1,v2), ccols = c(1,2))) 
} 

func_tabulate01 <- function(v1,v2){ 
    return(tabulate(1L + 1L*x + 2L*y)) 
} 

benchmark(replications = 1, order = "elapsed", 
    table = {res <- func_table(x,y)}, 
    find01 = {res <- func_find01(x,y)}, 
    find01B = {res <- func_find01B(x,y)}, 
    find02 = {res <- func_find01B(xb,yb)}, 
    logical = {res <- func_logical(x,y)}, 
    logicalB1 = {res <- func_logical(xb,yb)}, 
    logicalB2 = {res <- func_logicalB(x,y)}, 

    logical3B1 = {res <- func_logical3(xb,yb)}, 
    logical3B2 = {res <- func_logical3B(x,y)}, 

    tabulate = {res <- func_tabulate(x,y)}, 
    bigtabulate = {res <- func_bigtabulate(x,y)}, 
    bigtabulate2 = {res <- func_bigtabulate2(x1,y1)} 
) 

Answer 4

यहाँ, Rcpp साथ एक जवाब है केवल उन प्रविष्टियों कि दोनों 0 नहीं हैं सारणी । मुझे संदेह है कि इसमें सुधार करने के कई तरीके होंगे, क्योंकि यह असामान्य रूप से धीमा है; यह Rcpp के साथ मेरा पहला प्रयास भी है, इसलिए डेटा को स्थानांतरित करने से जुड़े कुछ स्पष्ट अक्षमताएं हो सकती हैं। मैंने एक उदाहरण लिखा है जो उद्देश्यपूर्ण रूप से सादा वेनिला है, जिसे दूसरों को यह प्रदर्शित करने देना चाहिए कि इसे कैसे सुधार किया जा सकता है।

library(Rcpp) 
library(inline) 
doCrossTab <- cxxfunction(signature(x="integer", y = "integer"), body=' 
    Rcpp::IntegerVector Vx(x); 
    Rcpp::IntegerVector Vy(y); 
    Rcpp::IntegerVector V(3); 
    for(int i = 0; i < Vx.length(); i++) { 
    if((Vx(i) == 1) & (Vy(i) == 1)){ V[0]++; } 
    else if((Vx(i) == 1) & (Vy(i) == 0)){ V[1]++; } 
    else if((Vx(i) == 0) & (Vy(i) == 1)){ V[2]++; } 
} 
    return(wrap(V)); 
    ', plugin="Rcpp") 

N = 3E8 के लिए समय परिणाम:

user system elapsed 
10.930 1.620 12.586 

यह मेरी 2nd जवाब में func_find01B जब तक 6X की तुलना में अधिक समय लगता है।

Answer 5

यहां आरसीपीपी चीनी का उपयोग कर एक उत्तर है।

N <- 1e8 
x <- sample(c(T,F),N,replace=T) 
y <- sample(c(T,F),N,replace=T) 

func_logical <- function(v1,v2){ 
    return(c(sum(v1 & v2), sum(v1 & !v2), sum(!v1 & v2), sum(!v1 & !v2))) 
} 


library(Rcpp) 
library(inline) 

doCrossTab1 <- cxxfunction(signature(x="integer", y = "integer"), body=' 
    Rcpp::LogicalVector Vx(x); 
    Rcpp::LogicalVector Vy(y); 
    Rcpp::IntegerVector V(4); 

    V[0] = sum(Vx*Vy); 
    V[1] = sum(Vx*!Vy); 
    V[2] = sum(!Vx*Vy); 
    V[3] = sum(!Vx*!Vy); 
    return(wrap(V)); 
    ' 
, plugin="Rcpp") 

system.time(doCrossTab1(x,y)) 

require(bit) 
system.time(
{ 
xb <- as.bit(x) 
yb <- as.bit(y) 
func_logical(xb,yb) 
}) 

जिसमें परिणाम:

> system.time(doCrossTab1(x,y)) 
    user system elapsed 
    1.067 0.002 1.069 
> system.time(
+ { 
+ xb <- as.bit(x) 
+ yb <- as.bit(y) 
+ func_logical(xb,yb) 
+ }) 
    user system elapsed 
    1.451 0.001 1.453 

तो, हम एक छोटे से गति को थोड़ा पैकेज से अधिक है, हालांकि मैं कितना प्रतिस्पर्धात्मक समय है पर हैरान हूँ मिल सकती है।

अद्यतन:

doCrossTab2 <- cxxfunction(signature(x="integer", y = "integer"), body=' 
    Rcpp::LogicalVector Vx(x); 
    Rcpp::LogicalVector Vy(y); 
    Rcpp::IntegerVector V(4); 
    V[0]=V[1]=V[2]=V[3]=0; 
    LogicalVector::iterator itx = Vx.begin(); 
    LogicalVector::iterator ity = Vy.begin(); 
    while(itx!=Vx.end()){ 
    V[0] += (*itx)*(*ity); 
    V[1] += (*itx)*(!*ity); 
    V[2] += (!*itx)*(*ity); 
    V[3] += (!*itx)*(!*ity);  
    itx++; 
    ity++; 
    } 
    return(wrap(V)); 
    ' 
, plugin="Rcpp") 

system.time(doCrossTab2(x,y)) 
# user system elapsed 
# 0.780 0.001 0.782