今天无意中发现Excel中除了Pearson，还有CorRel！二者究竟有何区别？Excel的Help含糊不清，没有发现有价值的提示信息。
google了一圈，终于在http://support.microsoft.com/kb/828129/zh-cn发现了权威的答复。
在>=2003版本中，Pearson与CorRel相同
在<2003版本中，Pearson与CorRel不同
在早于2003的Excel版本中，PEARSON可能会出现舍入错误，该缺陷在2003及更高版本的Excel中得到了改进。
CORREL总是通过这种改进的过程来实现。
2009-03-31
2009-03-28
简单相关系数并不简单
写这篇文章，源于pinggu.org上cathy_wxh关于“相关系数及其显著性水平”的求助贴，在与xmok77的讨论过程中，我发现自己对此仍有不少概念仍然不够清楚，需要澄清，因此萌发了写此文的念头。好，下面言归正传。
相关系数有很多种，如Pearson积差相关、Spearman等级相关、Kendall等级相关、Kendall和谐系数、Kappa一致性系数……，最常用的当然是Pearson积差相关系数，又称为简单相关系数，但它可不简单喔！
1、r的计算
首先复习一下有关知识：
总体参数协方差σxy=E{[X-E(X)][Y-E(Y)]}
总体参数相关系数ρxy=σxy/[sqrt(σxx)*sqrt(σyy)]
统计量样本协方差sxy=sum[(Xi-x)(Yi-y)]/(n-1)
统计量样本相关系数rxy=sxy/[sqrt(sxx)*sqrt(syy)]
当X、Y都服从正态分布时，ρxy常用rxy进行估计。
2、r的意义
相关关系不代表因果关系。
Pearson相关系数只反映线性相关程度，不反映曲性相关程度。
不少书上都谈到了相关系数的几何意义，但都没有谈透，我自己觉得Pearson相关系数反映了“散点向线段聚集”的程度。当然还有许多其它意义，我还想对此做进一步的探索。
3、r的分布
r究竟符合什么样的分布，我看到的教科书上都没有明确的答案，只是给出了两个检验方法，这在后面还会详细讲解。于是我试图自己推导，但高等数学早就忘得差不多了，实在无法进行。恰在此时，verycd.com上居然发布了Crystal Ball V11！简单的学习后我异想天开，利用Monte Carlo模拟对r的分布进行数值计算！其中遇到了不少困难，但都找到了解决办法，这不是本文的重点，以后我或许会写一篇关于Crystal Ball的介绍。下图是我的工作界面，相信对Crystall Ball熟悉的朋友一下就明白了。
这里简单解释一下大概的流程。
1)、对正态总体X抽样，得到容量为30的样本，将之作为Monte Carlo模拟的假设Assumption。
设定X01～X30共30个相互独立的随机变量，Xi~N(μx、σx^2)，这里μx、σx具体是多少并不重要。
2)、对正态总体Y抽样，得到容量为30的样本，将之作为Monte Carlo模拟的假设Assumption。
设定Y01～Y30共30个相互独立的随机变量，Yi~N(μy、σy^2)，这里μy、σy具体是多少并不重要。但注意对于每个Yi一定要设定它与Xi的相关系数为ρ。
3)、计算r，将之作为Crystal Ball的预测Forecast。
为了方便，可以同时计算r(n=2)、r(n=3)、……、r(n=30)
4)、Run the Model!
经过反复试验，我发现r的分布应该取决于总体参数ρ和抽样次数n。下图是进行100000次Monte Carlo模拟后，再用ImageMagick拼接后得到的结果：
ImageMagick网上有不少介绍，在此不再赘述，有兴趣的朋友google一下吧。
可惜的是Crystal Ball的绘图功能有限，只能定制X轴，无法定制Y轴，因此需要注意这些子图的X轴范围是一样的，都是-1～+1，但Y轴的范围都不太一样。但这无关大局。
从图中有几个有趣的地方：
1)、当n=2时，无论ρ是多少，r只能取±1，随着ρ由0→1，r取得1的概率逐步增加，取得-1的概率逐步减少。
2)、当ρ=0时：
当n=3，则r以很大的概率取得接近±1的值，取得0的概率反而很小。
当n=4，则r取得[-1,+1]之间任何数的概率都相等。
当n由5→+∞，r取得0的概率逐步增加，r也越来越象正态分布
3)、当n<10时，r的分布距离ρ实在太散，即如果用计算得到的r来估计ρ的话风险太大。>10时，r的分布距离ρ较为集中，即如果用计算得到的r来估计ρ的话风险较小。
5)、在同样的n下，当ρ由0→1，相应的r的分布也会越来越集中。
6)、当ρ=0时，r的分布左右对称；
当ρ>0时，r的分布左偏；
当ρ<0时，r的分布右偏。 r="1，因此图中没有画出ρ=" n="400！（我坚决不愿做Pearson的助手:)）。计算结果从最乐观的角度看也只是证明r的分布实在离正态分布太远。">3时就可以较好的接近正态分布！（看样子我也做不了Fisher的助手:)）。
3、r的检验
从前面可以看出，r似乎不遵从常见的概率分布，其偏离期望值ρ的程度很高，这就很好的解释了为什么教材书上格外强调对r需要进行显著性检验。
回想一下，当用样本平均值、样本方差这些统计量去估计总体均值、总体方差这些总体参数的时候可没有特别强调进行显著性检验。理由很简单：
如果总体X~N(μ、σ^2/n)，则样本平均值~N(μ、σ^2/n)；
如果总体X~N(μ、σ^2/n)，则样本方差有(n-1)s^2/σ^2~χ^2(n)；
这些分布的偏差都不是太大，因此常常被忽略，当然严格意义上讲这是不正确的，对于每个估计，都应该进行显著性检验。
那么，如何对其进行检验？
1)、t分布法
有人（我没有查到到底是哪位大牛，希望知道的朋友告诉我）发现当ρ=0时，r/sqrt[(1-r^2)/(n-2)]服从t(n-2)分布。
对此同样利用Crystal Ball进行100000次Monte Carlo模拟，再利用ImageMagick拼接，得到下图：
从直觉看似乎的确符合t的特点。
2)、z法/f法
Fisher发现若令z=ln(sqrt((1+x)/(1-x)))，则z(r)~N(z(ρ),1/n-3)。即[z(r)-z(ρ)]/sqrt[1/(n-3)]~N(0,1)。Fisher将之称为z分布，后人为表彰Fisher，将之称为Fisher分布！
对此同样利用Crystal Ball进行100000次Monte Carlo模拟，从直觉看，似乎的确如此。大家对正态图实在太熟悉了，这里就不再贴拼接图了。:)好好，我坦白，其实是我想偷懒还不行嘛。
这里再多说两句，绝大多数教材上都说，对于Pearson相关系数：一般先利用t方法进行“线性相关显著性检验”，如果相关不显著，则不管其相关系数r多大，都认为其不可接受，即认为ρ=0。
我想，应该也可以利用Fisher方法对其进行检验？即利用[z(r)-z(0)]/sqrt[1/(n-3)]~N(0,1)，不知是否正确，还请达人对此指正。
相关系数有很多种，如Pearson积差相关、Spearman等级相关、Kendall等级相关、Kendall和谐系数、Kappa一致性系数……，最常用的当然是Pearson积差相关系数，又称为简单相关系数，但它可不简单喔！
1、r的计算
首先复习一下有关知识：
总体参数协方差σxy=E{[X-E(X)][Y-E(Y)]}
总体参数相关系数ρxy=σxy/[sqrt(σxx)*sqrt(σyy)]
统计量样本协方差sxy=sum[(Xi-x)(Yi-y)]/(n-1)
统计量样本相关系数rxy=sxy/[sqrt(sxx)*sqrt(syy)]
当X、Y都服从正态分布时，ρxy常用rxy进行估计。
2、r的意义
相关关系不代表因果关系。
Pearson相关系数只反映线性相关程度，不反映曲性相关程度。
不少书上都谈到了相关系数的几何意义，但都没有谈透，我自己觉得Pearson相关系数反映了“散点向线段聚集”的程度。当然还有许多其它意义，我还想对此做进一步的探索。
3、r的分布
r究竟符合什么样的分布，我看到的教科书上都没有明确的答案，只是给出了两个检验方法，这在后面还会详细讲解。于是我试图自己推导，但高等数学早就忘得差不多了，实在无法进行。恰在此时，verycd.com上居然发布了Crystal Ball V11！简单的学习后我异想天开，利用Monte Carlo模拟对r的分布进行数值计算！其中遇到了不少困难，但都找到了解决办法，这不是本文的重点，以后我或许会写一篇关于Crystal Ball的介绍。下图是我的工作界面，相信对Crystall Ball熟悉的朋友一下就明白了。
这里简单解释一下大概的流程。
1)、对正态总体X抽样，得到容量为30的样本，将之作为Monte Carlo模拟的假设Assumption。
设定X01～X30共30个相互独立的随机变量，Xi~N(μx、σx^2)，这里μx、σx具体是多少并不重要。
2)、对正态总体Y抽样，得到容量为30的样本，将之作为Monte Carlo模拟的假设Assumption。
设定Y01～Y30共30个相互独立的随机变量，Yi~N(μy、σy^2)，这里μy、σy具体是多少并不重要。但注意对于每个Yi一定要设定它与Xi的相关系数为ρ。
3)、计算r，将之作为Crystal Ball的预测Forecast。
为了方便，可以同时计算r(n=2)、r(n=3)、……、r(n=30)
4)、Run the Model!
经过反复试验，我发现r的分布应该取决于总体参数ρ和抽样次数n。下图是进行100000次Monte Carlo模拟后，再用ImageMagick拼接后得到的结果：
ImageMagick网上有不少介绍，在此不再赘述，有兴趣的朋友google一下吧。
可惜的是Crystal Ball的绘图功能有限，只能定制X轴，无法定制Y轴，因此需要注意这些子图的X轴范围是一样的，都是-1～+1，但Y轴的范围都不太一样。但这无关大局。
从图中有几个有趣的地方：
1)、当n=2时，无论ρ是多少，r只能取±1，随着ρ由0→1，r取得1的概率逐步增加，取得-1的概率逐步减少。
2)、当ρ=0时：
当n=3，则r以很大的概率取得接近±1的值，取得0的概率反而很小。
当n=4，则r取得[-1,+1]之间任何数的概率都相等。
当n由5→+∞，r取得0的概率逐步增加，r也越来越象正态分布
3)、当n<10时，r的分布距离ρ实在太散，即如果用计算得到的r来估计ρ的话风险太大。>10时，r的分布距离ρ较为集中，即如果用计算得到的r来估计ρ的话风险较小。
5)、在同样的n下，当ρ由0→1，相应的r的分布也会越来越集中。
6)、当ρ=0时，r的分布左右对称；
当ρ>0时，r的分布左偏；
当ρ<0时，r的分布右偏。 r="1，因此图中没有画出ρ=" n="400！（我坚决不愿做Pearson的助手:)）。计算结果从最乐观的角度看也只是证明r的分布实在离正态分布太远。">3时就可以较好的接近正态分布！（看样子我也做不了Fisher的助手:)）。
3、r的检验
从前面可以看出，r似乎不遵从常见的概率分布，其偏离期望值ρ的程度很高，这就很好的解释了为什么教材书上格外强调对r需要进行显著性检验。
回想一下，当用样本平均值、样本方差这些统计量去估计总体均值、总体方差这些总体参数的时候可没有特别强调进行显著性检验。理由很简单：
如果总体X~N(μ、σ^2/n)，则样本平均值~N(μ、σ^2/n)；
如果总体X~N(μ、σ^2/n)，则样本方差有(n-1)s^2/σ^2~χ^2(n)；
这些分布的偏差都不是太大，因此常常被忽略，当然严格意义上讲这是不正确的，对于每个估计，都应该进行显著性检验。
那么，如何对其进行检验？
1)、t分布法
有人（我没有查到到底是哪位大牛，希望知道的朋友告诉我）发现当ρ=0时，r/sqrt[(1-r^2)/(n-2)]服从t(n-2)分布。
对此同样利用Crystal Ball进行100000次Monte Carlo模拟，再利用ImageMagick拼接，得到下图：
从直觉看似乎的确符合t的特点。
2)、z法/f法
Fisher发现若令z=ln(sqrt((1+x)/(1-x)))，则z(r)~N(z(ρ),1/n-3)。即[z(r)-z(ρ)]/sqrt[1/(n-3)]~N(0,1)。Fisher将之称为z分布，后人为表彰Fisher，将之称为Fisher分布！
对此同样利用Crystal Ball进行100000次Monte Carlo模拟，从直觉看，似乎的确如此。大家对正态图实在太熟悉了，这里就不再贴拼接图了。:)好好，我坦白，其实是我想偷懒还不行嘛。
这里再多说两句，绝大多数教材上都说，对于Pearson相关系数：一般先利用t方法进行“线性相关显著性检验”，如果相关不显著，则不管其相关系数r多大，都认为其不可接受，即认为ρ=0。
我想，应该也可以利用Fisher方法对其进行检验？即利用[z(r)-z(0)]/sqrt[1/(n-3)]~N(0,1)，不知是否正确，还请达人对此指正。
2009-03-22
用Excel生成二维函数曲面图
今天看到有人询问如何用Excel生成二维函数曲面图，或者用什么软件生成
用Octave/Matlab、Mathematica、Maple这些数学软件当然可以方便的画出二维概率密度函数图。
用Origin、Sigmlot这些半专业的科学图形软件当然也可以方便的画出二维概率密度函数图。
绝大多数人认为Excel无法生成二维函数曲面图，其实Excel当然也可以，只是繁琐一些，效果差一些而已。
步骤如下，
1、必须有一系列(x, y)向量值
2、将之code到区间上，比如[0.0～0.1)→0，[0.1～0.2)→1……[0.9～1.0)→9
3、利用透视表得到一张二维表格
4、根据这张二维表格，插入“曲面图”
在Excel2007中试验成功，效果可以见下图：
2009-03-18
如何提出Ho真是一门学问
今天上cos.name看到2008年统计学考研真题第四题“食品厂家说：净含量是每袋不低于250g。但有消费者向消协反映不是250g，消协据此要求厂家自检，同时消协也从中随机抽取20袋检验”
（1）如果厂家自己检验，你认为提出什么样的原假设和备则假设？并说明理由。
（2）如果从消费者利益出发，你认为应该提出什么样的原假设和备则假设？并说明理由。
我想：拒绝Ho是不容易的，因为必须p＜α。所以厂家基于自己的利益应该提出：
Ho：weight ≥ 250g；
Ha：weight ＜ 250g。
这样的话，即使实际测量得到的Weight＜250，如果短斤少两不是很厉害的话，也无法拒绝Ho。
当然，消费者基于自己的利益就应该提出：
Ho：weight ＜ 250g；
Ha：weight ≥ 250g；
但郑冰认为：
假设检验的主要目的是为了拒绝,所以厂家要说明自己的产品合格得拒绝原假设
因此应该提出：
H0：weight ＜ 250g；
谢益辉认为：
1、厂家会说，不拒绝又不表示接受，因此消费者如果提出H0：weight ＜ 250g则力度不够。
2、从是否容易拒绝零假设的角度来看，厂家提出weight ≥ 250g的确占点便宜。
（1）如果厂家自己检验，你认为提出什么样的原假设和备则假设？并说明理由。
（2）如果从消费者利益出发，你认为应该提出什么样的原假设和备则假设？并说明理由。
我想：拒绝Ho是不容易的，因为必须p＜α。所以厂家基于自己的利益应该提出：
Ho：weight ≥ 250g；
Ha：weight ＜ 250g。
这样的话，即使实际测量得到的Weight＜250，如果短斤少两不是很厉害的话，也无法拒绝Ho。
当然，消费者基于自己的利益就应该提出：
Ho：weight ＜ 250g；
Ha：weight ≥ 250g；
但郑冰认为：
假设检验的主要目的是为了拒绝,所以厂家要说明自己的产品合格得拒绝原假设
因此应该提出：
H0：weight ＜ 250g；
谢益辉认为：
1、厂家会说，不拒绝又不表示接受，因此消费者如果提出H0：weight ＜ 250g则力度不够。
2、从是否容易拒绝零假设的角度来看，厂家提出weight ≥ 250g的确占点便宜。
仔细想想应该这样说：
食品厂家对于计量是否准确，心里是有数的。
如果确实不短斤少两，食品厂家则按郑冰的做法，提出以下假设：
Ho：weight = 250g；
Ha：weight ＞ 250g。
这样可以非常有力的消除消费者的怀疑。
如果确实少上一点，食品厂家则按我的做法，提出以下假设：
Ho：weight = 250g；
Ha：weight ＜ 250g。
如果少得比较厉害，食品厂家还是及早承认错误吧。
不知这样解决如何？
MindManager 2008终于可以将.mmap直接转换成.swf文件了
MindManager 2008终于可以将.mmap直接转换成.swf文件了！
搞笑的是，无论是Save As的菜单栏中还是Export的菜单栏都没有直接的选项，必须到了Save As的对话框中，在保存类型中才有.swf文件类型，真让人好笑。真可谓“千户万换始出来，犹抱琵琶半遮面”。
不管如何，终究是一件好事，可惜blogger.com还是不支持上传.swf
搞笑的是，无论是Save As的菜单栏中还是Export的菜单栏都没有直接的选项，必须到了Save As的对话框中，在保存类型中才有.swf文件类型，真让人好笑。真可谓“千户万换始出来，犹抱琵琶半遮面”。
不管如何，终究是一件好事，可惜blogger.com还是不支持上传.swf
Subscribe to:
Posts (Atom)