机器学习:异常检测算法Seasonal Hybrid ESD及R语言实现

时间:2022-04-11 17:12:06

  Twritters的异常检测算法(Anomaly Detection)做的比较好,Seasonal Hybrid ESD算法是先用STL把序列分解,考察残差项。假定这一项符合正态分布,然后就可以用Generalized ESD提取离群点。

  目标是检测出时间序列数据集的异常点,如图所示,蓝色线是时间序列数据集,红色是圈是异常点。

  机器学习:异常检测算法Seasonal Hybrid ESD及R语言实现

  

  R语言实现如下,一些依赖包需要install.packages("")或者手动在cran社区下载(注意依赖包的下载)。本人github下载源码。

  1 主函数是,包含了主要逻辑,加载数据集,IMSHESD算法检测异常点和画出数据集和异常点。IMSHESD算法是主要功能,下面详细介绍。

library(zoo)
path_data="E:/Develop/Rstudio/IMS-H-ESDS/3151.csv"
path_sear="E:/Develop/Rstudio/IMS-H-ESDS/IS-H-ESD.R" source(path_sear)
data <-read.table(path_data,sep=",",skip=)
print(data)
test_data=IMSHESD(data)
plot(x=data[[]], y=data[[]],xlab="time",ylab="value",col="blue",type="l")
lines(x=test_data[[]], y=test_data[[]],col="red",type="p")

  2 IMSHESD算法是主要逻辑如下,通过Fourier转换自动求得时间序列的季节周期peri(必须满足数据集的长度>2*peri+1才可以应用时间序列分析),按照季节周期对数据集做划分,然后应用anmodetection异常检测算法探测异常。

IMSHESD<-function(data,group_peri=)
{
path_Fourier="E:/Develop/Rstudio/IMS-H-ESDS/Fourier.R"
path_data_group="E:/Develop/Rstudio/IMS-H-ESDS/data_group.R"
path_anmo_detection="E:/Develop/Rstudio/IMS-H-ESDS/anmo_detection.R"
source(path_Fourier)
source(path_data_group)
source(path_anmo_detection)
#data <-data_group(data,mode="median",group_period=)
peri=Fourier_trans(data)
print(peri)
if(ncol(data)!=)
{
print("The col of data must be two!")
stop()
}
if((*peri+)<length(data[[]])){
data_sep_length=ceiling(*peri+)
print(data_sep_length)
all_data <- vector(mode="list", length=ceiling(length(data[[]])/(data_sep_length)))
for(j in seq(,length(data[[]]), by=data_sep_length)){
start_data <- data[[]][j]
end_data <- data[[]][min(j + data_sep_length, length(data[[]]))]
if(j+data_sep_length<length(data[[]])){
all_data[[ceiling(j/(data_sep_length))]] <- subset(data, data[[]] >= start_data & data[[]] < end_data)
}else{
all_data[[ceiling(j/(data_sep_length))]] <- subset(data,data[[]] >= data[[]][length(data[[]])-data_sep_length] & data[[]] < end_data)
}
}
res=c()
for(i in :length(all_data))
{
res_temp=anmodetection(all_data[[i]],anoms_per=0.1,period=peri,alpha=0.05)
res=c(res,res_temp)
}
data_plot=rep(c(),length(res))
for(i in :length(res))
{
data_plot[i]=data[[]][which(data[[]]==res[i])]
}
anmo_point<-data.frame(res,data_plot)
}else{
print("This is not a seasonal time series")
stop()
}
}

  3 Fourier转换自动求得时间序列的季节周期peri。

Fourier_trans<-function(data)
{
install.packages("TSA")
library(TSA)
p=periodogram(data[])
dd=data.frame(freq=p$freq,spec=p$spec)
order=dd[order(-dd$spec),]
top2=head(order,)
time=min(1.00/top2$f)
}

  4 异常检测主要逻辑anmodetection函数,需要规定异常点的上限10%,STL分解数据集:周期+趋势+随机噪声=原始时间序列(分解方法有Twitters的Decompose和STL),残差项根据正态分布(方差未知使用学生t分布),提取离散点。假设要检测k个离群点,就对数据重复使用k次ESD检验,如果发现离群点就从数据里剔出去,然后在剩下的数据上重新检测(Generalized ESD)。

anmodetection<-function(data,anoms_per=0.10,period=,alpha=0.05,mode="addi")
{
num <- length(data[[]])
#cat("num",num)
num_anmo=trunc(anoms_per*length(data[[]]))
R_idx=rep(seq(),length(data[[]]))
if(ncol(data)!=)
{
print("The col of data must be two!")
stop()
}
data_decompose <- stl(ts(data[[2L]], frequency =period),"periodic",robust = TRUE)
seasonal_data <- data_decompose$time.series[,]
#cat("seasonal_data",seasonal_data)
trend_data <- data_decompose$time.series[,]
random_decomp <- data_decompose$time.series[,]
data<-data.frame(times=data[[]],count=(data[[]]-seasonal_data-median(data[[]])))
func_med <- match.fun(median)
func_mad <- match.fun(mad)
numb_anom=
for(n in :num_anmo)
{
data_norm<-abs(data[[]]-func_med(data[[]]))
data_mad<-func_mad(data[[]])
data_res<-data_norm/data_mad
R_res<-max(data_res)
max_temp_idx <- which(data_res == R_res)[]
R_idx_out=data[[]][max_temp_idx]
R_idx[n] <- R_idx_out
data <- data[-which(data[[]] == R_idx[n]), ]
p <- - alpha/(*(num-n+))
t <- qt(p,(num-n-))
thres <- t*(num-n) / sqrt((num-n-+t**)*(num-n+))
if(R_res>thres)
{
numb_anom <- n
}
}
if(numb_anom>)
{
R_idx <- R_idx[:numb_anom]
}
else
{
R_idx = NULL
}
return(R_idx)
}

  实验结果:红色圈标出了异常点,异常点是10%的比例。

机器学习:异常检测算法Seasonal Hybrid ESD及R语言实现

 参考

  Twitters异常检测方法,https://anomaly.io/blog/