Patent ID: 11927949
Assignee: ZHEJIANG UNIVERSITY
Field: Control (Instruments)
Classification: CPC G  H | IPC G

Claim 2:
3. The method for anomaly classification of an ICS communication network according to claim 2, wherein step 2) comprises:
2.1) based on a selected traffic aggregation scale and a short-cycle analysis scale, generating a SARIMA(p,d,q)×(P,D,Q)s time series by using a SARIMA(p,d,q)×(P,D,Q)s sequence definition method, comprising:
separately performing d-order difference calculation and D-order seasonal difference calculation on an auto regressive moving average ARMA(p,q) model, to obtain a SARIMA(p,d,q)×(P,D,Q)s model, and using AR(p) and MA(q) models to obtain the ARMA(p,q) model, wherein
the ARMA(p,q) model comprises:

Xt=ϕ1Xt−1+ϕ2Xt−2+ . . . +ϕpXt−p+εt−θ1εt−1− . . . −θqεt−q

in the foregoing formula, Xt represents a short-cycle stationary time series after averaging, with a relatively short length; ϕp represents a coefficient of an auto regressive term AR; θq represents a coefficient of a moving average term MA; εt represents a random error term; p represents an order number of AR; and q represents an order number of MA;
substituting a delay operator B into the ARIMA(p, d, q) model, to obtain an AR coefficient polynomial Φ(B)=1−ϕ1B− . . . −ϕp(B)p and an MA coefficient polynomial Θ(B)=1−θ1B− . . . −θq(B)q, wherein BXt=Xt−1;
introducing a difference operator Δd=(1−B)d, wherein the ARIMA(p,d,q) model comprises:

Φ(B)ΔdXt=Θ(B)εt

performing seasonal difference calculation on the ARIMA model, to obtain the SARIMA model, wherein the SARIMA model comprises:

Φp(B)Φp(Bs)ΔdΔsDXt=Θq(B)ΘQ(Bs)εt

wherein εt represents a white noise sequence, d represents an order number of a trend difference, D represents an order number of a seasonal difference compensated based on cycle s, Bs represents an s-order delay operator, ΔsD represents a seasonal difference operator, BsXt=Xt−s , ΔsD=1−Bs, ΦP(Bs) is a Q-order polynomial of Bs, and ΦP(Bs) is a P-order polynomial of Bs;
2.2) using Bayesian Information Criterion (BIC) to supervise, analyze, and determine p, d, q, P, D, and Q orders of the SARIMA(p, d, q)×(P, D, Q)s model;
2.3) using a least square method to estimate a p-order coefficient ϕk(k=1,2, . . . , p), a q-order coefficient θk(k=1,2, . . . , q), a seasonal P-order coefficient ϕks(k=s,2·s, . . . , P·s), and a seasonal Q-order coefficient θks(k=s,2·s, . . . , Q·s) of SARIMA(p, d, q)×(P, D, Q)s;
2.4) using the SARIMA(p, d, q)×(P, D, Q)s model under the optimal BIC to perform fitting analysis on the original sequence, and performing residual testing; and if a residual is white noise, performing inverse filtering on a fitted sequence to obtain a fitted value or a forecast value of the original sequence; or if the residual is not white noise, using the BIC to determine orders of the ARMA(p, q) model;
2.5) obtaining a mathematical expression of the short-cycle SARIMA(p, d, q)×(P, D, Q)s model;
2.6) collecting real-time traffic data from the ICS industrial switch, and generating a time series based on a specified sampling frequency γsamp and an aggregation scale, with a short cycle as an iterative cycle;
2.7) performing training based on the collected real-time traffic data, performing SARIMA(p, d, q)×(P, D, Q)s model training and adaptation in a short cycle, and then outputting an optimal model and parameters adapted to the optimal model, wherein a model in the i-th short cycle comprises:

{circumflex over (X)}Tfore(i)=fSARIMA(XTtrai(i), Tfore, s, ‘BIC’)

wherein fSARIMA( ) is a functional expression of the SARIMA(p, d, q)×(P, D, Q)s model, XTtrai(i) is a short-cycle training set for the i-th iteration, Tfore is the number of forecast sequences in the short cycle, s is a periodic parameter, ‘BIC’ is a criterion for selecting optimal (p, d, q, P, D, Q) parameters based on econometrics, and {circumflex over (X)}Tfore(i) is a time series forecast at the i-th iteration of the SARIMA model;
calculating a forecasting average {circumflex over (μ)}(i) of the i-th iteration:, μ
    ^
   
   
    (
    i
    )
   
  
  =
  
   
    1
    
     T
     fore
    
   
   ·
   
    
     ∑
     
      k
      =
      1
     
     
      T
      fore
     
    
    
     
      X
      ^
     
     k
     
      (
      i
      )
     
    
   
  
 

2.8) running the short-cycle SARIMA(p, d, q)×(P, D, Q)s model in a distributed manner, and performing real-time rolling modeling based on ICS traffic data collected in real time in the i-th short cycle; using real-time traffic data as a validation set, and comparing the validation set with upper and lower bounds of a confidence interval-based traffic threshold generated in the forecasting process, wherein upper and lower thresholds for the i-th short cycle comprise:, U
   
    T
    fore
   
   
    (
    i
    )
   
  
  =
  
   
    
     X
     ^
    
    
     T
     fore
    
    
     (
     i
     )
    
   
   +
   
    
     z
     
      (
      
       1
       -
       
        α
        
         P
         .
         I
        
       
      
      )
     
    
    ⁢
    
     
      
       1
       
        T
        fore
       
      
      ·
      
       
        ∑
        
         k
         =
         1
        
        
         T
         fore
        
       
       
        
         [
         
          
           
            X
            ^
           
           k
           
            (
            i
            )
           
          
          -
          
           
            μ
            ^
           
           
            (
            i
            )
           
          
         
         ]
        
        2
       
      
     
    
   
  
 

 
  
   L
   
    T
    fore
   
   
    (
    i
    )
   
  
  =
  
   
    
     X
     ^
    
    
     T
     fore
    
    
     (
     i
     )
    
   
   -
   
    
     z
     
      (
      
       1
       -
       
        α
        
         P
         .
         I
        
       
      
      )
     
    
    ⁢
    
     
      
       1
       
        T
        fore
       
      
      ·
      
       
        ∑
        
         k
         =
         1
        
        
         T
         fore
        
       
       
        
         [
         
          
           
            X
            ^
           
           k
           
            (
            i
            )
           
          
          -
          
           
            μ
            ^
           
           
            (
            i
            )
           
          
         
         ]
        
        2
       
      
     
    
   
  
 

wherein z(1−αP.I) is z distribution at a significance level of 1−αP.I, UTfore(i) is an upper bound of an online traffic threshold interval of the i-th short cycle, LTfore(i) is a lower bound of the online traffic threshold interval, UTfore(i) and LTfore(i) each are a time series with a length of Tfore, and αP.I is a confidence level;
normal traffic data of the ICS communication network during the i-th iteration comprises:, {
   
    
     
      
       
        l
        j
        
         (
         i
         )
        
       
       ≤
       
        x
        j
        
         (
         i
         )
        
       
       ≤
       
        u
        j
        
         (
         i
         )
        
       
      
     
     
      
       (
       
        
         j
         =
         1
        
        ,
        2
        ,
        3
        ,
        …
           
        ,
        
         T
         fore
        
       
       )
      
     
    
    
     
      
       
        U
        
         T
         fore
        
        
         (
         i
         )
        
       
       =
       
        {
        
         
          u
          1
          
           (
           i
           )
          
         
         ,
         
          u
          2
          
           (
           i
           )
          
         
         ,
         
          u
          3
          
           (
           i
           )
          
         
         ,
         …
            
         ,
         
          u
          
           T
           fore
          
          
           (
           i
           )
          
         
        
        }
       
      
     
     
       
     
    
    
     
      
       
        L
        
         T
         fore
        
        
         (
         i
         )
        
       
       =
       
        {
        
         
          l
          1
          
           (
           i
           )
          
         
         ,
         
          l
          2
          
           (
           i
           )
          
         
         ,
         
          l
          3
          
           (
           i
           )
          
         
         ,
         …
            
         ,
         
          l
          
           T
           fore
          
          
           (
           i
           )
          
         
        
        }
       
      
     
     
       
     
    
   
  
 

wherein XTfore(i) is a traffic sequence for the ICS communication network collected in real time in the i-th short cycle, and Tfore is a sample size;
based on a real-time correlation between an ICS forecast sequence XTfore(i) and a training sequence XTtrin(i) in the i-th short cycle, estimating the forecast sequence XTfore(i) of the i-th short cycle:

{circumflex over (X)}Tfore(i)≈XTtrai(i+1)−(XTtrai(i+1)∩XTtrai(i))

wherein a function ∩ is used to get an intersection of two time series sets;
2.9) after traffic determination ends, proceeding with training iteration of a next short cycle, outputting a new optimal model and parameters adapted to the model, and determining newly input real-time traffic data; and
2.10) repeating the whole process until a specified number of iterations is reached.