1

我有一个数据文件,其中包含每日 xy 位置和表示该位置是否为异常值的逻辑向量。这是一些(我知道创建得很差)示例数据:

x=seq(3,10,length.out=30)
y=seq(42,45,length.out=30)
outlier=c(F,F,F,F,F,F,F,F,T,T,T,F,F,F,F,F,F,F,F,F,F,T,F,T,F,F,F,F,F,F)
data=cbind(x,y,outlier)
> data
             x           y outlier
 [1,]  3.000000000 42.00000000       0
 [2,]  3.241379310 42.10344828       0
 [3,]  3.482758621 42.20689655       0
 [4,]  3.724137931 42.31034483       0
 [5,]  3.965517241 42.41379310       0
 [6,]  4.206896552 42.51724138       0
 [7,]  4.448275862 42.62068966       0
 [8,]  4.689655172 42.72413793       0
 [9,]  4.931034483 42.82758621       1
[10,]  5.172413793 42.93103448       1
[11,]  5.413793103 43.03448276       1
[12,]  5.655172414 43.13793103       0
[13,]  5.896551724 43.24137931       0
[14,]  6.137931034 43.34482759       0
[15,]  6.379310345 43.44827586       0
[16,]  6.620689655 43.55172414       0
[17,]  6.862068966 43.65517241       0
[18,]  7.103448276 43.75862069       0
[19,]  7.344827586 43.86206897       0
[20,]  7.586206897 43.96551724       0
[21,]  7.827586207 44.06896552       0
[22,]  8.068965517 44.17241379       1
[23,]  8.310344828 44.27586207       0
[24,]  8.551724138 44.37931034       1
[25,]  8.793103448 44.48275862       0
[26,]  9.034482759 44.58620690       0
[27,]  9.275862069 44.68965517       0
[28,]  9.517241379 44.79310345       0
[29,]  9.758620690 44.89655172       0
[30,] 10.000000000 45.00000000       0

我需要的是对 x 和 y 列取一个不重叠的 6 天平均值。这很容易使用rollapply(). 但是,我不希望outlier=1将值包含在 6 天平均值中;我也不希望 6 天窗口通过删除所有行 where 来“跨越”留下的空白outlier=T。相反,我想对“非重叠规则”做一个例外。

我认为最好使用上面的示例数据来解释这一点:第一个值应该是 1:6 行的平均值,而不是第二个值是 7:12 行(包括outlier=1值)或 c(7: 8,12:15)(跳过outlier=1值)我希望它与第一个窗口重叠并取 3:8 行的平均值。

因此,对于上述长度为 30 的样本数据,最终结果的长度应为 5,显示第 1:6、3:8、12:17、16:21 和 25:30 行的平均值(理想情况下,所有值都来自重叠窗口应该这样标记;即值 1:4 重叠,而最终值是唯一的)

4

1 回答 1

2

这是一个函数,它将为您提供所需平均值的端点索引:

findIndices<-function(outlier,window=6){
  r<-rle(outlier)
  rends<-cumsum(r$lengths)
  segs<-cbind(rends-r$lengths+1,rends)
  segs<-segs[with(r,lengths>=window & values==0),]

  indices<-unlist(apply(segs,1,function(x) seq(x[1]+window-1,x[2],by=window)))
  sort(unique(c(indices,segs[,2])))     
}

findIndices(data[,3])
## [1]  6  8 17 21 30

然后,您可以像这样获得所需的平均值:

id<-findIndices(data[,3])
require(zoo)
cbind(index=id,rollmean(data[,1:2],6)[id-5,])
##     index        x        y
## [1,]     6 3.603448 42.25862
## [2,]     8 4.086207 42.46552
## [3,]    17 6.258621 43.39655
## [4,]    21 7.224138 43.81034
## [5,]    30 9.396552 44.74138

您可以将它们放在一个函数中,如下所示:

maWithOutliers<-function(x,outlier,window){
  id<-findIndices(outlier,window)
  cbind(index=id,rollmean(x,window)[id-window+1,])
}

> maWithOutliers(data[,1:2],data[,3],6)
     index        x        y
[1,]     6 3.603448 42.25862
[2,]     8 4.086207 42.46552
[3,]    17 6.258621 43.39655
[4,]    21 7.224138 43.81034
[5,]    30 9.396552 44.74138
> maWithOutliers(data[,1:2],data[,3],4)
     index        x        y
[1,]     4 3.362069 42.15517
[2,]     8 4.327586 42.56897
[3,]    15 6.017241 43.29310
[4,]    19 6.982759 43.70690
[5,]    21 7.465517 43.91379
[6,]    28 9.155172 44.63793
[7,]    30 9.637931 44.84483
> 
于 2013-10-05T23:00:49.937 回答