[發(fā)明專利]一種基于樣本熵和貝葉斯的時間序列突變檢測方法

有效

• 申請?zhí)枺?/td>: 202010511526.1
• 申請日: 2020-06-08
• 公開（公告）號: CN111667009B
• 公開（公告）日: 2023-04-04
• 發(fā)明（設(shè)計）人: 孫東永;張洪波;徐明珠;孔令魁;李楊津;李振欣;王琪
• 申請（專利權(quán)）人: 長安大學(xué)
• 主分類號: G06N7/01
• 分類號: G06N7/01;G06F17/18;G06F16/2458
• 代理公司: 北京東方盛凡知識產(chǎn)權(quán)代理有限公司 11562
• 代理人: 賈耀淇

1.一種基于樣本熵和貝葉斯的時間序列突變檢測方法，其特征在于，包括如下步驟：

S1.將長度為n的徑流原始序列{x₁,x₂,…,x_n}進(jìn)行相空間重構(gòu)，構(gòu)造如下一組m維向量：

X(i)＝[x(i),x(i+1),…,x(i+m-1)]，i＝1,2,…,n-m+1

S2.計算向量X(i)和X(j)之間的歐氏距離d[X(i),X(j)]：

d[X(i),X(j)]＝max[|X(i+k)-X(j+k)|],k＝0,1,2,…,m-1

S3.給定允許偏差r，計算每一個i所對應(yīng)的歐式距離d[X(i),X(j)]小于r的數(shù)目N_r與n-m+1的比值

S4.對于數(shù)值取i的所有值的平均值，得到C^m(r)：

S5.將維數(shù)m增加1并重復(fù)步驟S1到S4以獲得C^m+1(r)；

S6.用數(shù)值m、r和n計算徑流時間序列的SampEn(m,r,n)值：

SampEn(m,r,n)＝-ln[C^m+1(r)/C^m(r)]

S7.根據(jù)待分析的徑流時間序列{x₁,x₂,…,x_n}，選擇數(shù)據(jù)滑動窗口長度h以及滑動步長L，確定樣本熵SampEn的維數(shù)m和容許偏差r；

S8.從徑流時間序列{x₁,x₂,…,x_n}的第i個數(shù)據(jù)開始以滑動窗口長度h選取徑流子序列，其中i＝1,2,…,n-h+1；

S9.通過樣本熵SampEn計算子序列的熵值；

S10.保持?jǐn)?shù)據(jù)滑動窗口長度h不變，以滑動步長L逐步移動窗口，重復(fù)S8～S9步驟，直至原流時間序列{x₁,x₂,…,x_n}結(jié)束；

S11.通過步驟S7～S10，得到一個長度為int[(n-h+1)/L]的樣本熵SampEn序列；

S12.繪制樣本熵SampEn的值隨時間的變化圖，初步確定序列的突變點；

S13.假設(shè)徑流時間序列{x₁,x₂,…,x_n}服從正態(tài)分布，序列k點處發(fā)生突變，計k點前的序列為X^k＝{x₁,x₂,…,x_k}，k點后的序列為X^k+1＝{x_k+1,x_k+2,…,x_n}；當(dāng)X^k和X^k+1的概率分布的統(tǒng)計參數(shù)不相同時，X^k的分布密度函數(shù)表示如下：

X^k+1的分布密度函數(shù)表示如下：

假設(shè)X^k和X^k+1的概率分布的分布函數(shù)均服從正態(tài)分布：

上式中，σ²由徑流時間序列{x₁,x₂,…,x_n}估算；

后驗分布如下：

2.根據(jù)權(quán)利要求1所述的基于樣本熵和貝葉斯的時間序列突變檢測方法，其特征在于，分布參數(shù)μ_a、μ_b的后驗密度函數(shù)的推導(dǎo)計算如下：

Step1.給定μ_a和μ_b，X^k和X^k+1的聯(lián)合分布函數(shù)導(dǎo)出如下：

Step2.根據(jù)Bayes準(zhǔn)則推出變異點k的后驗分布密度函數(shù)：

其中，p(j)是改變點k的位置的先驗分布；

假定p(j)為均勻分布：p(j)＝1/n,j＝1,2,…,n；后驗密度函數(shù)為：