简单相关系数并不简单

写这篇文章，源于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)，不知是否正确，还请达人对此指正。