2023-08-10发表2023-12-06更新生信

因素之间的相关性分析

最近在和师姐一起开展的品质性状可塑性的论文课题中, 希望讲述有关环境和杂交组合的相互作用; 但是由于实验没有重复, 因此没办法按照传统的做法对G, E, 和G×E进行拟合计算, 自然也不能很漂亮地画出图证明”xxx性状中, 基因型作用占比56%, 环境因素占比23%, 剩余部分是互作bulabula”. 因此需要一种分析, 表明”我们的表型是和环境有关系滴”. 于是, 就想到了对其相关性进行检验的方法.

相关性分析的适用条件

显而易见: 两个连续变量, 并且可能存在线性关系.

其实还有很多细节, 比如: 分析的两个变量需要一一配对并且来自同一个个体(例如, 研究小学生玩手机时长和成绩的关系, 那就是每个小孩对应一个时长+一个成绩), 群体中不存在异常值等等.

如果两个变量都符合正态分布, 那使用Pearson相关系数进行检验; 如果样本采样结果不符合正态分布, 则使用Spearman相关系数.

Pearson相关系数

Pearson相关系数的定义是:两个变量之间的协方差与二者标准差的商.

针对总体的Pearson相关系数用字母$\rho$表示:
$$
\rho_{X,Y} = \frac{\mbox{cov}(X,Y)}{\sigma_X\sigma_Y}=\frac{E[(X-\mu_X)(Y-\mu_Y)]}{\sigma_X\sigma_Y},
$$
其中协方差公式里面的$\mu_X$是指样本X的加权平均, 因此也只有总体才能用上协方差cov().

但更多的, 使用的是样本的Pearson相关系数, 用字母$r$表示:
$$
r=\frac{\sum_i(X_i-\bar X)(Y_i-\bar Y)}{\sqrt{\sum_i(X_i-\bar X)^2\sum_i(Y_i-\bar Y)^2}},
$$
由于是抽取样本, 就不存在”加权平均”这一说法, 即使用正常抽样所得的数据的平均值计算即可.

算出来的Pearson $r$是一个[-1, 1]范围的值, 代表从”完全线性负相关”到”完全线性正相关”.

但是也可以从例子看出, 一些复杂的非线性关系无法使用Pearson进行刻画.

Pearson相关性也有p值置信度, 但是一般python都能直接导出, 原理就不赘述了.

Spearman相关系数

Spearman相关系数的作用和Pearson类似, 是将数据的秩次代替数据本身. 也就是说, 不用样本值本身, 而是使用样本值在总体中的排序代替其数值.

Spearman vs Pearson

如图, 只要排序不变, 也就是后一个数据点的X坐标与Y坐标都一定严格比前一个点大, Spearman相关性就会认定为1.

这样的变换使得Spearman具有更强大的适用性, 体现在:

允许某一变量是分立的等级变量;
不服从正态分布(分布类型未知).

如何检验样本数据是不是服从正态分布?

有的建议画出Q-Q图, 判断分位点与数据分布的值的大小, 通过判断是否偏离中间线.

另外可以使用Shapiro-Wilk正态分布检验一组数据是否服从正态分布:

from scipy import stats
data = [0.873, 2.817, 0.121, -0.945, -0.055, -1.436, 0.360, -1.478, -1.637, -1.869]
stat, p = stats.shapiro(data) # no need for sorted
# if p<0.05, then may NOT normal distribution

Python代码实现

先看参数定义:

def pearsonr(x: Any, y: Any) -> Any
def spearmanr(a: int,
              b: int = None,
              axis: Optional[int] = 0,
              nan_policy: Optional[str] = 'propagate',
              alternative: Optional[str] = 'two-sided') -> Any

Pearson函数就是传入两个list就可以, 而Spearman有更多参数; nan_policy是指Spearman遇到NaN值的时候应当如何处理, 设定为'omit'即可使其忽略NaN值继续计算.

两者函数返回值都是tuple. 下面放两个我分别使用的例子:

from scipy import stats
#Pearson
psr, p_value = stats.pearsonr(list(w_df[pheno]), list(w_df[weather]))
#Spearman
psr, p_value = stats.spearmanr(list(p_df[pheno]), list(p_df[p2]), nan_policy='omit')

就到这里, 吃饭!

#相关性数理统计 scipy

评论