时间序列分类方法：SAX-VSM学习
前言
SAX是一种通过将时间序列离散化后投影到字符的算法，简单来说，SAX利用分段聚合近似PAA将时间序列平均分为n个段，再将该段用平均值表示，根据正态分布的分位数将这些值映射为一个字符。
VSM(Vector Space Model)，它是信息检索和自然语言处理领域中常用的一种模型。VSM 将文档表示为高维向量，每个维度代表一个特定的词语或术语，从而将文档转换为数值化的形式，便于计算机进行处理和分析。
SAX-VSM利用了SAX将时间序列转换为字符的特点，将自然语言处理的TF-IDF方法应用其中，实现了时间序列分类的效果。
本文主要为个人的学习记录以及实际测试，和大家一起学习qwq
一、理论部分
长度为n的时间序列 T T T： { t 1 , t 2 , . . . , t n } \{t_1,t_2,...,t_n\} {t1,t2,...,tn}
1.SAX
-
分段聚合近似PAA
- 标准化：将时间序列 T T T进行z标准化： T = T − μ σ T=\frac{T-\mu}{\sigma} T=σT−μ，这样的目的是为了将时间序列划分到统一尺度下，消除尺度差异；
- 分段：将时间序列 T T T平均划分为 m m m个段（这里假设可以n可以被m整除，在实际应用中会对不能整除的情况进行处理）；
- 聚合：计算每个段的平均值 t ˉ i = 1 m ∑ j = n m ( i − 1 ) + 1 n m i t j \bar t_i = \frac{1}{m}\sum_{j=\frac{n}{m}(i-1)+1}^{\frac{n}{m}i}t_j tˉi=m1∑j=mn(i−1)+1mnitj；
- 最终得到的离散序列：
T
ˉ
=
{
t
ˉ
1
,
t
ˉ
2
,
.
.
.
,
t
ˉ
w
}
\bar T=\{\bar t_1,\bar t_2,...,\bar t_w\}
Tˉ={tˉ1,tˉ2,...,tˉw}。
如上图中，蓝色曲线为原始数据，红色为每一段的平均值。
-
字符映射
原始的SAX在选择字符数量上限后根据正态分布的分位数进行分配，如下：
即:
s j = a l p h a j , i f f β j − 1 ≤ t ˉ i ≤ β j s_j = alpha_j,iff \beta_{j-1}\le \bar t_i \le \beta_j sj=alphaj,iffβj−1≤tˉi≤βj
这里主要讲的是SAX-VSM实现分类，就不详细去讲如何对字符进行距离比较了。
原始的时间序列经过SAX进行降维，减少了数据量，并且在分段聚合近似的过程中也实现了一定的降噪效果，但是需要注意的是，SAX是依靠均值来分配字符的，这导致人直观感觉到的和实际计算的有些出入，且会因为一些细小的变化导致分配的字符不同。
如上图所示，这是两个具有一定“相位差”的子序列，由于幅值的差别，导致离散化字符的差别很大。这是由于划分字符的规则是由分段的平均值决定，虽然第二段有显著的尖峰，但是平均下来不如总体略高的第四段（上）和第一段（下）。
2.词袋(Bag-of-Word)
对于一个长时间序列，往往需要利用滑窗的方法来将其中每一个子序列转换为SAX字符串（单词），为了避免相邻子序列因高度相似导致的相同字符串，当出现连续相同的字符串时仅保留第一个单词。
这样我们将一个长时间序列转换为了一个单词序列：
T
=
{
t
1
,
t
2
,
.
.
.
,
t
n
}
⇒
S
=
{
w
o
r
d
1
,
w
o
r
d
2
,
.
.
.
,
w
o
r
d
j
}
T=\{t_1,t_2,...,t_n\}\Rightarrow S=\{word_1,word_2,...,word_j\}
T={t1,t2,...,tn}⇒S={word1,word2,...,wordj}
此时，我们可以利用直方图的方式将所有的单词表示出来：
利用这种词袋方法的一个好处是可以实现时间序列的相位不变性，即将时间序列转换为无序的词袋后进行比较。利用直方图也可以进行相似或者分类的任务，但不是本文的重点。
3.词频-逆文档频率 TF-IDF
术语频率指的是一个词语在文档中出现的次数，如果一个词在一个文档中出现次数越多，对文档重要性越大。
t
f
t
,
d
=
{
l
o
g
(
1
+
f
t
,
d
)
,
i
f
f
t
,
d
>
0
0
,
o
t
h
e
r
w
i
s
e
tf_{t,d}=\left\{ \begin{aligned} log(1+f_{t,d}),if\ f_{t,d}>0\\ 0,otherwise \end{aligned} \right.
tft,d={log(1+ft,d),if ft,d>00,otherwise
简单来说，
t
t
t就是每一种单词，
d
d
d就是一类数据的单词集合，
f
t
,
d
f_{t,d}
ft,d是单词的词频。
逆文档频率的思想是从全局（所有类别的数据）来看，如果一个单词在很少的类别中出现，则它对于该类来说很重要：
i
d
f
t
,
D
=
l
o
g
∣
D
∣
∣
d
∈
D
:
t
∈
d
∣
=
l
o
g
N
d
f
t
idf_{t,D}=log\frac{|D|}{|d\in D:t\in d|}=log\frac{N}{df_t}
idft,D=log∣d∈D:t∈d∣∣D∣=logdftN
其中
N
N
N是数据的类别数，
d
f
t
df_t
dft是单词
t
t
t出现过的词袋数。
则最终
t
f
−
i
d
f
tf-idf
tf−idf：
t
f
∗
i
d
f
(
t
,
d
,
D
)
=
t
f
t
,
d
×
i
d
f
t
,
D
=
l
o
g
(
1
+
f
t
,
d
)
×
l
o
g
N
d
f
t
tf*idf(t,d,D)=tf_{t,d}\times idf_{t,D}=log(1+f_{t,d})\times log\frac{N}{df_t}
tf∗idf(t,d,D)=tft,d×idft,D=log(1+ft,d)×logdftN
一旦计算出所有频率值，单词频率矩阵就成为词频-逆文档频率权重矩阵，每一行表示对应的单词，每一列对应一类数据的权重向量，将其归一化处理。
将未知类别的时间序列转换为词袋后，在转换为词频向量，通过余弦相似度计算来划分类别：
s
i
m
i
l
a
r
i
t
y
(
a
,
b
)
=
c
o
s
(
θ
)
=
a
⋅
b
∣
∣
a
∣
∣
⋅
∣
∣
b
∣
∣
similarity(a,b)=cos(\theta)=\frac{a\cdot b}{||a||\cdot ||b||}
similarity(a,b)=cos(θ)=∣∣a∣∣⋅∣∣b∣∣a⋅b
如果未知类别的词袋中出现了TF-IDF矩阵中不存在的单词则忽略。
二、实际应用
个人学习的习惯是学到一种新的算法后，按照伪代码编出来，再跟一些标准的函数库pyts进行对比，SAX算法其原创团队已经提供了代码，但是由于其正态分位数是直接设定了几位有效数字，因此可能和理论上分配的单词稍微有所偏差，但是问题不大。
1.自编代码
1.1 Matlab部分代码：
function [saxseries,indexs] = SAX(x,winSize,wordSize,num,stepSize)
n = length(x);
beta = [-inf,norminv((1:num-1)/num,0,1)]; % 正态分布的分位数
temp = '';
alpha = ['a','b','c','d','e','f','g','h','i','j','k','l','m','n'];
tempPAA = PAA(x,winSize,wordSize,stepSize); % 预先计算好PAA得到的均值
saxseries = "";
indexs = [];
for i = 1:size(tempPAA,1)
for j = 1:wordSize
temp(j) = alpha(sum(beta<=tempPAA(i,j),2)); % 以当前均值大于分位数的个数作为索引，从字母表中取出对应的字母
end
if length(saxseries)==1||~all(temp==saxseries(end)) % 连续相同的单词仅保存第一个
saxseries(end+1) = temp;
indexs(end+1) = i;
end
temp = '';
end
saxseries(1) = []; % 第一个字符是""
end
上述代码为Matlab实现SAX的函数，其中PAA为预先计算各个子序列均值的函数，由于对每一个子序列都进行标准化非常消耗时间，因此利用了名为MASS的提前计算统计量的思想，可以提前计算许多标准化需要的均值和标准差，从而减少运行时间。
第一部分提到如果分段数不能整除，则可以按照如下的方式进行分配：
if m/N~=ceil(m/N)
temp = zeros(N,m);
for j = 1 : N
temp(j, :) = sub;
end
sub = reshape(temp, 1, N*m);
sub = reshape(sub,m,N);
else
sub = reshape(sub,m/N,N);
end
这里的思想就是将数据重复后进行叠加分配计算均值：
T
=
[
1
,
2
,
3
,
4
,
5
,
6
,
7
,
8
,
9
,
10
]
T=[1,2,3,4,5,6,7,8,9,10]
T=[1,2,3,4,5,6,7,8,9,10]共十个数据，如果要划分为3个段则将数据重复为：
T
′
=
[
1
2
3
4
5
6
7
8
9
10
1
2
3
4
5
6
7
8
9
10
1
2
3
4
5
6
7
8
9
10
]
T'=\begin{bmatrix}1&2&3&4&5&6&7&8&9&10\\1&2&3&4&5&6&7&8&9&10\\1&2&3&4&5&6&7&8&9&10\\\end{bmatrix}
T′=
111222333444555666777888999101010
再将其重新排列：
T
′
′
=
[
1
4
7
1
4
8
1
5
8
2
5
8
2
5
9
2
6
9
3
6
9
3
6
10
3
7
10
4
7
10
]
T''=\begin{bmatrix}1&4&7\\1&4&8\\1&5&8\\2&5&8\\2&5&9\\2&6&9\\3&6&9\\3&6&10\\3&7&10\\4&7&10\end{bmatrix}
T′′=
111222333444555666777888999101010
然后计算这三列数据的平均值即可。
后续需要计算单词数量、每个类别的单词数量等等，我比较习惯python的字典，因此SAX-VSM在python实现。
1.2 Python部分代码：
转换代码：
def trans(self, x, y):
label = sorted(tuple(set(y))) # 标签
allclass = {key: {} for key in label} # 每一个单词出现数量的字典
count_class = {key: 0 for key in label} # 每一个单词在不同类中出现数量的字典
tf = {key: {} for key in label}
for i in range(np.size(x, 0)):
data = x[i, :]
sax = self.sax(data, self.win_size, self.word_len)
un_sax, counts = np.unique(sax, return_counts=True)
for key, count in zip(un_sax, counts):
self.wordlist[key] = self.wordlist.get(key, 0) + count
allclass[y[i]][key] = allclass[y[i]].get(key, 0) + count
class_list[y[i]][key] = class_list[y[i]].get(key, 0) + count
N = len(label)
# 计算IDF
self.dft = {}
for key, values in self.wordlist.items():
n = 0
for i in range(len(self.label)):
if key in class_list[self.label[i]]:
n += 1
self.dft[key] = np.log(N / n) + 1
for i in range(len(self.label)):
dic = class_list[self.label[i]]
if key in dic:
tf[self.label[i]][key] = tf[self.label[i]].get(key, 0) + dic[key]
# 计算tf-idf矩阵
self.allwords = list(self.wordlist.keys())
label_len = len(self.label)
allwords_len = len(self.allwords)
# 预计算tf中每个单词的log值
tf_log_values = {label: {word: np.log(tf_count) + 1 for word, tf_count in tf[label].items()}
for label in self.label}
self.tf_idf = np.zeros((allwords_len, label_len))
for j, label in enumerate(self.label):
dic = class_list[label]
for key in dic.keys():
if key in self.allwords:
i = self.allwords.index(key)
self.tf_idf[i, j] = self.dft[key] * tf_log_values[label][key]
和pyts的SAX-VSM库对比了结果，完全相同，但是官方代码似乎是用了稀疏矩阵等方式，运行速度受数据量大小影响更小，我这种方法对于样本点数少的速度更快。
2.Pyts库函数
代码如下：
from pyts.classification import SAXVSM
from pyts.datasets import load_gunpoint,load_pig_central_venous_pressure,load_coffee
# 导入数据
X_train, X_test, y_train, y_test = load_gunpoint(return_X_y=True)
# 滑窗宽度为70，单词长度为10，单词表长度为10，分箱策略为正态分布，这里threshold_std是指当子序列标准差小于阈值则在标准化过程不除以标准差了
clf = SAXVSM(window_size=70,word_size=10, n_bins=10, strategy='normal', numerosity_reduction=True,threshold_std=0.01)
clf.fit(X_train, y_train)
score = clf.score(X_test, y_test)
print(score)
测试结果：
总结
本文简述了SAX-VSM的基本内容，并且展示了部分自己编的代码和pyts库函数的应用，个人体验来说SAX-VSM还是挺有趣的，因为SAX转换允许将很多相似形状的时间序列用一类单词表示，再根据它的索引可以一次性筛选出这些子序列，又由于转换为了词袋，使用了TF-IDF来进行分类任务，在一些不是特别复杂的数据里，其效果还是比较好的。但是需要注意，随着单词长度和单词表范围的增长，潜在的单词类别呈指数上升，如我自己编的代码就会大大增加计算量，但是用库函数的话还好。
这是我第一次写文章，可能内容差一些，还是希望能和大家一起学习。
参考文献
标签：10,self,label,tf,SAX,key,序列,VSM From： https://blog.csdn.net/weixin_48841772/article/details/137543395Lin J, Keogh E, Lonardi S, et al. A symbolic representation of time series, with implications for streaming algorithms[C]//Proceedings of the 8th ACM SIGMOD workshop on Research issues in data mining and knowledge discovery. 2003: 2-11.
Senin P, Malinchik S. Sax-vsm: Interpretable time series classification using sax and vector space model[C]//2013 IEEE 13th international conference on data mining. IEEE, 2013: 1175-1180.