数模5：相关系数

Tokirin.

已于 2023-08-30 16:08:40 修改

阅读量273

点赞数

文章标签：人工智能算法机器学习

于 2023-08-30 16:06:56 首次发布

本文链接：https://blog.csdn.net/TOK__/article/details/132456412

版权

Jarque-Bera（大样本 n>30）

Shapiro-wilk检验（小样本3≤n≤50）

1.简介

相关系数可用来衡量两个变量之间的相关性的大小，根据数据满足的不同条件，我们要选择不同的相关系数进行计算和分析（建模论文中最容易用错的方法）。

2.模型

（1）Pearson相关系数

A.总体Pearson相关系数

两组数据 $X:\{X_{1},X_{2},\cdots,X_{n}\}$ 和 $Y:\{Y_{1},Y_{2},\cdots,Y_{n}\}$ 是总体数据：

总体均值： $E\left(X\right)=\frac{\sum_{i=1}^{n}X_{i}}{n},E\left(Y\right)=\frac{\sum_{i=1}^{n}Y_{i}}{n}$

总体协方差： $\operatorname{Cov}(X,Y)=\frac{\sum\limits_{i=1}^n\left(X_i-E\left(X\right)\right)\left(Y_i-E\left(Y\right)\right)}n$

注：协方差的大小和两个变量的量纲有关，因此不适合做比较。

总体Pearson相关系数： $\color{red}{\rho_{XY}=\frac{\operatorname{Cov}(X,Y)}{\sigma_X\sigma_Y}\color{black}=\frac{\sum_{i=1}^n\frac{(X_i-E(X))}{\sigma_X}\frac{(Y_i-E(Y))}{\sigma_Y}}n}$

注：Pearson相关系数剔除变量量纲影响，即将X和Y标准化后的协方差。

B.样本Pearson相关系数

两组数据 $X:\{X_{1},X_{2},\cdots,X_{n}\}$ 和 $Y:\{Y_{1},Y_{2},\cdots,Y_{n}\}$ 是样本数据（一般查到的数据均是）：

样本均值： $\overline{X}=\frac{\sum_{i=1}^nX_i}{n},\overline{Y}=\frac{\sum_{i=1}^nY_i}{n}$

总体协方差： $\mathrm{Cov}(X,Y)=\frac{\sum_{i=1}^n(X_i-\overline{X})(Y_i-\overline{Y})}{n-1}$

样本Pearson相关系数： $\color{red}{r_{XY}=\frac{\operatorname{Cov}\left(X,Y\right)}{S_XS_Y}}$

其中，标准差 $S_X=\sqrt{\frac{\sum_{i=1}^n(X_i-\overline{X})^2}{n-1}}, \quad S_Y=\sqrt{\frac{\sum_{i=1}^n(Y_i-\overline{Y})^2}{n-1}}$

C.一些注意点

（1）非线性相关也会导致线性相关系数很大，例如图2。

（2）离群点对相关系数的影响很大，例如图3，去掉离群点后，相关系数为0.98。

（3）如果两个变量的相关系数很大也不能说明两者相关，例如图4，可能是受到了异常值的影响。

（4）相关系数计算结果为0，只能说不是线性相关，但说不定会有更复杂的相关关系（非线性相关），例如图5。

注：若两个变量是线性的关系，那么Pearson相关系数绝对值大的就是相关性强，反之亦然；在不确定两个变量关系的情况下，算出Pearson相关系数很大，并不能说明那两个变量线性相关，甚至不能认为他们相关，一定要画出散点图来看才行。

事实上，相关系数大小反应相关性必须在确定关系的情况下才有一定的意义；因此通常情况，我们更关注显著性（假设检验）。

D.Matlab常用的统计函数

E.假设检验

在假设检验之前，要先验证是否满足条件：

第一， 实验数据通常假设是成对的来自于正态分布的总体。因为我们在求Pearson相关性系数后，通常还会用t检验之类的方法来进行皮尔逊Pearson相关性系数检验，而t检验是基于数据呈正态分布的假设的。验证数据正态分布，具体见G部分内容。

第二，实验数据之间的差距不能太大。Pearson相关性系数受异常值的影响比较大。

第三：每组样本之间是独立抽样的。（构造t统计量时需要用到。)

1.提出假设：原设与备设，以是否显著与0为例：

$H_0:\rho=0; \quad H_1:\rho \ne 0;$

2.构造某一分布的统计量

(注：不含其他统计量，常见分布为标准正态分布，t分布， $\chi^{2}$ 分布，F分布)

$t=\rho\sqrt{\frac{n-2}{1-\rho^2}}$

t为自由度为n-2的t分布，证明参考文献[2]。

3.代入检验值

以相关系数0.5，样本数30为例：

$t^{*}=0.5\sqrt{\frac{30-2}{1-0.5^{2}}}=3.05505$

4.根据给定置信区间，确定拒绝域

一般是画出分布函数pdf，查表或电脑计算临界值，表参考[3]:

5.定论

如果落在拒绝域，则拒绝原设，接收备受；反之亦然。

F.p值检验

本质上和前面的步骤一致，但是能直观的数据显示是否满足条件，判断接收与拒绝。

$t^{*}=0.5\sqrt{\frac{30-2}{1-0.5^{2}}}=3.05505$ ，求出该值下的概率：

注：单双侧是否乘2

适当加入显著性标记，能够更直观。

G.验证数据正态分布

Jarque-Bera（大样本 n>30）

对于随机变量 $\{X_i\}$ ，偏度 $S$ ，峰值 $K$ ，构造统计量JB：

$JB=\frac n6[S^2+\frac{(K-3)^2}4]$

若 $\{X_i\}$ 正态分布，则大样本 $JB{\sim}\chi^{2}(2)$

在Matlab中，偏度，峰度函数分别为 skewness，kurtosis。

JB检验函数：[h,p] = jbtest(x,alpha)，当输出h等于1时，表示拒绝原假设；h等于0则代表不能拒绝原假设。 alpha就是显著性水平，一般取0.05，此时置信水平为1‐0.05=0.95。

Shapiro-wilk检验（小样本3≤n≤50）

因为一般的样本数据不会这么小，所以用到不多。Matlab没有具体的函数，可以使用SPSS。

具体过程参考[4]。

Q-Q图

在统计学中，Q‐Q图（Q代表分位数Quantile）是一种通过比较两个概率分布的分位数对这两个概率分布进行比较的概率图方法。首先选定分位数的对应概率区间集合，在此概率区间上，点(x,y)对应于第一个分布的一个分位数x和第二个分布在和x相同概率区间上相同的分位数。这里，我们选择正态分布和要检验的随机变量，并对其做出QQ图，可想而知，如果要检验的随机变量是正态分布，那么QQ图就是一条直线。要利用Q‐Q图鉴别样本数据是否近似于正态分布,只需看Q‐Q图上的点是否近似地在一条直线附近。（要求数据量非常大）

Matlab中函数 qqplot（data）。

（2）Spearman相关系数

A.Definition

两组数据X,YSpearman相关系数：

$r_s=1-\frac{6\sum_{i=1}^nd_i^2}{n\left(n^2-1\right)}$

其中， $d_i$ 为 $X_i,Y_i$ 等级差，具体实例如下下：

可得 $r_s=0.875$

注：可证明 $r_s$ 在-1与1之间。若数值相等，去所在位置算术平均。

Another definition

等级之间的Person系数：

Matlab中 corr函数可求Spearman相关系数：

（1）corr(X , Y , 'type' , 'Spearman')——X和Y必须是列向量

（2）corr(X , 'type' , 'Spearman') ——计算X矩阵各列之间的Spearman相关系数

B.假设检验

小样本（n<=30）

可以直接选择查表法：

见文件：

大样本

构造统计量： $r_s\sqrt{n-1}{\sim}N(0,1)$

检验： $H_0{:}r_s{=}0\text{,}H_1{:}r_s{\neq}0$

同理，也可以使用p值法，直接对比：
若p值大于0.05,因此我们无法拒绝原假设。（和0没有显著的差异）

注：Matlab：% 直接给出相关系数和p值 [R,P]=corr(Test, 'type' , 'Spearman')

（3）Comparison

1.连续数据，正态分布，线性关系，用Pearson相关系数是最恰当，当然用 Spearman相关系数也可以，就是效率没有pearson相关系数高。

2.上述任一条件不满足，就用Spearman相关系数，不能用Pearson相关系数。

3.两个定序数据之间也用spearman相关系数，不能用pearson相关系数。

(定序数据是指仅仅反映观测对象等级、顺序关系的数据，是由定序尺度计量形成的，表现为类别，可以进行排序，属于品质数据。例如：优、良、差；我们可以用1表示差、2表示良、3表示优，但请注意，用2除以1得出的2并不代表任何含义。定序数据最重要的意义代表了一组数据中的某种逻辑顺序。)

注：Spearman相关系数的适用条件比Pearson相关系数要广，只要数据满足单调关系（例如线性函数、指数函数、对数函数等）就能够使用。

3.代码

clear;clc;
% % 颜色定义
% % C = TheColor('sci',813);

%% Load data file
% If the filename contains spaces, use single quotes around it
load 'physical fitness test.mat';


%% Descriptive Statistics
MIN = min(Test);           % Minimum value for each column
MAX = max(Test);           % Maximum value for each column
MEAN = mean(Test);         % Mean for each column
MEDIAN = median(Test);     % Median for each column
SKEWNESS = skewness(Test); % Skewness for each column
KURTOSIS = kurtosis(Test); % Kurtosis for each column
STD = std(Test);           % Standard deviation for each column
RESULT = [MIN; MAX; MEAN; MEDIAN; SKEWNESS; KURTOSIS; STD]; % Store statistics in a matrix


%% Prior to calculating Pearson's correlation coefficient, it's important to visualize scatter plots to observe linear relationships between variable pairs.
% % Generate random m*n data matrix (using a 4x3 matrix as an example)
% m = 4;
% n = 3;
% data = rand(m, n);

% Create a scatter matrix plot using "parula" color map
[S,AX,BigAx,H,HAx] = plotmatrix(Test);
colormap('parula'); % Set colormap to "parula"


% Detail
% Ajust Axis
set(AX, 'linewidth',1, ...
        'FontName', 'Arial', ...
        'FontSize', 5)

set(HAx,'linewidth',1, ...
        'FontName', 'Arial', ...
        'FontSize', 5)

% Add label
n = size(Test, 2);
variable_names = {'Variable 1', 'Variable 2', 'Variable 3','Variable 4', 'Variable 5', 'Variable 6'};
for i = 1:n
    for j = 1:n
        ax = AX(i, j);
 
        if i == n
            xlabel(ax, variable_names{j}, 'FontWeight', 'bold');
        end
        
        if j == 1
            ylabel(ax, variable_names{i}, 'FontWeight', 'bold');
        end
    end
end

% Background color
set(gcf,'Color',[1 1 1])

% Set the title and axis labels
title('Scatter Matrix of Matrix Variables');


% Save the image to a specified path
save_path = 'scatter_matrix.png'; % Replace with your desired path and filename
saveas(gcf, save_path);

% Display the image
imshow(save_path);


%% Calculate correlation coefficients between columns
% Tools like SPSS can help create scatter plots.
R = corrcoef(Test); % Correlation coefficient matrix

% Draw Heatmap
xname = {'Variable 1', 'Variable 2', 'Variable 3','Variable 4', 'Variable 5', 'Variable 6'};
yname = {'Variable 1', 'Variable 2', 'Variable 3','Variable 4', 'Variable 5', 'Variable 6'};
h = heatmap(xname, yname, R);
h.CellLabelFormat = '%0.2f';
% h.CellLabelColor = 'none';
% colormap(gca, 'parula');


%% Hypothesis Testing
x = -4:0.1:4;
y = tpdf(x, 28);   % Compute probability density of t-distribution with 28 degrees of freedom
figure(1);
plot(x, y, '-');
grid on;
hold on;
critical_value = tinv(0.975, 28); % Calculate critical value (2.0484)
plot([-critical_value, -critical_value], [0, tpdf(-critical_value, 28)], 'r--');
plot([critical_value, critical_value], [0, tpdf(critical_value, 28)], 'r--');

% Calculate p-value
x = -4:0.1:4;
y = tpdf(x, 28);
figure(2)
plot(x, y, '-')
grid on
hold on
plot([-3.055, -3.055], [0, tpdf(-3.055, 28)], 'r-');
plot([3.055, 3.055], [0, tpdf(3.055, 28)], 'r-');
disp('The p-value corresponding to this test statistic is:\n');
disp((1 - tcdf(3.055, 28)) * 2); % Multiply by 2 for two-tailed test

% Calculate correlation coefficients and p-values between columns
[R, P] = corrcoef(Test);
% Mark significance levels in an Excel sheet
stars_3 = P < 0.01; % Mark positions with 3 stars
stars_2 = (P < 0.05) .* (P > 0.01); % Mark positions with 2 stars
stars_1 = (P < 0.1) .* (P > 0.05); % Mark positions with 1 star


%% Normal Distribution Test
% JB Test
n_c = size(Test, 2); % Number of columns
H = zeros(1, n_c); % Initialize matrices for efficiency
P = zeros(1, n_c);
for i = 1:n_c
    [h, p] = jbtest(Test(:, i), 0.05);
    H(i) = h;
    P(i) = p;
end
disp(H);
disp(P);

% Q-Q plot
qqplot(Test(:, 1));



%% Spearman's Rank Correlation
X = [3 8 4 7 2]'; % Ensure column vectors
Y = [5 10 9 10 6]';
% First method of calculation
correlation_coefficient = 1 - 6 * (1 + 0.25 + 0.25 + 1) / 5 / 24;
% Second method of calculation
coeff = corr(X, Y, 'type', 'Spearman');
% Equivalent to:
RX = [2 5 3 4 1];
RY = [1 4.5 3 4.5 2];
R = corrcoef(RX, RY);
% Calculate Spearman correlation for matrix columns
R = corr(Test, 'type', 'Spearman');

% Large Sample Hypothesis Testing
disp(sqrt(590) * 0.0301); % Calculate test statistic
disp((1 - normcdf(0.7311)) * 2); % Calculate p-value for two-tailed test

% Calculate Spearman's Rank Correlation coefficients and p-values between columns
[R, P] = corr(Test, 'type', 'Spearman');