Суммирование по индексации в Rcpp

Question

Пример данных

library(dplyr)
library(data.table)
library(Rcpp)

df <- data.frame(year = rep(2003:2004, each = 4), 
             mg = rep(rep(c("a", "b"), each = 2), times = 2),
             plant_date = c(20, 30, 20, 30, 33, 40, 33, 40),
             stage1 = c(40, 50, 42, 52, 43, 55, 48, 57),
             stage2 = c(55, 65, 57, 66, 58, 68, 59, 65),
             stage3 = c(61, 75, 63, 76, 66, 77, 68, 79))

set.seed(123)
dat <- data.frame(year = rep(2003:2004, each = 365), doy = rep(1:365, times = 2), 
              rainfall = sample(0:20, 730, replace = T))

final.dat <- dat %>% dplyr::left_join(df)

Что я хочу сделать, это для каждого года, взять комбинацию plant.date а также mgи сумма и среднее количество осадков от plant.date в stage1, stage1 вstage2а такжеstage2 в stage3, Я сделал это:

final.dat %>% dplyr::group_by(year, plant_date, mg) %>% 
dplyr::summarise(
sum_rain_stage1 = sum(rainfall[doy >= plant_date & doy <= stage1]),
sum_rain_stage2 = sum(rainfall[doy >= plant_date & doy <= stage2]),
sum_rain_stage3 = sum(rainfall[doy >= plant_date & doy <= stage3]),
mean_rain_stage1 = mean(rainfall[doy >= plant_date & doy <= stage1]),
mean_rain_stage2 = mean(rainfall[doy >= plant_date & doy <= stage2]),
mean_rain_stage3 = mean(rainfall[doy >= plant_date & doy <= stage3]))

Я хочу добиться аналогичного использования Rcpp, Но из-за нехватки знакомства мне удалось развить скелет этого следующим образом:

final.dat <- data.table(final.dat)   
setkey(x, year, plant_date, mg, doy)

Определите функцию (ниже) и создайте ее источник:

Rcpp::sourceCpp("rainfall_sum.cpp")) 

и запустите функцию следующим образом:

 final.dat[, c("sum_rain_stage1","sum_rain_stage2","sum_rain_stage3" , "mean_rain_stage1", , "mean_rain_stage2", "mean_rain_stage3") := rainfall_sum(doy, rainfall, plant_date, stage1, stage2 , stage3), keyby = .(year, plant_date, mg)]

Моя (неполная) функция:

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]

List rainfall_sum(NumericVector doy,
                  NumericVector rainfall,
                  double plant_date,
                  double stage1,
                  double stage2,
                  double stage3) {

 double sum_rain_stage1(n);
 double sum_rain_stage2(n);
 double sum_rain_stage3(n);
 double mean_rain_stage1(n);
 double mean_rain_stage2(n);
 double mean_rain_stage3(n);

  # This part I do not know how to sum/average rainfall from plant_date to stage1, stage1 to stage2 and stage2 to stage3
}

 return Rcpp::List::create(Rcpp::Named("sum_rain_stage1") = sum_rain_stage1,
                      Rcpp::Named("sum_rain_stage2") = sum_rain_stage2,
                      Rcpp::Named("sum_rain_stage3") = sum_rain_stage3,
                      Rcpp::Named("mean_rain_stage1") = mean_rain_stage1,
                      Rcpp::Named("mean_rain_stage2") = mean_rain_stage2,
                      Rcpp::Named("mean_rain_stage3") = mean_rain_stage3);
 }

Answer 1

Чтобы ответить на вопрос, здесь есть функция суммирования элементов, оставьте в качестве упражнения для вычисления среднего значения в C++.

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
double sumRcpp(NumericVector x){
  int n = x.size(); // Size of vector
  double sum = 0; // Sum value
  // For loop, note cpp index shift to 0
  for(int i = 0; i < n; i++){
    // Shorthand for sum = sum + x[i]
    sum += x[i];
  }
  return sum; // Obtain and return the Sum
}

// You can include R code blocks in C++ files processed with sourceCpp
// (useful for testing and development). The R code will be automatically 
// run after the compilation.
//

/*** R
sumRcpp(c(42,43))
*/