आर - एक परिवर्तनीय अंतराल पर रोलिंग आंकड़ों की गणना करने के लिए तेज़ तरीका

Question

मुझे उत्सुकता है कि अगर कोई वहां घुमावदार आंकड़ों (रोलिंग माध्य, औसत, प्रतिशत, आदि) की गणना करने के लिए एक (तेज़) तरीके से आ सकता है समय के अंतराल (खिड़की)।

ऐसा लगता है कि किसी को यादृच्छिक रूप से समय के अवलोकन दिए जाते हैं (यानी दैनिक, या साप्ताहिक डेटा नहीं, अवलोकनों में समय की टिकटें होती है, जैसे कि टिक डेटा में), और मान लीजिए कि आप केंद्र और फैलाव आंकड़े देखना चाहते हैं आप समय के अंतराल को चौड़ा और कसने में सक्षम हैं जिन पर इन आंकड़ों की गणना की जाती है।

मैंने लूप के लिए एक सरल बनाया जो ऐसा करता है। लेकिन यह स्पष्ट रूप से बहुत धीमा चलता है (वास्तव में मुझे लगता है कि मेरा लूप अभी भी डेटा की एक छोटी नमूना पर चल रहा है जिसे मैंने अपनी गति का परीक्षण करने के लिए स्थापित किया है)। मैं ऐसा करने के लिए ddply की तरह कुछ पाने की कोशिश कर रहा हूं - जो दैनिक आंकड़ों के लिए दौड़ने के लिए कठोर लगता है - लेकिन मैं इसके बाहर अपना रास्ता काम नहीं कर सकता।

उदाहरण:

नमूना सेट अप:!

df <- data.frame(Date = runif(1000,0,30)) 
df$Price <- I((df$Date)^0.5 * (rnorm(1000,30,4))) 
df$Date <- as.Date(df$Date, origin = "1970-01-01") 

उदाहरण समारोह (जो वास्तव में कई टिप्पणियों के साथ धीमी गति से चलाता है

SummaryStats <- function(dataframe, interval){ 
    # Returns daily simple summary stats, 
    # at varying intervals 
    # dataframe is the data frame in question, with Date and Price obs 
    # interval is the width of time to be treated as a day 

    firstDay <- min(dataframe$Date) 
    lastDay <- max(dataframe$Date) 
    result <- data.frame(Date = NULL, 
         Average = NULL, Median = NULL, 
         Count = NULL, 
         Percentile25 = NULL, Percentile75 = NULL) 

    for (Day in firstDay:lastDay){ 

    dataframe.sub = subset(dataframe, 
       Date > (Day - (interval/2)) 
       & Date < (Day + (interval/2))) 

    nu = data.frame(Date = Day, 
        Average = mean(dataframe.sub$Price), 
        Median = median(dataframe.sub$Price), 
        Count = length(dataframe.sub$Price), 
        P25 = quantile(dataframe.sub$Price, 0.25), 
        P75 = quantile(dataframe.sub$Price, 0.75)) 

    result = rbind(result,nu) 

    } 

    return(result) 

} 

आपका सलाह स्वागत किया जाएगा

Answer 1

Rcpp यदि गति अपनी प्राथमिक चिंता का विषय है एक अच्छा तरीका है। उदाहरण के द्वारा व्याख्या करने के लिए मैं रोलिंग माध्य आंकड़े का उपयोग करूंगा।

मानक: आर

x = sort(runif(25000,0,4*pi)) 
y = sin(x) + rnorm(length(x),0.5,0.5) 
system.time(rollmean_r(x,y,xout=x,width=1.1)) # ~60 seconds 
system.time(rollmean_cpp(x,y,xout=x,width=1.1)) # ~0.0007 seconds 

Rcpp के लिए कोड और आर समारोह

cppFunction(' 
    NumericVector rollmean_cpp(NumericVector x, NumericVector y, 
           NumericVector xout, double width) { 
    double total=0; 
    unsigned int n=x.size(), nout=xout.size(), i, ledge=0, redge=0; 
    NumericVector out(nout); 

    for(i=0; i<nout; i++) { 
     while(x[ redge ] - xout[i] <= width && redge<n) 
     total += y[redge++]; 
     while(xout[i] - x[ ledge ] > width && ledge<n) 
     total -= y[ledge++]; 
     if(ledge==redge) { out[i]=NAN; total=0; continue; } 
     out[i] = total/(redge-ledge); 
    } 
    return out; 
    }') 

rollmean_r = function(x,y,xout,width) { 
    out = numeric(length(xout)) 
    for(i in seq_along(xout)) { 
    window = x >= (xout[i]-width) & x <= (xout[i]+width) 
    out[i] = .Internal(mean(y[window])) 
    } 
    return(out) 
} 

rollmean_cpp के explantion के लिए अब बनाम Rcpp। x और y डेटा हैं। xout उन बिंदुओं का एक वेक्टर है जिस पर रोलिंग आंकड़े का अनुरोध किया जाता है। width रोलिंग विंडो की चौड़ाई * 2 है। ध्यान दें कि स्लाइडिंग विंडो के सिरों के लिए indeces ledge और redge में संग्रहीत हैं। ये अनिवार्य रूप से x और y में संबंधित तत्वों के लिए संकेतक हैं। ये इंडेस अन्य सी ++ फ़ंक्शंस (उदाहरण के लिए, औसत और पसंद) को कॉल करने के लिए बहुत फायदेमंद हो सकते हैं जो एक वेक्टर लेते हैं और इनपुट के रूप में इंडेस शुरू करते हैं और समाप्त करते हैं।

जो लोग डिबगिंग (लंबा) के लिए rollmean_cpp की एक "वर्बोज़" संस्करण चाहते हैं के लिए:

cppFunction(' 
    NumericVector rollmean_cpp(NumericVector x, NumericVector y, 
           NumericVector xout, double width) { 

    double total=0, oldtotal=0; 
    unsigned int n=x.size(), nout=xout.size(), i, ledge=0, redge=0; 
    NumericVector out(nout); 


    for(i=0; i<nout; i++) { 
     Rcout << "Finding window "<< i << " for x=" << xout[i] << "..." << std::endl; 
     total = 0; 

     // numbers to push into window 
     while(x[ redge ] - xout[i] <= width && redge<n) { 
     Rcout << "Adding (x,y) = (" << x[redge] << "," << y[redge] << ")" ; 
     Rcout << "; edges=[" << ledge << "," << redge << "]" << std::endl; 
     total += y[redge++]; 
     } 

     // numbers to pop off window 
     while(xout[i] - x[ ledge ] > width && ledge<n) { 
     Rcout << "Removing (x,y) = (" << x[ledge] << "," << y[ledge] << ")"; 
     Rcout << "; edges=[" << ledge+1 << "," << redge-1 << "]" << std::endl; 
     total -= y[ledge++]; 
     } 
     if(ledge==n) Rcout << " OVER "; 
     if(ledge==redge) { 
     Rcout<<" NO DATA IN INTERVAL " << std::endl << std::endl; 
     oldtotal=total=0; out[i]=NAN; continue;} 

     Rcout << "For interval [" << xout[i]-width << "," << 
       xout[i]+width << "], all points in interval [" << x[ledge] << 
       ", " << x[redge-1] << "]" << std::endl ; 
     Rcout << std::endl; 

     out[i] = (oldtotal + total)/(redge-ledge); 
     oldtotal=total+oldtotal; 
    } 
    return out; 
    }') 

x = c(1,2,3,6,90,91) 
y = c(9,8,7,5.2,2,1) 
xout = c(1,2,2,3,6,6.1,13,90,100) 
a = rollmean_cpp(x,y,xout=xout,2) 
# Finding window 0 for x=1... 
# Adding (x,y) = (1,9); edges=[0,0] 
# Adding (x,y) = (2,8); edges=[0,1] 
# Adding (x,y) = (3,7); edges=[0,2] 
# For interval [-1,3], all points in interval [1, 3] 
# 
# Finding window 1 for x=2... 
# For interval [0,4], all points in interval [1, 3] 
# 
# Finding window 2 for x=2... 
# For interval [0,4], all points in interval [1, 3] 
# 
# Finding window 3 for x=3... 
# For interval [1,5], all points in interval [1, 3] 
# 
# Finding window 4 for x=6... 
# Adding (x,y) = (6,5.2); edges=[0,3] 
# Removing (x,y) = (1,9); edges=[1,3] 
# Removing (x,y) = (2,8); edges=[2,3] 
# Removing (x,y) = (3,7); edges=[3,3] 
# For interval [4,8], all points in interval [6, 6] 
# 
# Finding window 5 for x=6.1... 
# For interval [4.1,8.1], all points in interval [6, 6] 
# 
# Finding window 6 for x=13... 
# Removing (x,y) = (6,5.2); edges=[4,3] 
# NO DATA IN INTERVAL 
# 
# Finding window 7 for x=90... 
# Adding (x,y) = (90,2); edges=[4,4] 
# Adding (x,y) = (91,1); edges=[4,5] 
# For interval [88,92], all points in interval [90, 91] 
# 
# Finding window 8 for x=100... 
# Removing (x,y) = (90,2); edges=[5,5] 
# Removing (x,y) = (91,1); edges=[6,5] 
# OVER NO DATA IN INTERVAL 

print(a) 
# [1] 8.0 8.0 8.0 8.0 5.2 5.2 NaN 1.5 NaN 

Answer 2

चलो देखते हैं ... आप लूप (आर में बहुत धीमी) कर रहे हैं, जिससे सबसेट बनाने में डेटा की अनावश्यक प्रतियां और डेटा सेट जमा करने के लिए rbind का उपयोग किया जा रहा है। यदि आप उनसे बचते हैं, तो चीजें काफी तेज हो जाएंगी। इस प्रयास करें ...

Summary_Stats <- function(Day, dataframe, interval){ 
    c1 <- dataframe$Date > Day - interval/2 & 
     dataframe$Date < Day + interval/2 
    c(
     as.numeric(Day), 
     mean(dataframe$Price[c1]), 
     median(dataframe$Price[c1]), 
     sum(c1), 
     quantile(dataframe$Price[c1], 0.25), 
     quantile(dataframe$Price[c1], 0.75) 
    ) 
} 
Summary_Stats(df$Date[2],dataframe=df, interval=20) 
firstDay <- min(df$Date) 
lastDay <- max(df$Date) 
system.time({ 
    x <- sapply(firstDay:lastDay, Summary_Stats, dataframe=df, interval=20) 
    x <- as.data.frame(t(x)) 
    names(x) <- c("Date","Average","Median","Count","P25","P75") 
    x$Date <- as.Date(x$Date) 
}) 
dim(x) 
head(x) 

Answer 3

जवाब में मेरे सवाल का "केविन" करने के लिए ऊपर, मुझे लगता है कि मैं बाहर से नीचे कुछ लगा।

यह फ़ंक्शन डेटा को लेता है (टाइम अवलोकन यादृच्छिक अंतराल पर आते हैं और एक समय टिकट द्वारा इंगित होते हैं) और एक अंतराल पर माध्य की गणना करता है।

library(Rcpp) 

cppFunction(' 
    NumericVector rollmean_c2(NumericVector x, NumericVector y, double width, 
           double Min, double Max) { 

double total = 0, redge,center; 
unsigned int n = (Max - Min) + 1, 
        i, j=0, k, ledge=0, redgeIndex; 
NumericVector out(n); 


for (i = 0; i < n; i++){ 
    center = Min + i + 0.5; 
    redge = center - width/2; 
    redgeIndex = 0; 
    total = 0; 

    while (x[redgeIndex] < redge){ 
    redgeIndex++; 
    } 
    j = redgeIndex; 

    while (x[j] < redge + width){ 
    total += y[j++]; 

    } 

    out[i] = total/(j - redgeIndex); 
} 
return out; 

    }') 

# Set up example data 
x = seq(0,4*pi,length.out=2500) 
y = sin(x) + rnorm(length(x),0.5,0.5) 
plot(x,y,pch=20,col="black", 
    main="Sliding window mean; width=1", 
    sub="rollmean_c in red  rollmean_r overlaid in white.") 


c.out = rollmean_c2(x,y,width=1,Min = min(x), Max = max(x)) 
lines(0.5:12.5,c.out,col="red",lwd=3) 

Answer 4

एक श्रृंखला के रूप में जुड़े बिंदुओं में से सब से Think। इस श्रृंखला को ग्राफ के रूप में सोचें, जहां प्रत्येक डेटा पॉइंट नोड होता है। फिर, प्रत्येक नोड के लिए, हम उन सभी अन्य नोड्स को ढूंढना चाहते हैं जो दूरी w या उससे कम दूरी पर हैं। ऐसा करने के लिए, मैं पहले एक मैट्रिक्स उत्पन्न करता हूं जो जोड़ों की दूरी देता है। n वें पंक्ति नोड्स n नोड्स के अलावा दूरी देता है।

# First, some data 
x = sort(runif(25000,0,4*pi)) 
y = sin(x) + rnorm(length(x),0,0.5) 

# calculate the rows of the matrix one by one 
# until the distance between the two closest nodes is greater than w 
# This algorithm is actually faster than `dist` because it usually stops 
# much sooner 
dl = list() 
dl[[1]] = diff(x) 
i = 1 
while(min(dl[[i]]) <= w) { 
    pdl = dl[[i]] 
    dl[[i+1]] = pdl[-length(pdl)] + dl[[1]][-(1:i)] 
    i = i+1 
} 

# turn the list of the rows into matrices 
rarray = do.call(rbind, lapply(dl,inf.pad,length(x))) 
larray = do.call(rbind, lapply(dl,inf.pad,length(x),"right")) 

# extra function 
inf.pad = function(x,size,side="left") { 
    if(side=="left") { 
    x = c(x, rep(Inf, size-length(x))) 
    } else { 
    x = c(rep(Inf, size-length(x)), x) 
    } 
    x 
} 

मैं प्रत्येक विंडो के किनारे को निर्धारित करने के लिए मैट्रिस का उपयोग करता हूं। इस उदाहरण के लिए, मैंने w=2 सेट किया है।

# How many data points to look left or right at each data point 
lookr = colSums(rarray <= w) 
lookl = colSums(larray <= w) 

# convert these "look" variables to indeces of the input vector 
ri = 1:length(x) + lookr 
li = 1:length(x) - lookl 

खिड़कियों परिभाषित के साथ, यह बहुत *apply कार्यों का उपयोग करने के लिए अंतिम जवाब पाने के लिए आसान है।

rolling.mean = vapply(mapply(':',li,ri), function(i) .Internal(mean(y[i])), 1) 

उपरोक्त सभी कोड मेरे कंप्यूटर पर लगभग 50 सेकंड लेते हैं। यह मेरे अन्य उत्तर में rollmean_r फ़ंक्शन से थोड़ा तेज़ है। हालांकि, यहां विशेष रूप से अच्छी बात यह है कि इंडेक्स प्रदान किए जाते हैं। फिर आप *apply फ़ंक्शंस के साथ जो भी आर फ़ंक्शन पसंद करते हैं उसका उपयोग कर सकते हैं। उदाहरण के लिए,

rolling.mean = vapply(mapply(':',li,ri), 
             function(i) .Internal(mean(y[i])), 1) 

में लगभग 5 सेकंड लगते हैं। और,

rolling.median = vapply(mapply(':',li,ri), 
             function(i) median(y[i]), 1) 

में लगभग 14 सेकंड लगते हैं। यदि आप चाहते थे, तो आप indeces प्राप्त करने के लिए मेरे अन्य उत्तर में आरसीपीपी फ़ंक्शन का उपयोग कर सकते हैं।