※ 引述《ppp1987 (ppp)》之铭言:
: [问题类型]:
: 程式咨询(我想用R 做某件事情,但是我不知道要怎么用R 写出来)
: [软件熟悉度]:
: 入门(写过其他程式,只是对语法不熟悉)
: [问题叙述]:
: 资料形式
: ID Day X
: 1 1 0.5
: 1 3 0.1
: 1 4 0.3
: 1 7 0.5
: 1 9 0.5
: 1 11 0.2
: 1 14 0.5
: 2 1 0.1
: 2 2 0.4
: 2 5 0.8
: 2 9 0.7
: 2 11 0.1
: 2 13 0.2
: 现在我的资料每个ID(有100多个ID)每天有一笔观测值(但不一定每天有)
: 我现在想要算每个ID当天的前七天有观测值的平均
: 例如
: (ID=1,Day=9) Xhat=(0.5+0.3+0.1)/3
: (ID=2,Day=11) Xhat=(0.7+0.8)/2
: 现在已经用for循环跑出结果(但是很慢 2万多笔约30分钟)
: 想请教各位大大有没有比较快的方法
: 谢谢
我后来想到好像其实不用这样转换
直接算其实会是最快的,用rollapply反而很慢
library(data.table)
library(pipeR)
DT <- fread('ID Day X
1 1 0.5
1 3 0.1
1 4 0.3
1 7 0.5
1 9 0.5
1 11 0.2
1 14 0.5
2 1 0.1
2 2 0.4
2 5 0.8
2 9 0.7
2 11 0.1
2 13 0.2')
mvAvg_f <- function(val, day, span = 7) {
sapply(day, function(y){
tmp <- day - y
mean(val[which(tmp < 0 & tmp >= -span)])
})
}
DT[ , x_mean := mvAvg_f(X, Day), by = .(ID)]
# ID Day X x_mean
# 1: 1 1 0.5 NaN
# 2: 1 3 0.1 0.5000000
# 3: 1 4 0.3 0.3000000
# 4: 1 7 0.5 0.3000000
# 5: 1 9 0.5 0.3000000
# 6: 1 11 0.2 0.4333333
# 7: 1 14 0.5 0.4000000
# 8: 2 1 0.1 NaN
# 9: 2 2 0.4 0.1000000
# 10: 2 5 0.8 0.2500000
# 11: 2 9 0.7 0.6000000
# 12: 2 11 0.1 0.7500000
# 13: 2 13 0.2 0.4000000
benchmark: http://pastebin.com/Msib1dEh
Unit: milliseconds
expr min lq mean median uq max neval
me_1 4630.6863 4675.9866 4760.6722 4724.280 4873.9444 4902.0621 10
cywhale 3923.2040 3946.4697 3980.4697 3964.108 4030.4913 4063.7901 10
me_2 301.1468 307.6433 313.9698 309.936 314.3256 349.7625 10
me_1是用melt.data.table + rollapply的程式
me_2是这个方法,只需要0.3秒就可以算完了
根本不需要去用CJ 或是 melt.data.table 再配上rollapply去处理
直接算是最简单省事的,想复杂反而增加运算复杂度
QQ... 这也教了我一课
12/07 18:55补充:
可以直接简化成下面这样:
DT[ , x_mean := sapply(Day, function(s) mean(X[between(Day, s-7, s-1)])),
by = .(ID)]
但是会比较慢,请看benchmark... 我也不懂why
12/07 19:07补充:
between效果不彰,如果换成下面这样就会很快了
DT[ , x_mean := sapply(Day, function(s) mean(X[Day >= s-7 & Day < s])),
by = .(ID)]
请看benchmark,三万列只要0.15秒就搞定了