纠正时间序列中的连续错误
Correcting consecutive errors in time series
我有非常大的连续数据集(>1M 行),由于传感器故障或其他外部因素,经常出现 "breaks" 或 "jumps"。这些中断对应于添加或删除的恒定值,并且仅持续有限的时间。我正在尝试将这些序列与其余数据重新对齐。
par(mfrow=c(2,1))
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
#plotting ideal vs real data
plot(dfe$date, dfe$valuer, main="real data", ylim=c(8,32))
plot(dfe$date, dfe$valueideal, main="ideal data", ylim=c(8,32))
所以我的数据看起来像 "real data",我希望他们将它们重新调整为像 "ideal data"。
到目前为止,我已经做了一个 for
循环,除了每个工件的第一个数据点外,它大部分都有效,并且它对其余数据有轻微影响。我不确定为什么或如何解决它:
#trying to solve it with a loop
dfe$valuel<-dfe$valuer
for (i in seq(2,length(dfe$valuel)-1,1)){
future<-diff(c(dfe$valuel[i],dfe$valuel[i+1]))
past<-diff(c(dfe$valuel[i-1],dfe$valuel[i]))
if (abs(future)>2*abs(past)){
dfe$valuel[i:length(dfe$valuel)]<-dfe$valuel[i:length(dfe$valuel)]-future
}
}
plot(dfe$date, dfe$valuel, main="loop corrected data", ylim=c(8,32))
我也担心在我非常大的数据集上使用这种方法,我不确定这需要多长时间。所以我也尝试过使用这种 方法,但效果不佳,可能是因为很难选择相关的 delta_max
值:
#trying to solve it with a vectorised function
remove_artifacts <- function(weights, delta_max) {
# calculate deltas, and set first delta to zero
dw <- c(0, diff(x))
# create vector of zeros and abs(observations) > delta_max
# dw * (logical vector) results in either:
# dw * 0 (if FALSE)
# dw * 1 (if TRUE)
dm <- dw * (abs(dw) > delta_max)
# subtract the cumulative sum of observations > delta_max
return(weights - cumsum(dm))
}
dfe$valuedm<-remove_artifacts(dfe$valuer, 10)
plot(dfe$date, dfe$valuedm, main="remove artifacts function", ylim=c(8,32))
所以我的问题是,我怎样才能有效地纠正这些连续的数据中断?
这是一个不完美的解决方案。首先,我用你的代码设置问题。
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
接下来,我使用 strucchange
包中的 breakpoints
来查找时间序列中的断点。
# Find breakpoints
bp <- strucchange::breakpoints(valuer ~ date, data = dfe)
# Get breakpoints plus start & end of time series
int <- c(1, bp$breakpoints + 1, nrow(dfe))
在这里,我根据断点为数据集添加标签。我绘制了包括颜色在内的图表,以查看我们的表现如何。 (一个散兵游勇,还算不错。)
# Create labels
dfe$label <- cut(1:nrow(dfe),
breaks = int,
include.lowest = TRUE,
right = FALSE)
# Plot "real" data coloured by label
with(dfe, plot(date, valuer, col = label, main="real data", ylim=c(8,32)))
然后我切换到 data.table,因为那是我的果酱。
# Load library
library(data.table)
# Convert to data.table
setDT(dfe)
我按 label
分组,然后使用每个区间的平均值进行校正。
# Offset by mean
dfe[, corrected := valuer - mean(valuer), by = label]
# Plot again
with(dfe, plot(date, corrected, col = label, main = "Corrected data", ylim = c(-10, 10)))
由 reprex package (v0.2.1.9000)
创建于 2019-12-02
掉队者稍微偏离了这个间隔,但更正后的解决方案并不糟糕。
TL;DR
# Find breakpoints
bp <- strucchange::breakpoints(valuer ~ date, data = dfe)$breakpoints
# Add start & end points
int <- c(1, bp + 1, nrow(dfe))
# Tag intervals
dfe$label <- cut(1:nrow(dfe),
breaks = int,
include.lowest = TRUE,
right = FALSE)
# Correct by subtracting mean from each interval
data.table::setDT(dfe)[, corrected := valuer - mean(valuer), by = label]
这是另一个解决方案,它使用 data.table 时速度很快,因此可以就地修改。首先,我设置问题。
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
接下来,我加载 data.table 并将数据帧转换为数据 table。
# Load data.table
library(data.table)
# Convert data frame into data.table
setDT(dfe)
我使用向量化方法计算连续值的差异,而不是问题中的循环。
# Calculate changes
dfe[, delta := c(abs(diff(valuer)), 0)]
这些差异用于将时间序列分成间隔:
# Labels intervals
dfe[, int := cut(.I,
breaks = c(0, which(delta > 3 * sd(delta)/mean(delta)), nrow(dfe)),
include.lowest = TRUE)]
我将所有间隔都集中在零上。
# Mean of zero
dfe[, value_new := valuer - mean(valuer), by = int]
然后,我添加一个偏移量作为第一组的平均值。
# Correct offset
dfe[, value_new := value_new + dfe[, mean(valuer), by = int][, first(V1)]]
最后,我绘制了结果。
# Plot result
with(dfe, plot(date, value_new, main="real data", ylim=c(8,32)))
由 reprex package (v0.3.0)
于 2019-12-11 创建
我有非常大的连续数据集(>1M 行),由于传感器故障或其他外部因素,经常出现 "breaks" 或 "jumps"。这些中断对应于添加或删除的恒定值,并且仅持续有限的时间。我正在尝试将这些序列与其余数据重新对齐。
par(mfrow=c(2,1))
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
#plotting ideal vs real data
plot(dfe$date, dfe$valuer, main="real data", ylim=c(8,32))
plot(dfe$date, dfe$valueideal, main="ideal data", ylim=c(8,32))
所以我的数据看起来像 "real data",我希望他们将它们重新调整为像 "ideal data"。
到目前为止,我已经做了一个 for
循环,除了每个工件的第一个数据点外,它大部分都有效,并且它对其余数据有轻微影响。我不确定为什么或如何解决它:
#trying to solve it with a loop
dfe$valuel<-dfe$valuer
for (i in seq(2,length(dfe$valuel)-1,1)){
future<-diff(c(dfe$valuel[i],dfe$valuel[i+1]))
past<-diff(c(dfe$valuel[i-1],dfe$valuel[i]))
if (abs(future)>2*abs(past)){
dfe$valuel[i:length(dfe$valuel)]<-dfe$valuel[i:length(dfe$valuel)]-future
}
}
plot(dfe$date, dfe$valuel, main="loop corrected data", ylim=c(8,32))
我也担心在我非常大的数据集上使用这种方法,我不确定这需要多长时间。所以我也尝试过使用这种 delta_max
值:
#trying to solve it with a vectorised function
remove_artifacts <- function(weights, delta_max) {
# calculate deltas, and set first delta to zero
dw <- c(0, diff(x))
# create vector of zeros and abs(observations) > delta_max
# dw * (logical vector) results in either:
# dw * 0 (if FALSE)
# dw * 1 (if TRUE)
dm <- dw * (abs(dw) > delta_max)
# subtract the cumulative sum of observations > delta_max
return(weights - cumsum(dm))
}
dfe$valuedm<-remove_artifacts(dfe$valuer, 10)
plot(dfe$date, dfe$valuedm, main="remove artifacts function", ylim=c(8,32))
所以我的问题是,我怎样才能有效地纠正这些连续的数据中断?
这是一个不完美的解决方案。首先,我用你的代码设置问题。
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
接下来,我使用 strucchange
包中的 breakpoints
来查找时间序列中的断点。
# Find breakpoints
bp <- strucchange::breakpoints(valuer ~ date, data = dfe)
# Get breakpoints plus start & end of time series
int <- c(1, bp$breakpoints + 1, nrow(dfe))
在这里,我根据断点为数据集添加标签。我绘制了包括颜色在内的图表,以查看我们的表现如何。 (一个散兵游勇,还算不错。)
# Create labels
dfe$label <- cut(1:nrow(dfe),
breaks = int,
include.lowest = TRUE,
right = FALSE)
# Plot "real" data coloured by label
with(dfe, plot(date, valuer, col = label, main="real data", ylim=c(8,32)))
然后我切换到 data.table,因为那是我的果酱。
# Load library
library(data.table)
# Convert to data.table
setDT(dfe)
我按 label
分组,然后使用每个区间的平均值进行校正。
# Offset by mean
dfe[, corrected := valuer - mean(valuer), by = label]
# Plot again
with(dfe, plot(date, corrected, col = label, main = "Corrected data", ylim = c(-10, 10)))
由 reprex package (v0.2.1.9000)
创建于 2019-12-02掉队者稍微偏离了这个间隔,但更正后的解决方案并不糟糕。
TL;DR
# Find breakpoints
bp <- strucchange::breakpoints(valuer ~ date, data = dfe)$breakpoints
# Add start & end points
int <- c(1, bp + 1, nrow(dfe))
# Tag intervals
dfe$label <- cut(1:nrow(dfe),
breaks = int,
include.lowest = TRUE,
right = FALSE)
# Correct by subtracting mean from each interval
data.table::setDT(dfe)[, corrected := valuer - mean(valuer), by = label]
这是另一个解决方案,它使用 data.table 时速度很快,因此可以就地修改。首先,我设置问题。
#simulating perfect dataset
dfe<-data.frame(
date=seq(as.Date('2015-07-12'),as.Date('2015-07-12')+49, by = '1 day'),
valueideal=round(sin(seq(1,50,1))+20)
)
#introducing artifacts
dfe$valuer<-dfe$valueideal
dfe$valuer[10:20]<-dfe$valueideal[10:20]+10
dfe$valuer[30:35]<-dfe$valueideal[30:35]-10
接下来,我加载 data.table 并将数据帧转换为数据 table。
# Load data.table
library(data.table)
# Convert data frame into data.table
setDT(dfe)
我使用向量化方法计算连续值的差异,而不是问题中的循环。
# Calculate changes
dfe[, delta := c(abs(diff(valuer)), 0)]
这些差异用于将时间序列分成间隔:
# Labels intervals
dfe[, int := cut(.I,
breaks = c(0, which(delta > 3 * sd(delta)/mean(delta)), nrow(dfe)),
include.lowest = TRUE)]
我将所有间隔都集中在零上。
# Mean of zero
dfe[, value_new := valuer - mean(valuer), by = int]
然后,我添加一个偏移量作为第一组的平均值。
# Correct offset
dfe[, value_new := value_new + dfe[, mean(valuer), by = int][, first(V1)]]
最后,我绘制了结果。
# Plot result
with(dfe, plot(date, value_new, main="real data", ylim=c(8,32)))
由 reprex package (v0.3.0)
于 2019-12-11 创建