爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 119661|回复: 377

[原创]Matlab气候突变检测程序合集(更新了)

  [复制链接]

新浪微博达人勋

发表于 2017-1-13 13:07:57 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
本帖最后由 sfhstcn2 于 2017-10-19 07:29 编辑

    mutation_test这是一包代码,是基于《现代气候统计诊断与预测技术》(魏凤英)书中的气候突变检测一章的7种检测方法改写+编写的用于Matlab计算的代码。里面一共有8个function,都可直接调用,正好对应了书中的8种检验方法:
1 Movet  滑动t检验
2 Cramert  克拉默法
3 Yamamoto  山本法
4 MannKendall  M-K法
5 Pettitt  佩蒂特法
6 LePage  勒帕热法
7 BGSA  B-G分割算法
8 wvlt_trans 小波变换法
每个函数的介绍和用法都有精心撰写的详细注释,help 函数名 即可看到,不过由于想让程序更国(zhuang)际(bi)一点,因而注释都写成了英文。导入原序列后,每个程序都能算出用于检验的序列,还能直接画出检验图(当然如果不想画图,只要算出来的序列,可以把程序里面figure前面注释掉的return注释取消)。这里附一张B-G分割算法画的图。

    用法:下载->解压->放在某个路径下->set path->可以直接调用了。    这里补了一个我测试用的序列,比如对于BGSA,导入数据后输入:[ FLAGS ] = BGSA( year,zscore(series),.95,20 );
    当然,突变检验远不止这些方法,如果大家还有什么其它好的检验方法可以提在下面,附上链接什么的,有程序当然更好,我争取整合进来,壮大这一包代码。


    mutation_test_V2.0修改了之前的一些错误;增进了函数调用的灵活性(可以不必对所有输入的变量赋值,比如对于BGSA,直接[ FLAGS ] = BGSA( series );都行,程序会给这些参数赋上缺省值default。当然,具体设置视情况而定);增加了小波变换突变点检测的程序。小波变换的函数是用的很经典的http://paos.colorado.edu/research/wavelets/上的。
    接下来再来个高级一点的,mutation_test_menu这包增加了一个可以做突变检验的下拉菜单(如图2所示)。设置好路径后,命令窗口输入mt_menu; 就出来了。直接点击操作就可以,很简单。弹出的对话框用于设置各种参数。(注意:在计算画图之前一定要将序列变量名赋成y,对应的时间变量名赋成x!
    另外,大家是否发现右上角还多了个Image Saving菜单?那是保存图片用的,详细见贴:http://bbs.06climate.com/forum.php?mod=viewthread&tid=56495








BGSA

BGSA
截图.jpg

mutation_test.zip

6 KB, 下载次数: 476, 下载积分: 金钱 -5

序列.xlsx

11.27 KB, 下载次数: 325, 下载积分: 金钱 -5

mutation_test_V2.0.zip

15.97 KB, 下载次数: 547, 下载积分: 金钱 -5

mutation_test_menu.zip

18.46 KB, 下载次数: 314, 下载积分: 金钱 -5

售价: 1 贡献  [记录]

评分

参与人数 13金钱 +80 贡献 +11 收起 理由
嘴吧大肚皮小 + 5 很给力!
xxxuue + 2 很给力!
huangbicheng + 20 + 2 赞一个!
漠小北 + 1 很给力!
wangzhength + 1
lishuoyang00 + 1 赞一个!
yibao3 + 1 很给力!
liuhj29 + 2 很给力!
茶叶蛋花 + 1 很给力!
青风海魅 + 1 很给力!
denghaijun + 10 + 2 很给力!
二爷名声在外 + 15 + 5 很给力!
雷子 + 20 + 2 赞一个!

查看全部评分

密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2018-8-4 21:47:59 | 显示全部楼层
sfhstcn2 发表于 2018-3-14 13:11
看了下,这个方法只寻找统计量“秩和”的最大值作为突变点,对于一些趋势变化不明显的序列确实是有可能出 ...

感谢楼主的无私奉献,有关pettitt的检验,发现用如下程序是可以的,但用您的程序是出现了如前两位所说的问题,即最后位是突变点而且通过了显著性检验,请明察,菜鸟小白,多拍砖 :)
https://www.mathworks.com/matlabcentral/fileexchange/60973-pettitt-change-point-test-for-univariate-time-series-data
密码修改失败请联系微信:mofangbao
回复 支持 3 反对 0

使用道具 举报

新浪微博达人勋

发表于 2017-4-7 11:06:42 | 显示全部楼层
楼主的pettitt检验最后一句写错了,应该是除以(N^3+N^2)
密码修改失败请联系微信:mofangbao
回复 支持 2 反对 0

使用道具 举报

新浪微博达人勋

发表于 2017-4-7 11:23:26 | 显示全部楼层
楼主的Pettitt test
t0 = find(s==kt0);


这句也写错了
应该这样:t0 = find(abs(s)==kt0);
密码修改失败请联系微信:mofangbao
回复 支持 2 反对 0

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2019-2-17 13:58:12 | 显示全部楼层
本帖最后由 sfhstcn2 于 2019-2-17 14:02 编辑
漠小北 发表于 2019-2-17 13:21
哇,谢谢前辈解答,那如果我是想找突变点,就是变化趋势方向发生转变(由增加变为减少或者由减少变为增加 ...

可以对序列做个低通滤波,然后做一阶差分(相当于求一阶导),找到一阶差分后序列的等于或者最接近于0的点(f'(x)=0极值点)
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2019-2-15 23:55:46 | 显示全部楼层
漠小北 发表于 2019-2-15 15:10
前辈加楼主:

         您好!

不敢当不敢当。是这样子,是关于“突变检验”和“突变点”含义理解的问题,“突变点”并不是“序列趋势变化点f'(x)=0甚至f''(x)=0”,而是指的“序列出现了某种明显均值或者方差变化的地方”并且这种变化通过了检验的显著水平,因而“突变点”可能不是一点,甚至可能是一段。Pettitt检验你给的这个序列结果显示是2005年,P=0.0042,显著的;你做的t检验结果中2005年这点突变也是显著的。当然t检验结果显示后面那一段都是显著的,而Pettitt检验只有2005年这一点,这是由于两者算法本身的差异,Pettitt检验只给出了检验所计算的统计量中最大值所在处,并不能说明其它点就不是显著突变点。
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

发表于 2018-8-10 21:55:58 | 显示全部楼层
zzcsunmoonstar 发表于 2018-7-25 21:04
构造的时间序列是一个很明显的先增后减突变,不知道为什么检测不出来?

M-K检验等属均值突变检验,建议采用滑动F检验或者改进的有序聚类分析等均方差突变检验方法,综合诊断突变点
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

发表于 2018-1-21 22:18:09 | 显示全部楼层
jinjin 发表于 2018-1-20 21:27
不好意思楼主 还想冒昧请教您  BG分割图中的阴影部分代表什么呢  还有Y轴代表统计量吗

不好意思楼主 自己重新看了程序,上面解释的很清楚  ,打扰了 ,y轴代表自己检测的原始数据,阴影代表突变密集区  
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

发表于 2017-2-9 14:00:03 | 显示全部楼层
我发过一个帖子,关于突变检测的,可惜是R语言的,搞得自己头大。附上链接
http://bbs.06climate.com/forum.php?mod=viewthread&tid=49384
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

发表于 2017-1-13 13:43:48 | 显示全部楼层
可以说非常牛逼了。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-13 13:48:03 | 显示全部楼层
很国际的  赞  拜服
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-13 15:39:33 | 显示全部楼层
好东西,感谢lz分享!
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-20 09:58:57 | 显示全部楼层
感谢楼主分享
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-23 02:36:41 | 显示全部楼层
感谢楼主分享
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-24 08:31:32 | 显示全部楼层
感谢楼主分享
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-1-24 20:39:40 | 显示全部楼层
非常牛逼了。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-2-11 10:27:56 | 显示全部楼层
赞!
密码修改失败请联系微信:mofangbao
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表