तेजी से आर

Question

में पैदा ~ एक यादृच्छिक प्रक्रिया के 10^9 चरणों मैं एक निम्न कार्य करने के लिए है:

X(0)=0 
X(t+1)=X(t)+Y(t) 

:

प्रक्रिया सूत्र द्वारा वर्णित के 10^9 चरणों उत्पन्न जहां Y(t) वितरण N(0,1) साथ स्वतंत्र यादृच्छिक चर रहे हैं। tX(t) का मान नकारात्मक था।

x<-c(0,0) 
    z<-0 
    loop<-10^9 
    for(i in 2:loop) { 
    x[1]<-x[2] 
    x[2]<-x[1]+rnorm(1, 0, 1) 
    if (x[2]<0) {z<-z+1} 
    } 

हालांकि, यह बहुत धीमी है:

मैं निम्नलिखित कोड की कोशिश की। मैं इसे कैसे बढ़ा सकता हूं?

Answer 1

एक समाधान vectorized @ G5W द्वारा प्रस्तावित के साथ जाना है, लेकिन यह छोटे टुकड़ों में तोड़ने के लिए किसी भी स्मृति अतिप्रवाह मुद्दों से बचने के है। यह आपको वेक्टरकृत समाधान की गति देता है, लेकिन खंड आकार के प्रबंधन से आप नियंत्रित कर सकते हैं कि प्रक्रिया कितनी मेमोरी का उपयोग करती है।

निम्नलिखित 1e + 07 के ब्लॉक में समस्या को तोड़ता है, और 100 बार लूप करके आपको कुल 1e + 09 मिलता है।

पहले ब्लॉक के अंत में, आप 0 से नीचे का समय और अंतिम बिंदु रिकॉर्ड करते हैं। अंत बिंदु को फिर अगले ब्लॉक में खिलाया जाता है, और आप 0 से नीचे का समय, और नया समापन बिंदु रिकॉर्ड करते हैं।

अंत में शून्य से कुल समय प्राप्त करने के लिए 100 रन औसत। जबकि लूप में cat पर कॉल प्रगति की निगरानी करने और प्रगति को देखने के लिए हैं, इस पर टिप्पणी की जा सकती है।

funky <- function(start, length = 1e+07) { 
    Y <- rnorm(length) 
    Z <- cumsum(Y) 
    c(sum(Z<(-start))/length, (tail(Z, 1) + start)) 
} 

starttime <- Sys.time() 
resvect <- vector(mode = "numeric", length = 100) 
result <- funky(0) 
resvect[1] <- result[1] 
i <- 2 
while (i < 101) { 
    cat(result, "\n") 
    result <- funky(result[2]) 
    resvect[i] <- result[1] 
    i <- i + 1 
} 
mean(resvect) 
# [1] 0.1880392 
endtime <- Sys.time() 
elapsed <- endtime - starttime 
elapsed 
# Time difference of 1.207566 mins 

Answer 2

यह बहुत तेजी से होना चाहिए, लेकिन कुछ भी की एक अरब कुछ समय लग सकता। लंबाई के छोटे मूल्यों के साथ इसका परीक्षण करना अच्छा हो सकता है - जैसे 10^6।

length = 10^9 
Y = rnorm(length) 
sum(cumsum(Y)<0)/length 

संपादित

@ की टिप्पणी user3666197 मैं इस परीक्षण किया है और वह सही था पर आधारित है। यह समाधान में छोटी संख्याओं के लिए अच्छी तरह से काम, लेकिन एक बार चरणों की संख्या बहुत बड़ी होने के लिए हो जाता है, यह विफल रहता है।

मैं ओपी के कोड के खिलाफ मेरी "vectorized" संस्करण का परीक्षण किया। जब यादृच्छिक चलने की लंबाई 10^8 थी, तो मेरे कोड में लगभग 7 सेकंड लग गए और ओपी के कोड में 131 सेकंड (मेरे लैपटॉप पर) लिया गया। हालांकि, जब मैं 10^9 लंबाई बढ़ (मूल प्रश्न के अनुसार), मेरी संस्करण डिस्क स्वैपिंग का एक बहुत का कारण बना और मैं इस प्रक्रिया को मारने के लिए किया था। यह समाधान ओपी द्वारा अनुरोधित पैमाने पर विफल रहता है।

Answer 3

यादृच्छिकता के स्रोत को देखते हुए तकनीकी रूप से जेनरेट की गई स्ट्रीम की पुनरावर्तनीयता के लिए एक आवश्यकता दोनों को पूरा करने के लिए अन्यथा निर्धारिती हार्डवेयर की क्षमता के रूप में बनाया गया है और "जेनरेट" के लिए सभी स्थितियां- एक निश्चित छद्म-यादृच्छिक जनरेटर द्वारा यादृच्छिकता एल्गोरिदम, यादृच्छिकता का स्रोत आसानी से शुद्ध-[SERIAL] से "बस" - [CONCURRENT] या एक सत्य- [PARALLEL] मोडस ऑपरंदी के किसी भी रूप में परिवर्तनीय नहीं है।

यह कहा, पीआरजी कदम एक pure- [SERIAL] कोड-निष्पादन को फिर से परिभाषित करने के लिए किसी भी प्रयास के लिए केंद्रीय बिंदु (अवरुद्ध) है।

यह (गैर) -नकारात्मक X(t) -values ​​का प्रतिशत परिवर्तन नहीं करता है, लेकिन सिर्फ एक दिया पीआरजी-हार्डवेयर कार्यान्वयन के लिए, वहाँ कोई छोटा रास्ता नहीं है, लेकिन परस्पर की पीढ़ी के pure- [SERIAL] अनुक्रम (क्रमानुसार कि निर्धारित करता है) निर्भर मूल्य।

"धीमी" पाश या अर्ध unrolling -vectorised procesing (आर-भाषा कार्यान्वयन सुविधाओं का दोहन लेकिन लगभग हार्डवेयर सीपीयू अनुदेश सेट स्तर चाल (मूल्यों के रूप में _{अभी भी क्रमानुसार निर्भर} कर रहे हैं) - तो नहीं एक भाषा खेल -चेंजर, लेकिन कुछ जानबूझकर धीमे कोड-निष्पादन कन्स्ट्रक्टरों को बाधित करने का थोड़ा सा) सबसे अधिक होने की उम्मीद कर सकते हैं।

Answer 4

सामान्य में, इस तरह की समस्याओं के लिए, आप एक-से-एक सी में ++ Rcpp पैकेज का उपयोग कर अपने समारोह अनुवाद कर सकते हैं। यह एक काफी गति देना चाहिए।

पहले, आर संस्करण:

random_sum <- function(loop = 1000) { 
    x<-c(0,0) 
    z<-0 
    for(i in 2:loop) { 
    x[1]<-x[2] 
    x[2]<-x[1]+rnorm(1, 0, 1) 
    if (x[2]<0) {z<-z+1} 
    } 
    z/loop 
} 
set.seed(123) 
random_sum() 
# [1] 0.134 

अब सी ++ संस्करण:

library("Rcpp") 
cppFunction(" 
    double random_sum_cpp(unsigned long loop = 1000) { 
    double x1 = 0; 
    double x2 = 0; 
    double z = 0; 
    for (unsigned long i = 2; i < loop; i++) { 
     x1 = x2; 
     x2 = x1 + Rcpp::rnorm(1)[0]; 
     if (x2 < 0) z = z+1; 
    } 
    return z/loop; 
    }") 

set.seed(123) 
random_sum_cpp() 
# [1] 0.134 

पूर्णता के लिए, हम भी vectorized संस्करण है कि प्रस्तावित किया गया था पर विचार करते हैं:

random_sum_vector <- function(loop = 1000) { 
    Y = rnorm(loop) 
    sum(cumsum(Y)<0)/loop 
} 
set.seed(123) 
random_sum_vector() 
# [1] 0.134 

हम देखें कि यह एक ही यादृच्छिक बीज के लिए एक ही परिणाम देता है, इसलिए यह एक व्यवहार्य दावेदार लगता है।

बेंचमार्क में, सी ++ संस्करण और vectorized संस्करण सी पर एक मामूली बढ़त प्रदर्शित vectorized संस्करण के साथ, इसी तरह प्रदर्शन ++ संस्करण:

> microbenchmark(random_sum(100000), 
       random_sum_vector(100000), 
       random_sum_cpp(100000)) 
Unit: milliseconds 
        expr  min   lq  mean  median   uq  max neval 
     random_sum(1e+05) 184.205588 199.859266 209.220232 205.137043 211.026740 274.47615 100 
random_sum_vector(1e+05) 6.320690 6.631704 7.273645 6.799093 7.334733 18.48649 100 
    random_sum_cpp(1e+05) 8.950091 9.362303 10.663295 9.956996 11.079513 21.30898 100 

हालांकि, vectorized संस्करण स्मृति और will blow up your memory for long loops. साथ गति बंद कारोबार सी ++ संस्करण वर्चुअल रूप से कोई स्मृति का उपयोग नहीं करता है।

10^9 चरणों के लिए, सी ++ संस्करण मेरी मशीन पर लगभग 2 मिनट (110 सेकंड) में चलता है। मैंने आर संस्करण की कोशिश नहीं की। छोटे बेंचमार्क के आधार पर, शायद इसमें लगभग 7 घंटे लगेंगे।

> microbenchmark(random_sum_cpp(10^9), times = 1) 
Unit: seconds 
       expr  min  lq  mean median  uq  max neval 
random_sum_cpp(10^9) 110.2182 110.2182 110.2182 110.2182 110.2182 110.2182  1 

Answer 5

का उपयोग वैक्टर आम तौर पर छोरों के लिए की तुलना में बेहतर प्रदर्शन निकलेगा। बहुत बड़ी संख्या के साथ यहां मुद्दा (यानी 10^9) मेमोरी सीमाएं हैं। चूंकि आप केवल ऋणात्मक सूचकांक के अंतिम प्रतिशत में रुचि रखते हैं जो नकारात्मक हैं, निम्नलिखित कार्य करेंगे (10^9 चरणों में कुछ मिनट लगते हैं)।

update_state <- function (curr_state, step_size) { 
    n <- min(curr_state$counter, step_size) 
    r <- rnorm(min(curr_state$counter, step_size)) 
    total <- curr_state$cum_sum + cumsum(r) 

    list('counter' = curr_state$counter - n, 
     'neg_count' = curr_state$neg_count + length(which(total < 0)), 
     'cum_sum' = curr_state$cum_sum + sum(r)) 
} 


n <- 10^9 
curr_state <- list('counter' = n, 'neg_count' = 0, 'cum_sum' = 0) 

step_size <- 10^8 
while (curr_state$counter > 0) { 
    curr_state <- update_state(curr_state = curr_state, step_size = step_size) 
} 

print(curr_state) 
print(curr_state$neg_count/ n)