您现在的位置是:网站首页> 内容页

算法之python实现近似熵、互近似熵算法

  • 菲律宾金佰利娱乐
  • 2019-03-29
  • 63人已阅读
简介理论基础近似熵?定义:近似熵是一个随机复杂度,反应序列相邻的m个点所连成折线段的模式的互相近似的概率与由m+1个点所连成的折线段的模式相互近似的概率之差。作用:用来描述复杂系统的不规则

理论基础

近似熵?

定义:近似熵是一个随机复杂度,反应序列相邻的m个点所连成折线段的模式的互相近似的概率与由m+1个点所连成的折线段的模式相互近似的概率之差。

作用:用来描述复杂系统的不规则性,越是不规则的时间序列对应的近似熵越大。反应维数改变时产生的新的模式的可能性的大小。

对于eeg信号来说,由于噪声存在、和信号的微弱性、多重信号源叠加,反映出来的是混沌属性,但是同一个人在大脑活动相对平稳的情况下,其eeg近似熵应该变化不大。

证明和对应几何意义可参考论文:https://wenku.baidu.com/view/4ec89e44b307e87101f696ef.html

互近似熵

从近似熵定义引申出来的,近似熵描述的是一段序列的自相似程度,互近似熵比较的是两段序列的复杂度接近程度;熵值越大越不相似,越小越相似;

近似熵算法分析

    设存在一个以等时间间隔采样获得的m维的时间序列u(1)u(2)...u(N).

    定义相关参数维数m一般取值为2,相似容限即阀值r,其中,维数表示向量的长度;r表示“相似度”的度量值.

    重构m维向量X(1)X(2)...X(N−m+1)其中X(i)=[u(i)u(i+1)...u(i+m−1)],X(j)=[u(j)u(j+1)...u(j+m−1)]计算X(i)和X(j)之间的距离,由对应元素的最大差值决定;d[XX∗]=maxa|u(a)−u∗(a)|d[XX∗]=maxa⁡|u(a)−u∗(a)|

    统计所有的d[XX∗]<=r的个数g则g/(N-M)就是本次的i取值对应的相似概率,计算所有i和j取值的概率对数的平均值,即熵值Φm(r);

    取m+1重复3、4过程,计算近似熵:

ApEn=Φm(r)−Φm+1(r)

参数选择:通常选择参数m=2或m=3;通常选择r=0.2∗std,其中std表示原时间序列的标准差.

互近似熵计算和近似熵的步骤一样,把计算X(i)和X(j)之间的距离改为计算序列a的向量X(i)和序列b的向量Y(j)的距离;相似容限r为两个原序列的0.2倍协方差;

python代码实现

使用Pincus提出的近似熵定义计算近似熵

class BaseApEn(object): """ 近似熵基础类 """ def __init__(self m r): """ 初始化 :param U:一个矩阵列表,for example: U = np.array([85 80 89] * 17) :param m: 子集的大小,int :param r: 阀值基数,0.1---0.2 """ self.m = m self.r = r @staticmethod def _maxdist(x_i x_j): """计算矢量之间的距离""" return np.max([np.abs(np.array(x_i) - np.array(x_j))]) @staticmethod def _biaozhuncha(U): """ 计算标准差的函数 :param U: :return: """ if not isinstance(U np.ndarray): U = np.array(U) return np.std(U ddof=1)class ApEn(BaseApEn): """ Pincus提出的算法,计算近似熵的类 """ def _biaozhunhua(self U): """ 将数据标准化, 获取平均值 所有值减去平均值除以标准差 """ self.me = np.mean(U) self.biao = self._biaozhuncha(U) return np.array([(x - self.me) / self.biao for x in U]) def _dazhi(self U): """ 获取阀值 :param U: :return: """ if not hasattr(self "f"): self.f = self._biaozhuncha(U) * self.r return self.f def _phi(self m U): """ 计算熵值 :param U: :param m: :return: """ # 获取矢量列表 x = [U[i:i + m] for i in range(len(U) - m + 1)] # 获取所有的比值列表 C = [len([1 for x_j in x if self._maxdist(x_i x_j) <= self._dazhi(U)]) / (len(U) - m + 1.0) for x_i in x] # 计算熵 return np.sum(np.log(list(filter(lambda a: a C)))) / (len(U) - m + 1.0) def _phi_b(self m U): """ 标准化数据计算熵值 :param m: :param U: :return: """ # 获取矢量列表 x = [U[i:i + m] for i in range(len(U) - m + 1)] # 获取所有的比值列表 C = [len([1 for x_j in x if self._maxdist(x_i x_j) <= self.r]) / (len(U) - m + 1.0) for x_i in x] # 计算熵 return np.sum(np.log(list(filter(lambda x: x C)))) / (len(U) - m + 1.0) def jinshishang(self U): """ 计算近似熵 :return: """ return np.abs(self._phi(self.m + 1 U) - self._phi(self.m U)) def jinshishangbiao(self U): """ 将原始数据标准化后的近似熵 :param U: :return: """ eeg = self._biaozhunhua(U) return np.abs(self._phi_b(self.m + 1 eeg) - self._phi_b(self.m eeg))if __name__ == "__main__": U = np.array([2 4 6 8 10] * 17) G = np.array([3 4 5 6 7] * 17) ap = ApEn(2 0.2) ap.jinshishang(U) # 计算近似熵

说明:

jinshishang函数直接计算近似熵

jinshishangbiao函数将原始数据标准化后计算近似熵

使用Pincus提出的近似熵定义计算互近似熵

class HuApEn(BaseApEn): def _xiefangcha(self U G): """ 计算协方差的函数 :param U: 序列1,矩阵 :param G: 序列2,矩阵 :return: 协方差,float """ if not isinstance(U np.ndarray): U = np.array(U) if not isinstance(G np.ndarray): G = np.array(G) if len(U) != len(G): raise AttributeError("参数错误!") return np.cov(U G ddof=1)[0 1] def _biaozhunhua(self U G): """ 对数据进行标准化 """ self.me_u = np.mean(U) self.me_g = np.mean(G) self.biao_u = self._biaozhuncha(U) self.biao_g = self._biaozhuncha(G) # self.biao_u = self._xiefangcha(U G) # self.biao_g = self._xiefangcha(U G) return np.array([(x - self.me_u) / self.biao_u for x in U]) np.array( [(x - self.me_g) / self.biao_g for x in U]) def _dazhi(self U G): """ 获取阀值 :param r: :return: """ if not hasattr(self "f"): self.f = self._xiefangcha(U G) * self.r return self.f def _phi(self m U G): """ 计算熵值 :param m: :return: """ # 获取X矢量列表 x = [U[i:i + m] for i in range(len(U) - m + 1)] # 获取y矢量列表 y = [G[g:g + m] for g in range(len(G) - m + 1)] # 获取所有的条件概率列表 C = [len([1 for y_k in y if self._maxdist(x_i y_k) <= self._dazhi(U G)]) / (len(U) - m + 1.0) for x_i in x] # 计算熵 return np.sum(np.log(list(filter(lambda x_1: x_1 C)))) / (len(U) - m + 1.0) def _phi_b(self m U G): """ 标准化数据计算熵值 :param m: :param U: :return: """ # 获取X矢量列表 x = [U[i:i + m] for i in range(len(U) - m + 1)] # 获取y矢量列表 y = [G[g:g + m] for g in range(len(G) - m + 1)] # 获取所有的条件概率列表 C = [len([1 for y_k in y if self._maxdist(x_i y_k) <= self.r]) / (len(U) - m + 1.0) for x_i in x] # 计算熵 return np.sum(np.log(list(filter(lambda x: x C)))) / (len(U) - m + 1.0) def hujinshishang(self U G): """ 计算互近似熵 :return: """ return np.abs(self._phi(self.m + 1 U G) - self._phi(self.m U G)) def hujinshishangbiao(self U G): """ 将原始数据标准化后的互近似熵 :param U: :param G: :return: """ u g = self._biaozhunhua(U G) return np.abs(self._phi_b(self.m + 1 u g) - self._phi_b(self.m u g))

使用洪波提出的快速实用算法计算近似熵

class NewBaseApen(object): """新算法基类""" @staticmethod def _get_array_zeros(x): """ 创建N*N的0矩阵 :param U: :return: """ N = np.size(x 0) return np.zeros((N N) dtype=int) @staticmethod def _get_c(z m): """ 计算熵值的算法 :param z: :param m: :return: """ N = len(z[0]) # 概率矩阵C计算 c = np.zeros((1 N - m + 1)) if m == 2: for j in range(N - m + 1): for i in range(N - m + 1): c[0 j] += z[j i] & z[j + 1 i + 1] if m == 3: for j in range(N - m + 1): for i in range(N - m + 1): c[0 j] += z[j i] & z[j + 1 i + 1] & z[j + 2 i + 2] if m != 2 and m != 3: raise AttributeError("m的取值不正确!") data = list(filter(lambda x:x c[0]/(N - m + 1.0))) if not all(data): return 0 return np.sum(np.log(data)) / (N - m + 1.0)class NewApEn(ApEn NewBaseApen): """ 洪波等人提出的快速实用算法计算近似熵 """ def _get_distance_array(self U): """ 获取距离矩阵 :param U: :return: """ z = self._get_array_zeros(U) fa = self._dazhi(U) for i in range(len(z[0])): z[i :] = (np.abs(U - U[i]) <= fa) + 0 return z def _get_shang(self m U): """ 计算熵值 :param U: :return: """ # 获取距离矩阵 Z = self._get_distance_array(U) return self._get_c(Z m) def hongbo_jinshishang(self U): """ 计算近似熵 :param U: :return: """ return np.abs(self._get_shang(self.m + 1 U) - self._get_shang(self.m U))

使用洪波提出的快速实用算法计算互近似熵

class NewHuApEn(HuApEn NewBaseApen): """ 洪波等人提出的快速实用算法计算互近似熵 """ def _get_distance_array(self U G): """ 获取距离矩阵 :param U:模板数据 :return:比较数据 """ z = self._get_array_zeros(U) fa = self._dazhi(U G) for i in range(len(z[0])): z[i :] = (np.abs(G - U[i]) <= fa) + 0 return z def _get_shang(self m U G): """ 计算熵值 :param U: :return: """ # 获取距离矩阵 Z = self._get_distance_array(U G) return self._get_c(Z m) def hongbo_hujinshishang(self U G): """ 对外的计算互近似熵的接口 :param U: :param G: :return: """ return np.abs(self._get_shang(self.m + 1 U G) - self._get_shang(self.m U G))

简单测试

if __name__ == "__main__": import time import random U = np.array([random.randint(0 100) for i in range(1000)]) G = np.array([random.randint(0 100) for i in range(1000)]) ap = NewApEn(2 0.2) ap1 = NewHuApEn(2 0.2) t = time.time() print(ap.jinshishang(U)) t1 = time.time() print(ap.hongbo_jinshishang(U)) t2 = time.time() print(ap1.hujinshishang(U G)) t3 = time.time() print(ap1.hongbo_hujinshishang(U G)) t4 = time.time() print(t1-t) print(t2-t1) print(t3-t2) print(t4-t3)

测试后发现使用快速算法比使用定义算法的计算效率提高了6倍以上。

参考:

https://wenku.baidu.com/view/4ec89e44b307e87101f696ef.html

http://blog.sina.com.cn/s/blog_6276ec79010118cx.html

https://blog.csdn.net/Cratial/article/details/79707169

1 0 9)

文章评论

Top