计算数据的经验分布函数与MATLAB作图

计算数据的经验分布函数与MATLAB作图

写在前面

因为某些原因,需要处理某些数据,比如说某项测量数据与理论值的误差,我们就需要检验误差是否符合正态分布。最直观的方法就是直接做出经验分布函数的图来进行观察。所以这里简单的写一下如何做出经验分布图。

第一步,做经验分布图

对数据按升序排列

首先我们将所有一共n个数据按照升序排列。之所以升序,是因为这样方便我们后续计算。

将区间[a, b]分成小段

a是所有数据的最小值,b是所有数据的最大值,我们将所有数据分成小区间,比如说将[a,b]分为

  1. [a,a+(b-a)/100]
  2. [a, a+2*(b-a)/100]
  3. ...
  4. [a, b]

注意,区间左值永远是a,然后计算每个区间内落入数据点的个数,这样就得到了一个类似于频数分布直方图的东西。我们记第i个区间内数据点的个数为f_i,令所有f_i除以n得到数据落入第i个区间的频率p_n,也就是

p_i = f_i/n

当n无穷大的时候,就可以认为p_i是数据在这个区间内的概率,并且有p_n=1,也就是具有归一性。

那么我们就直接做出区间右值与p_i的函数图象,就是这个数据的经验分布函数了。

写程序

假设我们现在的数据都在data里,是1xn大小的向量。

[~, n] = size(data);
%最大值与最小值
d_max = max(data);
d_min = min(data);
%划区间,step是小的区间长度
x = d_min:step:d_max;
len = length(x);
%提前给区间的频数f预留空间
func = zeros(1, len)
for i = 1:len
    %sum内的data<=(d_min+i*step)是逻辑判断语句
    %整个语句意思是找出data内小于等于d_min+i*step的总个数
    func(i) = sum(data<=(d_min+i*step))
end
plot(x, func/n, 'b-', 'LineWidth', 1)
title('经验分布函数')

第二步,做正态分布的图

要知道我们的经验函数虽然假设服从正态分布,所以我们要知道服从的是均值为多少,标准差为多少的正态分布,在这里我们就需要计算我们原始数据的均值mu与标准差sigma。

计算出mu和sigma后做出对应的正态分布图即可。

下面是程序

%均值
mu = sum(data)/n;
%标准差
sigma = sqrt(1/(n-1)*sum((data-mu).^2));
%正太分布数据点
func_norm = normcdf(x, mu, sigma);
plot(x, func_norm, 'k-', 'LineWidth', 1)

第三步,做置信带

实际上我也不是很清楚置信带具体的含义是什么,因为对于标准正态分布的95%置信带,意思就是样本点落在置信带内的概率为95%,但是如果是经验分布函数的话,画上一个置信带的意义就不是很明确了。

不过即便如此,这里还是画上比较好。

%f实际上就是经验分布函数,与上文中计算出来的func一模一样,x是横轴
%flo是置信下区间,fup是置信上区间
[f, x, flow, fup] = ecdf(data);
stairs(x, flow, 'r:', 'LineWidth', 1)
stairs(x, fup, 'r:', 'LineWidth', 1)

%写一下图例等
legend('经验分布函数', '正态分布函数', '置信下限', '置信上限')
xlabel('x')
ylabel('F(x)')

最后得出图像

如下图所示

可以看出来我们的经验分布函数和正态分布吻合的比较好。

感想

实际上我感觉置信带应该画正态分布的,而不是如果我们的经验分布函数所有点都在正态分布的置信带里面,则说明我们的数据服从正态分布,置信程度为95%。这样子还比较能够说得过去。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 159,015评论 4 362
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 67,262评论 1 292
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 108,727评论 0 243
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 43,986评论 0 205
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 52,363评论 3 287
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 40,610评论 1 219
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 31,871评论 2 312
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 30,582评论 0 198
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 34,297评论 1 242
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 30,551评论 2 246
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 32,053评论 1 260
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 28,385评论 2 253
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 33,035评论 3 236
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 26,079评论 0 8
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,841评论 0 195
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 35,648评论 2 274
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 35,550评论 2 270

推荐阅读更多精彩内容

  • 这个不错分享给大家,从扣上看到的,就转过来了 《电脑专业英语》 file [fail] n. 文件;v. 保存文...
    麦子先生R阅读 6,515评论 5 24
  • 水东方阅读 100评论 0 0
  • 我发现很多讲解编程的书籍,在前面都会详细地讲解相关的预备知识,而大多数读者却更希望马上进行实践。没错,人们总是对基...
    Allan要做活神仙阅读 647评论 0 0
  • 我很想做个小女人 不用自己努力拼搏 在累的时候会有个人摸摸我的头 跟我说不要太辛苦我来… 我很想做个小女人 不用自...
    Higanbana燕阅读 207评论 0 0
  • 小涛哥 2017.2.5 十几年前的事了 在一个很小的县城 人群中传来“嘭”的一声 吓的我哆嗦了片刻 尤其在不防备...
    宋大圣阅读 194评论 0 3