做EEG频谱分析,看这一篇文章就够了!

2023-11-01

所谓频谱分析,又称为功率谱分析或者功率谱密度(Power Spectral Density, PSD)分析,实际就是通过一定方法求解信号的功率power随着频率变化曲线。笔者在这里对目前常用的频谱分析方法做一个总结,并重点介绍目前EEG分析中最常用的频谱分析方法,并给出相应的Matlab程序。

1.频谱分析的方法有哪些?
目前来说,功率谱分析的方法大致可以分为两大类:第一类是经典的功率谱计算方法,第二类是现代功率谱计算方法,如图1所示。
其中第一类经典功率谱分析方法,又可以分为直接法、间接法和改进的直接法。直接法又称之为周期图法,简单地说,其直接利用信号的傅里叶变换系数的幅度平方来计算信号的功率谱。间接法又称为自相关函数法,其先估算出信号的自相关函数,然后对自相关函数求傅里叶变换从而得到信号的功率谱。改进的直接法,是针对直接法存在的缺点改进而来的方法,包括Barlett法、Welch法和Nuttall法。
第二类现代功率谱计算方法,又可以分为基于参数建模的功率谱计算和基于非参数建模的功率谱计算。基于参数建模的功率谱计算方法又分为基于AR模型、MA模型、ARMA模型等方法;基于非参数建模的功率谱计算方法主要基于矩阵特征分解的功率谱估计,主要包括基于MUSIC算法的功率谱估计和基于特征向量的功率谱估计。
本文,笔者主要对经典功率谱分析方法中的直接法(周期图法)以及在EEG频谱分析中最常用的改进直接法(Welch法)进行详细的介绍,并给出相应的Matlab程序。

在这里插入图片描述
图1

2.直接法计算PSD
直接法又称之为周期图法,是由Schuster于1899年提出,其方法很简单:对于长度为N的序列x(n),其傅里叶变换为X(k),那么x(n)在每个频率(或k)处的功率可以表示为
在这里插入图片描述
即傅里叶变换系数的幅值平方除以数据长度N(注意:上述公式是Power或者功率的计算方法,对于功率谱密度PSD来说,还需要除以采样率!)。
根据直接法求解PSD的定义,可以直接通过调用Matlab中的fft函数(fft函数是计算信号的傅里叶变换)进行计算;
此外,Matlab中有专门的函数periodogram实现直接法的PSD计算。
例1:按照直接法计算PSD的定义,利用Matlab中的fft函数直接计算信号y=5sin(2pif1t)+3cos(2pif2*t)+δ(其中f1=20Hz,f2=35Hz,δ为一随机噪声)的PSD。
结果如图2所示,Matlab程序可以在公众号后台输入“PSDcode”进行下载,下载后可以直接在Matlab中运行出以下结果。
Matlab中有专门的函数periodogram,也可以实现直接法的PSD计算,关于其用法,这里笔者就不再赘述,各位可以直接在Matlab中help一下这个函数即可,其使用方法很简单。

在这里插入图片描述图2

3.直接法计算PSD
上述直接法计算PSD,虽然可以直接FFT,计算速度快,但是频率分辨率比较低。因此,研究者提出了改进的直接法来计算PSD。**其中Welch方法就是其中的一种。Welch方法的思路是:**先把长度为N的信号分成L段,每段数据长为M,则N=LM;然后把窗函数w加到每段数据上,求出每段数据的功率谱;最后对每段数据的功率谱进行平均,得到整个信号的功率谱。各个数据段之间可以有重叠,窗函数w可以选择如Hanning、Hamming等任意一种窗口。
Welch方法是EEG的PSD计算中最常用的一种方法,在Matlab中直接采用pwelch函数可以实现Welch方法的功率谱估计,其一般调用格式如下:
[Pxx,F] = pwelch(X,WINDOW,NOVERLAP,NFFT,Fs)
[Pxx,F] = pwelch(X,WINDOW,NOVERLAP,NFFT,Fs,REQRANGE)
关于函数中的参数含义,各位可以在Matlab命令窗口中输入“help pwelch”即可查看其详细用法,这里笔者不再赘述。
例2:采用Welch方法计算信号y=5sin(2pif1t)+3cos(2pif2*t)+δ(其中f1=20Hz,f2=35Hz,δ为一随机噪声)的PSD。
计算结果如图3所示,与图2利用直接法求得的PSD基本一样;Matlab程序可以在公众号后台输入“PSDcode”进行下载,下载后可以直接在Matlab中运行出以下结果。

在这里插入图片描述
图3

4.总结
本文首先对目前进行PSD计算的不同方法进行了总结和简单介绍,重点详细介绍了如何利用直接法和改进的直接法(Welch法)来计算信号的PSD,并给出了Matlab程序。其中Welch法是目前计算EEG的PSD最常用的方法之一,理解和学会使用Welch法进行频率分析对于我们做EEG研究来说至关重要。

本文来自**河南悦影医药科技有限公司(简称“悦影科技”)**旗下以脑科学和神经科学为主题的公益性自媒体平台——“脑之说”微信公众号,欢迎各位朋友搜索关注“脑之说”公众号。本文为悦影科技原创技术文章,转载请联系赵老师(微信号:kervin_zhao;邮箱:yueyingkj@163.com)。

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

做EEG频谱分析,看这一篇文章就够了! 的相关文章

  • 如何在Matlab中自定义轮廓线?

    我正在准备一个等高线图 我应该在其中突出显示特定级别的等高线 例如 我的轮廓线值位于 1 和 1 之间 我想突出显示与值 0 相对应的线 我尝试使用以下过程来执行此操作 M c contourf longitude latitude del
  • MATLAB中如何画水平线和垂直线?

    我目前正在尝试在 MATLAB 中绘制简单的垂直线和水平线 例如 我想绘制线 y 245 我该怎么做呢 MATLAB 根据您提供的向量逐点进行绘图 因此 要创建一条水平线 您需要改变x同时保持y对于垂直线恒定 反之亦然 xh 0 10 yh
  • 两个 y 轴与相同的 x 轴[重复]

    这个问题在这里已经有答案了 可能的重复 在单个图中绘制 4 条曲线 具有 3 个 y 轴 https stackoverflow com questions 1719048 plotting 4 curves in a single plo
  • MATLAB 滚动图

    我有一个脑电图数据库 我想绘制它 数据库是一个19 1000 134的矩阵 其中 19 是通道数 在第一种方法中 我只使用一个渠道 1000 个样本大小 采样率为 500 Hz 时为 1000 个点 即 2 秒数据 134 epochs的数
  • 在 matlab 中求 3d 峰的体积

    现在我有一个带有峰值的 3D 散点图 我需要找到其体积 我的数据来自图像 因此 x 和 y 值表示 xy 平面上的像素位置 z 值是每个像素的像素值 这是我的散点图 scatter3 x y z 20 z filled 我试图找到数据峰值的
  • 使用mat2cell将MxN的矩阵划分为1xN大小的M矩阵

    我有一个大小为 MxN 的矩阵 比方说 1867x3 1867 行和 3 列 我想将其分成 1867 个大小为 1x3 的单元格 我使用了mat2cell X 1 1866 这里X是矩阵 1867x3 结果给出了两个单元格 一个单元格的大小
  • MATLAB - 如何将子图一起缩放?

    我在一张图中有多个子图 每个图的 X 轴是相同的变量 时间 每个图上的 Y 轴都不同 无论是它所代表的内容还是数据的大小 我想要一种同时放大所有图的时间尺度的方法 理想情况下 可以在其中一张图上使用矩形缩放工具 并让其他图相应地更改其 X
  • 在 MATLAB 中重命名文件

    我正在尝试以编程方式重命名工作目录中的文件a temp txt to b hello txt 您建议如何这样做 MATLAB中有一个简单的文件重命名函数吗 我认为您正在寻找 MOVEFILE
  • 如何在 Matlab 中将数组打印到 .txt 文件?

    我才刚刚开始学习Matlab 所以这个问题可能非常基本 我有一个变量 a 2 3 3 422 6 121 9 4 55 我希望将值输出到 txt 文件 如下所示 2 3 3 422 6 121 9 4 55 我怎样才能做到这一点 fid f
  • 在 Matlab 中将 datenum 转换为 datetime 的最快方法

    我在 Matlab 中将 datenum 转换为 datetime 时遇到问题 Given dnum floor now floor now 1 我尝试了以下方法 datenum dnum 但这没有用 我发现有效的方法是 datetime
  • Matlab:如何更改矩阵的存储方式?从 1x1x3 到 1x3?

    我目前有 val 1 0 7216 val 2 0 7216 val 3 0 7216 但我想要 0 7216 0 716 0 721 我可以做什么样的操作来做到这一点 The reshape函数将在这里解决问题 Arrange the e
  • 黑白随机着色的六角格子

    我正在尝试绘制一个 10 000 x 10 000 随机半黑半白的六边形格子 我不知道如何将该格子的六边形随机填充为黑色和白色 这是我真正想要从这段代码中得到的示例 但我无法做到 https i stack imgur com RkdCw
  • 在 Matlab 中高效获取像素坐标

    我想在 Matlab 中创建一个函数 给定一个图像 该函数将允许人们通过单击图像中的像素来选择该像素并返回该像素的坐标 理想情况下 人们能够连续单击图像中的多个像素 并且该函数会将所有相应的坐标存储在一个矩阵中 有没有办法在Matlab中做
  • matlab中的排列函数是如何工作的

    这是一个有点愚蠢的问题 但我似乎无法弄清楚排列在 matlab 中是如何工作的 以文档为例 A 1 2 3 4 permute A 2 1 ans 1 3 2 4 到底是怎么回事 这如何告诉 matlab 3 和 2 需要交换 哇 这是我迄
  • 访问图像的 Windows“标签”元数据字段

    我正在尝试进行一些图像处理 所以现在我正在尝试读取图像 exif 数据 有 2 个内置函数可用于读取图像的 exif 数据 问题是我想读取图像标签 exifread and imfinfo这两个函数都不显示图像标签 Is there any
  • “Desort”向量(撤消排序)

    在Matlab中 sort返回排序后的向量和索引向量 显示哪个向量元素已移动到以下位置 v ix sort u Here v是一个包含所有元素的向量u 但已排序 ix是一个向量 显示每个元素的原始位置v in u 使用 Matlab 的语法
  • 如何在MATLAB中显示由三个矩阵表示的图像?

    我有 3 个相同大小的 2D 矩阵 假设 200 行和 300 列 每个矩阵代表三种 基本 颜色 红色 绿色和蓝色 之一的值 矩阵的值可以在 0 到 255 之间 现在我想组合这些矩阵以将它们显示为彩色图像 200 x 300 像素 我怎样
  • 如何将复杂的 csv 文件导入到 Matlab 中的数值向量

    我想知道我们应该如何读取由字符串 双精度数和字符等组成的复杂 csv 文件 例如 您能否提供一个可以在此 csv 文件中提取数值的成功命令 Click here http www ecb europa eu stats money yc d
  • 如何在 MATLAB 中绘制纹理映射三角形?

    我有一个三角形 u v 图像中的坐标 我想在 3D 坐标处绘制这个三角形 X Y Z 与图像中的三角形进行纹理映射 Here u v X Y Z都是具有三个元素的向量 代表三角形的三个角 我有一个非常丑陋 缓慢且令人不满意的解决方案 其中我
  • MATLAB:在不使用循环的情况下提取矩阵的多个部分

    我有一个巨大的 2D 矩阵 我想从中提取 15 个不同的 100x100 部分 我有两个向量 x 和 y 其中保存了零件的左上角索引 我用过这样的东西 result cam1 x 1 end x 1 end 99 y 1 end y 1 e

随机推荐

  • jenkins - Manage and Assign Roles

    Role Strategy Plugin 插件 针对多个project进行权限控制 访问 上几张图 希望你能看明白 哈哈 1 png 710dba0dgy1fkgqp3cze1j219g0kmn24 jpg 710dba0dgy1fkgqp
  • MySQL查询语句in子查询的优化

    项目中有需要 使用MySQL的in子查询 查询符合in子查询集合中条件的数据 但是没想到的是 MySQL的in子查询会如此的慢 让人无法接收 于是上网搜索解决办法 下面记录下 一 原始in子查询 SELECT FROM basic zdjb
  • Ubuntu系统上安装WPS

    前言 在Ubuntu系统下 想使用WPS的功能 觉得用起来更加方便 所以在此记录一下安装的步骤 记录两种安装方法 方法一 Ubuntu Software中搜索WPS 如图所示 在Ubuntu Software中搜索WPS 可能需要稍等一会再
  • python使用局部敏感性哈希算法,在海量数据中查询相似序列

    文章目录 一 原生python实现 二 第三方库datasketch使用 1 官方示例 2 LSH算法 3 MinHashLSHForest 局部敏感性哈希是指 相似的哈希具有相似的原始序列 整体思路 首先将数据装在不同的桶里 通过桶之间的
  • 2023国赛数学建模思路 - 案例:随机森林

    文章目录 1 什么是随机森林 2 随机深林构造流程 3 随机森林的优缺点 3 1 优点 3 2 缺点 4 随机深林算法实现 建模资料 0 赛题思路 赛题出来以后第一时间在CSDN分享 https blog csdn net dc sinor
  • 隐私计算S2赛季-谁是真正的王者

    去年至今 隐私计算大约经历了如火如荼的一年 身为局中人 看穿居中事 道尽居中话 为的无非是让更多的来了解这个比较细分的AI领域 秋天本是硕果累累的丰收季 隐私计算这个行业算是金秋吗 一喜一悲 一喜为百花齐放 我所知道在布 挂 局 钩 隐私计
  • VL53L0X调试总结

    最近调VL53L0X花了不少时间 特总结下 https www st com content st com en search html q vl53l t products page 1 VL53L0X测距2m VL53L1X测距4m 支
  • networkx 中文学习手册

    文章目录 创建图表 节点 边 检查图的元素 从图中删除元素 使用图构造函数 什么用作节点和边 访问边和邻居 向图 节点和边添加属性 图形属性 节点属性 边缘属性 多图 图生成器和图操作 1 应用经典的图操作 例如 2 使用对经典小图之一的调
  • Harmony OS WiFi编程——连接热点、创建热点

    本节主要介绍如何在HiSpark WiFi IoT套件上使用Hamony OS的WiFi相关编程接口 相关知识点 WiFi的工作模式 AP模式 热点模式 提供无线接入服务 允许其它无线设备接入 提供数据访问 一般的无线路由 网桥工作在该模式
  • JavaNote 1.7final、finally、访问权限

    一 final 1 final的变量的值不能被改变 2 final的方法不能被重写 3 final的类不能被继承 二 finally finally 语句块 必须执行 通常在finally语句块中执行资源清除工作 如关闭打开的文件 删除临时
  • 基于Sklearn实现LDA算法

    文章目录 一 LDA算法 二 sklearn实现LDA 三 结果如图 四 总结 五 参考 一 LDA算法 1 线性判别分析 Linear Discriminant Analysis LDA 方法常被用于数据预处理中的降维 dimension
  • ArcMap布局视图的图例设置,如显示符号、标注、描述等

    转载请注明作者 独孤尚良dugushangliang 出处 https blog csdn net dugushangliang article details 81305762 如上图所示 为了达到以上的图例显示效果 鄙人上下求索 废了不
  • vue动态生成二维码,扫码登录

    首先找到对应的的三个接口 1 二维码获取key接口 接口说明 调用此接口可生成一个key 2 二维码生成接口 接口说明 调用此接口传入上一个接口生成的key可生成二维码图片的base64和二维码信息 可使用base64展示图片 或者使用二维
  • mysql自动启动设置用Systemctl start mysqld启动

    1 如果你是用yum安装的话就不需要进行设置了 用systemctl restart mysqld启动报如下错 2 查看系统服务有没有mysqld chkconfig list 3 MySQL启动关闭添加到 etc init d mysql
  • 解决Chrome浏览器中部分字体显示模糊的问题

    如果在Chrome浏览器中查看某些网页时 发现大部分字体显示清晰 但是另外部分字体显示模糊看不清的话 有可能是浏览器字体设置的问题 解决方式如下 1 点击Chrome浏览器右上角的 按钮 点击 设置 菜单 或直接在地址栏中输入 chrome
  • 解决java.net.BindException: Address already in use(Bind failed)端口占用问题

    sudo lsof i 8080 删掉图中两个进程 kill 9 2960 其中 9是九
  • 分别采用prim算法与kruskal算法构造最小生成树(第一次作业)

    分别采用prim算法与kruskal算法构造最小生成树 1 问题 举一个实例 画出采用Prim算法构造最小生成树的过程 并按实验报告模板编写算法 举一个实例 画出采用Kruskal算法构造最小生成树的过程 并按实验报告模板编写算法 有n个村
  • axios get请求特殊字符编码问题

    这几天在写一个项目 然后就遇到了请求的编码问题 然后在度娘上搜到了答案 请求拦截器配置处理 this axiosInstance interceptors request use config AxiosRequestConfig gt c
  • go换源国内并根据mod文件下载依赖

    go env w GO111MODULE on go env w GOPROXY https goproxy cn direct 根据mod文件下载依赖 此命令需要在go mod同级目录下执行 go mod tidy
  • 做EEG频谱分析,看这一篇文章就够了!

    所谓频谱分析 又称为功率谱分析或者功率谱密度 Power Spectral Density PSD 分析 实际就是通过一定方法求解信号的功率power随着频率变化曲线 笔者在这里对目前常用的频谱分析方法做一个总结 并重点介绍目前EEG分析中