4

馬が参加した一連のレースの時間 (日) と着順 (pos) を指定して、馬の移動最新加重平均着順を計算したいと思います。このような統計はハンディキャップに役立ちます。

現在、「ループ内ループ」アプローチを使用しています。この問題に対するより高速でエレガントな R 言語のアプローチはありますか?

#
# Test data
#

day <- c(0, 6, 10, 17, 21, 26, 29, 31, 34, 38, 41, 47, 48, 51, 61)
pos <- c(3, 5, 6, 1, 1, 3, 4, 1, 2, 2, 2, 6, 4, 5, 6)
testdata <- data.frame(id = 1, day = day, pos = pos, wt.pos = NA)

#
# No weight is given to observations earlier than cutoff
#

cutoff <- 30

#
# Rolling recency-weighted mean (wt.pos)
#

for(i in 2:nrow(testdata)) {
  wt <- numeric(i-1)
  for(j in 1:(i-1))
    wt[j] <- max(0, cutoff - day[i] + day[j] + 1)
    if (sum(wt) > 0)
      testdata$wt.pos[i] <- sum(pos[1:j] * wt) / sum(wt)
}

> testdata

   id day pos   wt.pos
1   1   0   3       NA
2   1   6   5 3.000000
3   1  10   6 4.125000
4   1  17   1 4.931034
5   1  21   1 3.520548
6   1  26   3 2.632911
7   1  29   4 2.652174
8   1  31   1 2.954128
9   1  34   2 2.436975
10  1  38   2 2.226891
11  1  41   2 2.119048
12  1  47   6 2.137615
13  1  48   4 3.030534
14  1  51   5 3.303704
15  1  61   6 4.075000
4

2 に答える 2

0

私は行きます

# Calculate `wt` for all values of `i` in one go
wt <- lapply(2:nrow(testdata), function(i)
    pmax(0, cutoff - day[i] + day[1:(i-1)] + 1))

# Fill in the column
testdata$wt.pos[-1] <- mapply(
    function(i, w) if(sum(w) > 0) sum(pos[1:i]*w)/sum(w) else NA,
    1:(nrow(testdata)-1), wt)

maxのすべての値に対して2番目の引数をj同時に計算することにより、計算をベクトル化したことに注意してください。これにより、速度が何桁も向上します。

外側のループとifケースをベクトル化する簡単な方法は見つかりませんでしたが(やり過ぎのように見えるCで書き直すことを除けば)、、、lapplyおよび同様のものはループmapplyよりも高速です。for

于 2012-09-05T16:14:50.760 に答える
0

このバージョンでは、1 つまたは複数の変数 (たとえば、最終順位、速度評価など) と 1 つまたは複数の被験者 (馬) について移動最新性加重平均を計算する方法を示します。

library(plyr)

day <- c(0, 6, 10, 17, 21, 26, 29, 31, 34, 38, 41, 47, 48, 51, 61)
pos <- c(3, 5, 6, 1, 1, 3, 4, 1, 2, 2, 2, 6, 4, 5, 6)
dis <- 100 + 0.5 * (pos - 1)
testdata1 <- data.frame(id = 1, day = day, pos = pos, dis = dis)
day <- c(0, 4, 7, 14, 22, 23, 31, 38, 42, 47, 52, 59, 68, 69, 79)
pos <- c(1, 3, 2, 6, 4, 5, 2, 1, 4, 5, 2, 1, 5, 5, 2)
dis <- 100 + 0.5 * (pos - 1)
testdata2 <- data.frame(id = 2, day = day, pos = pos, dis = dis)
testdata <- rbind(testdata1, testdata2)

# Moving recency-weighted mean
rollmean <- function(day, obs, cutoff = 90) {
  obs <- as.matrix(obs)
  wt <- lapply(2:nrow(obs), function(i)
    pmax(0, cutoff - day[i] + day[1:(i-1)] + 1))
  wt.obs <- lapply(1:(nrow(obs)-1), FUN =
    function(i)
      if(sum(wt[[i]]) > 0) {
        apply(obs[1:i, , drop = F] * wt[[i]], 2, sum) / sum(wt[[i]])
      } else {
        rep(NA, ncol(obs))
      }
  )
  answer <- rbind(rep(NA, ncol(obs)), do.call(rbind, wt.obs))
  if (!is.null(dimnames(answer)))
    dimnames(answer)[[2]] <- paste("wt", dimnames(answer)[[2]], sep = ".")
  return(answer)
}

x <- dlply(testdata, .(id), .fun =
  function(DF) rollmean(DF$day, DF[, c("pos", "dis"), drop = F])
)
y <- do.call(rbind, x)
于 2012-09-08T22:57:11.413 に答える