Phát hiện bất thường với Anomalize trong R

Anomalize: Mô tả ngắn

Khi nói về bất thường, chúng ta xét đến những điểm dữ liệu là ngoại lệ hoặc các sự kiện đặc biệt. Việc nhận diện các sự kiện này tương đối dễ trong những tập dữ liệu nhỏ và có thể thực hiện bằng một số biểu đồ phân tích đơn giản như boxplot. Tuy nhiên, mọi thứ sẽ trở nên phức tạp khi chuyển sang các tập dữ liệu lớn, đặc biệt là chuỗi thời gian. Chuỗi thời gian là dữ liệu được ghi nhận theo khoảng thời gian cố định trong một giai đoạn; khi phân tích sẽ cho thấy xu hướng hoặc tính mùa vụ. Việc xác định bất thường trong các trường hợp này là một khía cạnh khá khó.

Lúc này, gói anomalize dành cho phát hiện bất thường trong phân tích chuỗi thời gian phát huy tác dụng. Đây là thuật toán phát hiện bất thường theo phong cách tidy, dựa trên thời gian và có khả năng mở rộng từ một đến nhiều chuỗi thời gian.

Có một số gói và phương pháp sẵn có hỗ trợ cho quá trình phát triển của nó, hay có thể nói là sự kết hợp các nguồn lực hiện có với cách tiếp cận có khả năng mở rộng.

Các dự án mã nguồn mở đã hỗ trợ gồm:

• Gói AnomalyDetection của Twitter: Có trên Github (cran anaomalyDetection là một công việc khác).
• Hàm forecast::tsoutliers() của Rob Hyndman có trong gói forecast.
• Gói tsoutliers của Javier Lopez-de-lacalle trên CRAN.

Tất cả các gói và hàm này được dùng để tích hợp vào một quy trình làm việc có thể mở rộng.

Nói về quy trình làm việc của anomalize, nó được chia làm ba phần:

• Phân rã chuỗi thời gian với time_decompose().
• Phát hiện bất thường trên phần dư với anomalize().
• Chuyển đổi cận dưới và cận trên của bất thường với time_recompose().

Quy trình phát hiện bất thường

Phân rã chuỗi thời gian

Bước đầu tiên là phân rã chuỗi thời gian bằng time_decompose(). Giá trị đo được hay giá trị số cần phát hiện cho một nhóm cụ thể được phân rã thành bốn cột là observed, season, trend và remainder. Phương pháp mặc định dùng để phân rã là stl, tức phân rã theo mùa sử dụng bộ làm mượt Loess.

Hồi quy Loess là phương pháp phổ biến nhất dùng để làm mượt chuỗi thời gian biến động mạnh; nó khớp nhiều mô hình hồi quy trong các lân cận cục bộ, hay có thể nói dữ liệu được chia nhỏ và hồi quy được áp dụng cho từng phần. Điều này hữu ích trong chuỗi thời gian vì chúng ta biết biên của thời gian, vốn là biến X trong trường hợp này. Phương pháp này hoạt động tốt khi xu hướng chi phối tính mùa vụ của chuỗi thời gian.

Ở đây, xu hướng là sự tăng trưởng dài hạn diễn ra qua nhiều quan sát, còn tính mùa vụ là khuôn mẫu chu kỳ xảy ra theo chu kỳ ngày cho từng phút hoặc giờ, hoặc theo tuần.

Có một kỹ thuật thứ hai có thể dùng để phân rã theo mùa trong chuỗi thời gian dựa trên trung vị, đó là phương pháp Twitter cũng được dùng trong gói AnomalyDetection. Nó tương tự STL ở phần loại bỏ thành phần mùa vụ. Điểm khác là ở cách loại bỏ xu hướng: nó dùng trung vị từng đoạn của dữ liệu (một hoặc vài trung vị chia theo khoảng thời gian xác định) thay vì khớp bộ làm mượt. Phương pháp này hoạt động tốt khi tính mùa vụ chi phối xu hướng trong chuỗi thời gian.

Hãy nói về đầu ra của hàm time_decompose(), như đã nhắc ở trên, nó tạo ra 4 cột:

• observed: Giá trị thực tế.
• season: Thành phần mùa vụ hoặc chu kỳ. Mặc định là mùa vụ theo tuần.
• trend: Xu hướng dài hạn. Mặc định là khoảng 3 tháng.
• remainder: Dùng để phân tích ngoại lệ. Đơn giản là observed trừ đi cả season và trend.

Hàm time_decompose() có tham số merge; đặt TRUE để giữ dữ liệu gốc cùng với các cột được tạo ra.

# Using data package provided in the anomalize package and taking single time series of package purrr

purrr_package = tidyverse_cran_downloads%>%
  filter(package == "purrr")%>%
  ungroup()

purrr_anomaly  = purrr_package %>%
  time_decompose(count)

purrr_anomaly%>% glimpse()
Được hỗ trợ bởi

 ## Observations: 425
## Variables: 5
## $ date      <date> 2017-01-01, 2017-01-02, 2017-01-03, 2017-01-04, 201...
## $ observed  <dbl> 550, 1012, 1515, 1702, 1696, 1613, 860, 796, 2008, 2...
## $ season    <dbl> -2158.8893, 692.6571, 1087.5708, 1052.3294, 939.9377...
## $ trend     <dbl> 1496.712, 1511.009, 1525.307, 1539.604, 1553.901, 15...
## $ remainder <dbl> 1212.1777, -1191.6661, -1097.8773, -889.9334, -797.8...Được hỗ trợ bởi

Phát hiện bất thường trong phần dư

Sau khi phân tích chuỗi thời gian hoàn tất và phần dư đã có các đặc tính mong muốn để thực hiện phát hiện bất thường, bước này sẽ tạo thêm ba cột mới.

• remainder_l1: Giới hạn dưới của phần dư.
• remainder_l2: Giới hạn trên của phần dư.
• anomaly: Cột cho biết quan sát có phải bất thường hay không.

Các điểm bất thường là những điểm có đòn bẩy cao làm méo phân phối. Gói anomalize triển khai hai phương pháp chống chịu với các điểm có đòn bẩy cao:

• IQR: Khoảng tứ phân vị trong
• GESD: Kiểm định Sai lệch Chuẩn hóa Cực trị Tổng quát

IQR

Đây là phương pháp tương tự được dùng trong hàm tsoutliers() của gói forecast. Trong IQR, ta lấy phân phối và dùng khoảng tứ phân vị trong 25% và 75% để thiết lập phân phối của phần dư. Các giới hạn mặc định được đặt bằng hệ số 3 lần trên và dưới khoảng tứ phân vị trong; bất kỳ phần dư nào vượt quá giới hạn đều được xem là bất thường.

GESD

Trong GESD, các bất thường được đánh giá dần dần bằng cách loại bỏ các điểm vi phạm nặng nhất rồi tính lại thống kê kiểm định và các giá trị tới hạn; nói đơn giản là phạm vi được tính lại sau khi nhận diện bất thường theo cách lặp.

Cả IQR và GESD đều có ưu nhược điểm riêng. IQR nhanh hơn do không có vòng lặp, nhưng không chính xác bằng GESD vì các bất thường làm lệch trung vị, điều được loại bỏ trong GESD.

purrr_anomaly = purrr_anomaly%>%
  anomalize(remainder)

purrr_anomaly%>% glimpse()
Được hỗ trợ bởi

## Observations: 425
## Variables: 8
## $ date         <date> 2017-01-01, 2017-01-02, 2017-01-03, 2017-01-04, ...
## $ observed     <dbl> 550, 1012, 1515, 1702, 1696, 1613, 860, 796, 2008...
## $ season       <dbl> -2158.8893, 692.6571, 1087.5708, 1052.3294, 939.9...
## $ trend        <dbl> 1496.712, 1511.009, 1525.307, 1539.604, 1553.901,...
## $ remainder    <dbl> 1212.1777, -1191.6661, -1097.8773, -889.9334, -79...
## $ remainder_l1 <dbl> -4330.511, -4330.511, -4330.511, -4330.511, -4330...
## $ remainder_l2 <dbl> 4400.459, 4400.459, 4400.459, 4400.459, 4400.459,...
## $ anomaly      <chr> "No", "No", "No", "No", "No", "No", "No", "No", "...Được hỗ trợ bởi

Chuyển đổi cận dưới và cận trên của bất thường

Bước cuối của quy trình là tạo các cận dưới và cận trên xung quanh giá trị quan sát với time_recompose. Hàm này tái hợp season, trend, remainder_l1 và remainder_l2 thành các giới hạn mới gồm:

• recomposed_l1: Cận dưới của các ngoại lệ quanh giá trị quan sát.
• recomposed_l2: Cận trên của các ngoại lệ quanh giá trị quan sát.

purrr_anomaly = purrr_anomaly%>%
  time_recompose()

purrr_anomaly%>% glimpse()
Được hỗ trợ bởi

## Observations: 425
## Variables: 10
## $ date          <date> 2017-01-01, 2017-01-02, 2017-01-03, 2017-01-04,...
## $ observed      <dbl> 550, 1012, 1515, 1702, 1696, 1613, 860, 796, 200...
## $ season        <dbl> -2158.8893, 692.6571, 1087.5708, 1052.3294, 939....
## $ trend         <dbl> 1496.712, 1511.009, 1525.307, 1539.604, 1553.901...
## $ remainder     <dbl> 1212.1777, -1191.6661, -1097.8773, -889.9334, -7...
## $ remainder_l1  <dbl> -4330.511, -4330.511, -4330.511, -4330.511, -433...
## $ remainder_l2  <dbl> 4400.459, 4400.459, 4400.459, 4400.459, 4400.459...
## $ anomaly       <chr> "No", "No", "No", "No", "No", "No", "No", "No", ...
## $ recomposed_l1 <dbl> -4992.689, -2126.845, -1717.634, -1738.578, -183...
## $ recomposed_l2 <dbl> 3738.281, 6604.125, 7013.336, 6992.392, 6894.298...Được hỗ trợ bởi

Vẽ biểu đồ bất thường

purrr_anomaly%>%
  plot_anomaly_decomposition()+
  ggtitle("Plotting Anomalies")
Được hỗ trợ bởi

Điều chỉnh tham số

Mô hình hóa bộ phát hiện bất thường sẽ không thể hoàn chỉnh nếu không điều chỉnh các tham số, vốn phụ thuộc hoàn toàn vào dữ liệu.

Hãy cùng điều chỉnh tham số: các tham số ở mỗi tầng của quy trình là khác nhau vì mỗi tầng đảm nhiệm một tác vụ riêng.

Điều chỉnh tham số cho phân rã

Như đã trao đổi, xu hướng và tính mùa vụ là nền tảng để phân rã chuỗi thời gian, nên cũng cần điều chỉnh tần suất và xu hướng trong hàm time_decompose. Theo mặc định, các giá trị được gán tự động: 7 ngày cho tần suất ở cả hai phương pháp (STL, Twitter), và cho xu hướng là 91 ngày với STL và 85 ngày với Twitter.

Bạn có thể tinh chỉnh một hoặc cả hai tham số tùy nhu cầu, nhưng hãy quan sát kỹ trước khi điều chỉnh vì thay đổi thiếu thận trọng có thể dẫn đến quá khớp hoặc thiếu khớp trong quá trình phân rã.

purrr_package %>%
  time_decompose(count, frequency = "auto", trend = "2 weeks")%>%
  anomalize(remainder)%>%
  plot_anomaly_decomposition()+
  ggtitle("Trend = 2 Weeks / Frequency = auto ")
Được hỗ trợ bởi

Điều chỉnh tham số cho phát hiện bất thường

Như đã nêu, ở đây bất thường được quyết định theo các giá trị phần dư bằng cách tính các giới hạn để phân loại ngoại lệ. Hai tham số alpha và max_anoms điều khiển hàm anomalize(). Nói đơn giản, alpha điều khiển băng giới hạn; mặc định đặt 0,05, giảm giá trị này sẽ làm tăng độ rộng băng, từ đó khiến một điểm khó bị coi là bất thường hơn.

purrr_package%>%
  time_decompose(count)%>%
  anomalize(remainder, alpha = 0.05)%>%
  time_recompose()%>%
  plot_anomalies(time_recompose = T)+
  ggtitle("alpha = 0.05")
Được hỗ trợ bởi

purrr_package%>%
  time_decompose(count)%>%
  anomalize(remainder, alpha = 0.025)%>%
  time_recompose()%>%
  plot_anomalies(time_recompose = T)+
  ggtitle("alpha = 0.025")
Được hỗ trợ bởi

Tham số max_anom kiểm soát tỷ lệ phần trăm dữ liệu có thể là bất thường. Tham số này hữu ích khi alpha khó tinh chỉnh, và bạn muốn tập trung vào các bất thường mạnh nhất.

purrr_package%>%
  time_decompose(count)%>%
  anomalize(remainder, alpha = 0.2, max_anoms = 0.2)%>%
  time_recompose()%>%
  plot_anomalies(time_recompose = T)+
  ggtitle("20% anomaly Allowed")
Được hỗ trợ bởi

purrr_package%>%
  time_decompose(count)%>%
  anomalize(remainder, alpha = 0.2, max_anoms = 0.05)%>%
  time_recompose()%>%
  plot_anomalies(time_recompose = T)+
  ggtitle("5% anomaly Allowed")
Được hỗ trợ bởi

Tài liệu tham khảo

• Tham chiếu chi tiết tất cả hàm của Anomalize
• Bài blog thú vị với đầy đủ chi tiết

Nếu bạn muốn học thêm về R, hãy tham gia khóa Giới thiệu Phân tích Chuỗi thời gian của DataCamp.

Khám phá Hướng dẫn Phân tích Chuỗi thời gian bằng R của chúng tôi.