性状的pearson相关与表型相关的差异

pearson相关是指两个数量性状间的积差相关，其计算公式如下：
$${\rho }_{x,y}=\frac{Cov(x,y)}{{\sigma }_x\ast {\sigma }_y}$$
式中，$Cov(x,y)$是两性状数值间的协方差，${\sigma }_x$、${\sigma }_y$分别是两个性状的标准差。

实际中，可用下述方法计算：
$$r_{x,y}=\frac{{\sum }_{i=1}^n(x_i-\bar{x})\ast (y_i-\bar{y})}{{\sum }_{i=1}^n(x_i-\bar{x}{)}^2\ast {\sum }_{i=1}^n(y_i-\bar{y}{)}^2}$$

而遗传参数中的表型相关是指两个性状的表型协方差与其表型标准差积的比值，公式如下：
$$r_{p(x,y)}=\frac{Co{v}_{p(x,y)}}{{\sigma }_{px}\ast {\sigma }_{py}}$$
式中，$Co{v}_{p(x,y)}$是两性状间的表型协方差，${\sigma }_{px}$、${\sigma }_{py}$分别是两个性状的表型标准差。

因此，pearson相关与表型相关是不同的两个参数，对于固定的样本，指定性状间的pearson相关值是恒定的，但其表型相关会随着混合线性模型的不同而发生变化，即表型相关值是可以变化的。

下面举例说明：

// pearson correlation
Maize.df<-read.csv(file='MaizeRILs.csv')
cor(Maize.df[,c(6,7,9)])

指定3个性状间的pearson相关值如下：

// result
> cor(Maize.df[,c(6,7,9)])
            pollen    silking     height
pollen   1.0000000  0.9756044 -0.4193743
silking  0.9756044  1.0000000 -0.4335516
height  -0.4193743 -0.4335516  1.0000000

现在来看看，表型相关的情况，我们将指定3个模型，具体如下：

// 3 mixed model
m1<-echidna(cbind(height,pollen,silking)~Trait+Trait:Location,
            random=~us(Trait):RIL,
            residual=~units:us(Trait),
            mulT=TRUE,
            es0.file = 'MaizeRILs.es0')

m1a <- update(m1, random=~us(Trait):RIL
                         +us(Trait):Location:RIL)

m1b <- update(m1, random=~us(Trait):RIL
             +us(Trait):Location:RIL
             +us(Trait):Location:Rep)

表型相关的计算：

// result
> m1.pvc <- update(m1, vpredict=c('F covp1 us(Trait):RIL+units:us(Trait)',
+                                 'R rp1 covp1'))
> cat(m1.pvc$evp)
......
R rp1 covp1
     rp1 1   2  1 = covp1 15/SQR[covp1 14*covp1 16]=     0.4427        0.0757
     rp1 1   3  1 = covp1 17/SQR[covp1 14*covp1 19]=     0.3146        0.0782
     rp1 1   3  2 = covp1 18/SQR[covp1 16*covp1 19]=     0.7176        0.0380

同理可得另2个模型的表型相关值：

// An highlighted block
> cat(m1a.pvc$evp)
......
R rp1a covp1a
     rp1a 2  2  1 = covp1 21/SQR[covp1 20*covp1 22]=     0.4433        0.0761
     rp1a 2  3  1 = covp1 23/SQR[covp1 20*covp1 25]=     0.3151        0.0786
     rp1a 2  3  2 = covp1 24/SQR[covp1 22*covp1 25]=     0.7177        0.0381
> cat(m1b.pvc$evp)
......
R rp1b covp1b
     rp1b 2  2  1 = covp1 27/SQR[covp1 26*covp1 28]=     0.4715        0.0757
     rp1b 2  3  1 = covp1 29/SQR[covp1 26*covp1 31]=     0.3432        0.0787
     rp1b 2  3  2 = covp1 30/SQR[covp1 28*covp1 31]=     0.7079        0.0403

此时，大家就可以看出模型不同，表型相关值是不同的。
有意思的是，性状height与其它两个性状间的pearson相关是负相关，而其表型相关都是正的。
综上，pearson相关与表型相关不是一回事。