Подтвердить что ты не робот

Подсчет данных столбца в матрице с сбросами

Я собираю данные о том, насколько мои кошки попадают в матрицу:

m <- cbind(fluffy=c(1.1,1.2,1.3,1.4),misterCuddles=c(0.9,NA,1.1,1.0))
row.names(m) <- c("2013-01-01", "2013-01-02", "2013-01-03","2013-01-04")

Что дает мне это:

           fluffy misterCuddles
2013-01-01    1.1           0.9
2013-01-02    1.2            NA
2013-01-03    1.3           1.1
2013-01-04    1.4           1.0

В каждый день я хотел бы узнать, сколько дней подряд каждая кошка пошла на номер 2. Таким образом, результирующая матрица должна выглядеть так:

           fluffy misterCuddles
2013-01-01      1             1
2013-01-02      2             0
2013-01-03      3             1
2013-01-04      4             2

Есть ли способ сделать это эффективно? Функция cumsum делает что-то подобное, но это примитив, поэтому я не могу его модифицировать в соответствии с моими грязными, грязными потребностями.

Я мог бы запустить цикл for и сохранить счет так:

m.output <- matrix(nrow=nrow(m),ncol=ncol(m))
for (column in 1:ncol(m)) {
  sum <- 0
  for (row in 1:nrow(m)) {
    if (is.na(m[row,column])) sum <- 0
    else sum <- sum + 1

    m.output[row,column] <- sum
  }
}

Это самый эффективный способ сделать это? У меня много кошек, и я записал летние данные о корме. Могу ли я как-то параллелизировать это по столбцу?

4b9b3361

ОТВЕТЫ

Ответ 1

Это должно сработать. Обратите внимание: каждая из ваших кошек является независимым человеком, поэтому вы можете превратить свой фрейм данных в список и использовать mclapply, который использует параллельный подход.

count <- function(y,x){
  if(is.na(x)) return(0)
  return (y + 1)
}

oneCat = m[,1]

Reduce(count,oneCat,init=0,accumulate=TRUE)[-1]

EDIT: вот полный ответ

count <- function(x,y){
 if(is.na(y)) return(0)
 return (x + 1)
}

mclapply(as.data.frame(m),Reduce,f=count,init=0,accumulate=TRUE)

EDIT2: Основная проблема заключается в том, что я получаю дополнительные 0 в начале, поэтому...

result = mclapply(as.data.frame(m),Reduce,f=count,init=0,accumulate=TRUE)
finalResult = do.call('cbind',result)[-1,]
rownames(finalResult) = rownames(m)

выполняет задание.

Ответ 2

Все ответы здесь на самом деле слишком сложны (в том числе мои собственные, из ранее скопированных ниже). Семейство ответов Reduce просто маскирует цикл for в одном вызове функции. Мне нравятся Роланд и Ананда, но я думаю, что это слишком много происходит.

Таким образом, здесь простое векторное решение:

reset <- function(x) {
    s <- seq_along(x)
    s[!is.na(x)] <- 0
    seq_along(x) - cummax(s)
}

> apply(m, 2, reset)
     fluffy misterCuddles
[1,]      1             1
[2,]      2             0
[3,]      3             1
[4,]      4             2

Он также работает на примере Роланда:

m2 <- cbind(fluffy=c(NA,1.1,1.2,1.3,1.4,1.0,2),
           misterCuddles=c(NA,1.3,2,NA,NA,1.1,NA))

> apply(m2, 2, reset)
     fluffy misterCuddles
[1,]      0             0
[2,]      1             1
[3,]      2             2
[4,]      3             0
[5,]      4             0
[6,]      5             1
[7,]      6             0

Раньше: это не векторизация, но также работает:

pooprun <- function(x){
    z <- numeric(length=length(x))
    count <- 0
    for(i in 1:length(x)){
        if(is.na(x[i]))
            count <- 0
        else
            count <- + count + 1
        z[i] <- count
    }
    return(z)
}
apply(m, 2, pooprun)

> apply(m, 2, pooprun)
     fluffy misterCuddles
[1,]      1             1
[2,]      2             0
[3,]      3             1
[4,]      4             2

СЪЕМКА

Здесь я просто обернуваю все ответы в вызове функции (на основе их имени).

> library(microbenchmark)
> microbenchmark(alexis(), hadley(), thomas(), matthew(), thomasloop(), usobi(), ananda(), times=1000)
Unit: microseconds
         expr     min       lq   median       uq       max neval
     alexis()   1.540   4.6200   5.3890   6.1590   372.185  1000
     hadley()  87.755   92.758   94.298  96.6075  1767.012  1000
     thomas()  92.373  99.6860 102.7655 106.6140   315.223  1000
    matthew() 128.168 136.2505 139.7150 145.4880  5196.344  1000
 thomasloop() 133.556 141.6390 145.1030 150.4920 84131.427  1000
      usobi() 148.182 159.9210 164.7320 174.1620  5010.445  1000
     ananda() 720.507 742.4460 763.6140 801.3335  5858.733  1000

И вот результаты для данных примера Роланда:

> microbenchmark(alexis(), hadley(), thomas(), matthew(), thomasloop(), usobi(), ananda(), times=1000)
Unit: microseconds
         expr     min       lq   median       uq      max neval
     alexis()   2.310   5.3890   6.1590   6.9290   75.438  1000
     hadley()  75.053   78.902   80.058   83.136 1747.767  1000
     thomas()  90.834  97.3770 100.2640 104.3050  358.329  1000
    matthew() 139.715 149.7210 154.3405 161.2680 5084.728  1000
 thomasloop() 144.718 155.4950 159.7280 167.4260 5182.103  1000
      usobi() 177.048 188.5945 194.3680 210.9180 5360.306  1000
     ananda() 705.881 729.9370 753.4150 778.8175 8226.936  1000

Примечание. Решения Alexis и Hadley заняли довольно много времени, чтобы фактически определить как функции на моей машине, в то время как другие работают из коробки, но в противном случае Алексис станет явным победителем.

Ответ 3

Другой вариант, аналогичный @Usobi в том, что он использует Reduce, но с немного другим подходом:

apply(!is.na(m), 2, Reduce, f=function(x,y) if (y) x + y else y, accumulate=TRUE)
#      fluffy misterCuddles
# [1,]      1             1
# [2,]      2             0
# [3,]      3             1
# [4,]      4             2

Ответ 4

Я сохранил фрагмент из здесь, который почти точно подходит для такой проблемы:

countReset <- function(x) {
  x[!is.na(x)] <- 1
  y <- ave(x, rev(cumsum(rev(is.na(x)))), FUN=cumsum)
  y[is.na(y)] <- 0
  y
}
apply(m, 2, countReset)
#            fluffy misterCuddles
# 2013-01-01      1             1
# 2013-01-02      2             0
# 2013-01-03      3             1
# 2013-01-04      4             2

Ответ 5

Так как я нахожусь в периоде, когда пытаюсь привыкнуть к .Call, здесь появляется другая идея, которая, кажется, работает и, возможно, быстро. (Не верьте мне на слово, хотя мои навыки не заслуживают доверия!!):

library(inline)  #use "inline" package for convenience

f <- cfunction(sig = c(R_mat = "numeric", R_dims = "integer"), body = '
 R_len_t *dims = INTEGER(R_dims);
 R_len_t rows = dims[0], cols = dims[1];
 double *mat = REAL(R_mat);

 SEXP ans;
 PROTECT(ans = allocMatrix(INTSXP, rows, cols));
 R_len_t *pans = INTEGER(ans);

 for(int ic = 0; ic < cols; ic++)
  {
   pans[0 + ic*rows] = ISNA(mat[0 + ic*rows]) ? 0 : 1;

   for(int ir = 1; ir < rows; ir++)
    {
     if(ISNA(mat[ir + ic*rows]))
      {
       pans[ir + ic*rows] = 0;
      }else
      {
       if(!ISNA(mat[(ir - 1) + ic*rows]))
        {
         pans[ir + ic*rows] = pans[(ir - 1) + ic*rows] + 1;
        }else
        {
         pans[ir + ic*rows] = 1;
        }
      }
    }
  }

 UNPROTECT(1);

 return(ans);
')

f(m, dim(m))
#     [,1] [,2]
#[1,]    1    1
#[2,]    2    0
#[3,]    3    1
#[4,]    4    2
f(mm, dim(mm))   #I named Roland matrix, mm ; I felt that I had to pass this test!
#     [,1] [,2]
#[1,]    0    0
#[2,]    1    1
#[3,]    2    2
#[4,]    3    0
#[5,]    4    0
#[6,]    5    1
#[7,]    6    0

Ответ 6

Итак, решение этой проблемы состоит из двух частей:

  • Функция, которая принимает вектор за кошку и возвращает вектор, указывающий мне на каждую дату, сколько дней с момента последнего NA
  • Функция, которая принимает матрицу NxM и возвращает матрицу NxM, применяя функцию (1) к каждому столбцу

Для (2) я адаптировал это из ответа @Usobi:

daysSinceLastNA <- function(matrix, vectorFunction, cores=1) {
  listResult <- mclapply(as.data.frame(matrix), vectorFunction, mc.cores=cores)
  result <- do.call('cbind', listResult)
  rownames(result) <- rownames(matrix)
  result
}

Для (1) у меня есть два решения:

@решение ананда-махто:

daysSinceLastNA_1 <- function(vector) {
  vector[!is.na(vector)] <- 1
  result <- ave(vector, rev(cumsum(rev(is.na(vector)))), FUN=cumsum)
  result[is.na(result)] <- 0
  result
}

@Usobi:

daysSinceLastNA_2 <- function(vector) {
  reduction <- function(total, additional) ifelse(is.na(additional), 0, total + 1)
  Reduce(reduction, vector, init=0, accumulate=TRUE)[-1]
}

Затем я называю их следующим образом:

> system.time(result1 <- daysSinceLastNA (test, daysSinceLastNA_1 ))
   user  system elapsed 
   5.40    0.01    5.42 
> system.time(result2 <- daysSinceLastNA (test, daysSinceLastNA_2 ))
   user  system elapsed 
  58.02    0.00   58.03

В моем тестовом наборе данных, который представляет собой примерно матрицу 2500x2500, первый подход на порядок быстрее.

Если я запускаю linux с 64 ядрами, решение (1) запускается через 2 секунды, а решение (2) запускается через 6 секунд.

Ответ 7

Для такого рода проблем, который легко решается с помощью цикла for, я нахожу Rcpp очень естественным ответом.

library(Rcpp)

cppFunction("NumericVector cumsum2(NumericVector x) {
  int n = x.length();
  NumericVector out(x);

  for(int i = 0; i < n; ++i) {
    if (NumericVector::is_na(x[i]) || i == 0) {
      x[i] = 0;
    } else {
      x[i] = x[i - 1] + 1;
    }
  }

  return out;
}")

Код требует немного больше бухгалтерского учета, чем эквивалентный R-код, но основная часть функции - очень простой цикл.

Затем вы можете применить в R как любую другую векторную функцию:

m2 <- cbind(
  fluffy=c(NA,1.1,1.2,1.3,1.4,1.0,2),
  misterCuddles=c(NA,1.3,2,NA,NA,1.1,NA)
)

apply(m2, 2, cumsum2)

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