（數據科學學習手札19）R中基本統計分析技巧總結

時間 2019-11-21

標籤數據科學學習手札基本統計分析技巧總結简体版

原文原文鏈接

在獲取數據，而且完成數據的清洗以後，首要的事就是對整個數據集進行探索性的研究，這個過程當中會利用到各類描述性統計量和推斷性統計量來初探變量間和變量內部的基本關係，本篇筆者便基於R，對一些經常使用的數據探索方法進行總結：api

1.描述性統計量部分數組

1.1 計算描述性統計量的常規方法app

summary()
less

summary()函數提供了最小值、最大值、四分位數和數值型變量的均值，以及因子向量和邏輯型向量的頻數統計：ide

> #掛載鳶尾花數據
> data(iris)
> #計算鳶尾花各變量的基本描述統計量
> summary(iris)
  Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
 Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
 1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
 Median :5.800   Median :3.000   Median :4.350   Median :1.300  
 Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
 3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
 Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
       Species  
 setosa    :50  
 versicolor:50  
 virginica :50

利用sapply()經過自編函數進行偏度峯度的計算：函數

> #自定義偏度峯度計算函數
> skew_kurt <- function(x,na.omit=FALSE){
+   if(na.omit)
+     x <- x[!is.na(x)]
+   m <- mean(x)
+   n <- length(x)
+   skew <- sum((x-mean(x))^3/sd(x)^3)/n
+   kurt <- sum((x-m)^4/sd(x)^4)/n-3
+   return(c(skew=skew,kurtosis=kurt))
+ }
> #利用sapply將偏度峯度計算擴展到數據框型變量
> sapply(iris[,-5],skew_kurt)
         Sepal.Length Sepal.Width Petal.Length Petal.Width
skew        0.3086407   0.3126147   -0.2694109  -0.1009166
kurtosis   -0.6058125   0.1387047   -1.4168574  -1.3581792

Hmisc包中的describe()函數：spa

Hmisc包中的describe()函數可返回變量和觀測的數量、缺失值和惟一值的數目、平均值、分位數，以及五個最大的值和五個最小的值：code

> library(Hmisc)
> Hmisc::describe(iris)
iris 

 5  Variables      150  Observations
-------------------------------------------------------------------------------
Sepal.Length 
       n  missing distinct     Info     Mean      Gmd      .05      .10 
     150        0       35    0.998    5.843   0.9462    4.600    4.800 
     .25      .50      .75      .90      .95 
   5.100    5.800    6.400    6.900    7.255 

lowest : 4.3 4.4 4.5 4.6 4.7, highest: 7.3 7.4 7.6 7.7 7.9
-------------------------------------------------------------------------------
Sepal.Width 
       n  missing distinct     Info     Mean      Gmd      .05      .10 
     150        0       23    0.992    3.057   0.4872    2.345    2.500 
     .25      .50      .75      .90      .95 
   2.800    3.000    3.300    3.610    3.800 

lowest : 2.0 2.2 2.3 2.4 2.5, highest: 3.9 4.0 4.1 4.2 4.4
-------------------------------------------------------------------------------
Petal.Length 
       n  missing distinct     Info     Mean      Gmd      .05      .10 
     150        0       43    0.998    3.758    1.979     1.30     1.40 
     .25      .50      .75      .90      .95 
    1.60     4.35     5.10     5.80     6.10 

lowest : 1.0 1.1 1.2 1.3 1.4, highest: 6.3 6.4 6.6 6.7 6.9
-------------------------------------------------------------------------------
Petal.Width 
       n  missing distinct     Info     Mean      Gmd      .05      .10 
     150        0       22     0.99    1.199   0.8676      0.2      0.2 
     .25      .50      .75      .90      .95 
     0.3      1.3      1.8      2.2      2.3 

lowest : 0.1 0.2 0.3 0.4 0.5, highest: 2.1 2.2 2.3 2.4 2.5
-------------------------------------------------------------------------------
Species 
       n  missing distinct 
     150        0        3 
                                           
Value          setosa versicolor  virginica
Frequency          50         50         50
Proportion      0.333      0.333      0.333
-------------------------------------------------------------------------------

pastecs包中的stat.desc()函數計算描述性統計量：orm

stat.desc(x,basic=TRUE,desc=TRUE,norm=FALSE,p=0.95),其中x是一個數據框或時間序列（非數值型的部分只會返回空值），若basic=TRUE（默認值）,則計算x中全部變量的值、空值、缺失值的數量，以及最小值、最大值、值域，還有總和。若desc=TRUE（默認值），則計算中位數、平均數、平均數的標準偏差、平均數置信度爲95%的置信區間、方差、標準差以及變異係數。若norm=TRUE（非默認），則返回正態分佈統計量，包括偏度和峯度（以及他們的統計顯著程度）和Shapiro-Wilk（正態性檢驗）結果：blog

> library(pastecs)
> stat.desc(iris)
             Sepal.Length  Sepal.Width Petal.Length  Petal.Width Species
nbr.val      150.00000000 150.00000000  150.0000000 150.00000000      NA
nbr.null       0.00000000   0.00000000    0.0000000   0.00000000      NA
nbr.na         0.00000000   0.00000000    0.0000000   0.00000000      NA
min            4.30000000   2.00000000    1.0000000   0.10000000      NA
max            7.90000000   4.40000000    6.9000000   2.50000000      NA
range          3.60000000   2.40000000    5.9000000   2.40000000      NA
sum          876.50000000 458.60000000  563.7000000 179.90000000      NA
median         5.80000000   3.00000000    4.3500000   1.30000000      NA
mean           5.84333333   3.05733333    3.7580000   1.19933333      NA
SE.mean        0.06761132   0.03558833    0.1441360   0.06223645      NA
CI.mean.0.95   0.13360085   0.07032302    0.2848146   0.12298004      NA
var            0.68569351   0.18997942    3.1162779   0.58100626      NA
std.dev        0.82806613   0.43586628    1.7652982   0.76223767      NA
coef.var       0.14171126   0.14256420    0.4697441   0.63555114      NA

psych包中的describe()計算描述性統計量：

在psych包中，能夠經過describe()來計算非缺失值的數量、平均數、標準差、中位數、截尾均值、絕對中位數、最小值、最大值、值域、偏度、峯度和平均值的標準偏差（注意，此方法會將因子型變量當成數值型來計算）：

> library(psych)
> psych::describe(iris)
             vars   n mean   sd median trimmed  mad min max range  skew
Sepal.Length    1 150 5.84 0.83   5.80    5.81 1.04 4.3 7.9   3.6  0.31
Sepal.Width     2 150 3.06 0.44   3.00    3.04 0.44 2.0 4.4   2.4  0.31
Petal.Length    3 150 3.76 1.77   4.35    3.76 1.85 1.0 6.9   5.9 -0.27
Petal.Width     4 150 1.20 0.76   1.30    1.18 1.04 0.1 2.5   2.4 -0.10
Species*        5 150 2.00 0.82   2.00    2.00 1.48 1.0 3.0   2.0  0.00
             kurtosis   se
Sepal.Length    -0.61 0.07
Sepal.Width      0.14 0.04
Petal.Length    -1.42 0.14
Petal.Width     -1.36 0.06
Species*        -1.52 0.07

1.2 分組計算描述性統計量

在比較多組個體或觀測時，關注的焦點常常是各組的描述性統計信息，而不是樣本總體的描述性統計信息。一樣地，在R中完成這個任務有若干種方法：

利用aggregate()來對數據數組求描述性統計量：

aggregate(data,by,fun)，其中data爲待考察的數據所在的數據框，至少有兩列，其中一列爲分組依據的類別型數據，即by，fun爲對分組後數據採起的運算，且fun只能爲返回單值的簡單函數（很遺憾）：

> aggregate(iris[,1:4],by=list(Species=iris$Species),mean)
     Species Sepal.Length Sepal.Width Petal.Length Petal.Width
1     setosa        5.006       3.428        1.462       0.246
2 versicolor        5.936       2.770        4.260       1.326
3  virginica        6.588       2.974        5.552       2.026
> aggregate(iris[,1:4],by=list(Species=iris$Species),sd)
     Species Sepal.Length Sepal.Width Petal.Length Petal.Width
1     setosa    0.3524897   0.3790644    0.1736640   0.1053856
2 versicolor    0.5161711   0.3137983    0.4699110   0.1977527
3  virginica    0.6358796   0.3224966    0.5518947   0.2746501

利用by()函數來計算更復雜的分組數據描述性統計量：

by(data,INDICES,FUN)，其中data爲輸入的數據框，INDICES是一個因子型或幾個因子型變量組成的list，是分組依據，FUN是任意函數（包括自編函數）：

> by(iris[,-5],iris$Species,psych::describe)
iris$Species: setosa
             vars  n mean   sd median trimmed  mad min max range skew kurtosis
Sepal.Length    1 50 5.01 0.35    5.0    5.00 0.30 4.3 5.8   1.5 0.11    -0.45
Sepal.Width     2 50 3.43 0.38    3.4    3.42 0.37 2.3 4.4   2.1 0.04     0.60
Petal.Length    3 50 1.46 0.17    1.5    1.46 0.15 1.0 1.9   0.9 0.10     0.65
Petal.Width     4 50 0.25 0.11    0.2    0.24 0.00 0.1 0.6   0.5 1.18     1.26
               se
Sepal.Length 0.05
Sepal.Width  0.05
Petal.Length 0.02
Petal.Width  0.01
----------------------------------------------------------- 
iris$Species: versicolor
             vars  n mean   sd median trimmed  mad min max range  skew
Sepal.Length    1 50 5.94 0.52   5.90    5.94 0.52 4.9 7.0   2.1  0.10
Sepal.Width     2 50 2.77 0.31   2.80    2.78 0.30 2.0 3.4   1.4 -0.34
Petal.Length    3 50 4.26 0.47   4.35    4.29 0.52 3.0 5.1   2.1 -0.57
Petal.Width     4 50 1.33 0.20   1.30    1.32 0.22 1.0 1.8   0.8 -0.03
             kurtosis   se
Sepal.Length    -0.69 0.07
Sepal.Width     -0.55 0.04
Petal.Length    -0.19 0.07
Petal.Width     -0.59 0.03
----------------------------------------------------------- 
iris$Species: virginica
             vars  n mean   sd median trimmed  mad min max range  skew
Sepal.Length    1 50 6.59 0.64   6.50    6.57 0.59 4.9 7.9   3.0  0.11
Sepal.Width     2 50 2.97 0.32   3.00    2.96 0.30 2.2 3.8   1.6  0.34
Petal.Length    3 50 5.55 0.55   5.55    5.51 0.67 4.5 6.9   2.4  0.52
Petal.Width     4 50 2.03 0.27   2.00    2.03 0.30 1.4 2.5   1.1 -0.12
             kurtosis   se
Sepal.Length    -0.20 0.09
Sepal.Width      0.38 0.05
Petal.Length    -0.37 0.08
Petal.Width     -0.75 0.04
> 
> #前面定義過的偏度峯度自編函數做爲sapply的參數跟在FUN位置的後面
> by(iris[,-5],iris$Species,sapply,skew_kurt)
iris$Species: setosa
         Sepal.Length Sepal.Width Petal.Length Petal.Width
skew        0.1129778  0.03872946    0.1000954    1.179633
kurtosis   -0.4508724  0.59595073    0.6539303    1.258718
----------------------------------------------------------- 
iris$Species: versicolor
         Sepal.Length Sepal.Width Petal.Length Petal.Width
skew       0.09913926  -0.3413644   -0.5706024 -0.02933377
kurtosis  -0.69391378  -0.5493203   -0.1902555 -0.58731442
----------------------------------------------------------- 
iris$Species: virginica
         Sepal.Length Sepal.Width Petal.Length Petal.Width
skew        0.1110286   0.3442849    0.5169175  -0.1218119
kurtosis   -0.2032597   0.3803832   -0.3651161  -0.7539586

利用doBy中的summaryBy()來分組計算描述性統計量：

summaryBy(x1+x2+x3~by,data,FUN)，其中data爲目標數據框，x1,x2,x3,by均爲data中的變量，且by爲分組依據變量，FUN爲任意函數：

> library(doBy)
> data(mtcars)
> summaryBy(mpg+wt~cyl,data=mtcars,FUN=summary)
  cyl mpg.Min. mpg.1st Qu. mpg.Median mpg.Mean mpg.3rd Qu. mpg.Max. wt.Min.
1   4     21.4       22.80       26.0 26.66364       30.40     33.9   1.513
2   6     17.8       18.65       19.7 19.74286       21.00     21.4   2.620
3   8     10.4       14.40       15.2 15.10000       16.25     19.2   3.170
  wt.1st Qu. wt.Median  wt.Mean wt.3rd Qu. wt.Max.
1     1.8850     2.200 2.285727    2.62250   3.190
2     2.8225     3.215 3.117143    3.44000   3.460
3     3.5325     3.755 3.999214    4.01375   5.424
> summaryBy(mpg+wt~cyl,data=mtcars,FUN=skew_kurt)
  cyl   mpg.skew mpg.kurtosis    wt.skew wt.kurtosis
1   4  0.2591965   -1.6450119  0.3001262  -1.3559552
2   6 -0.1583137   -1.9069714 -0.2221115  -1.9839377
3   8 -0.3628186   -0.5655154  0.9877797  -0.7126583

利用psych包中的describeBy()分組計算概述統計量：

psych中的describeBy(data,by)，其中data爲目標數據框，by爲分組依據的數據：

> library(psych)
> describeBy(iris[,-5],iris$Species)

 Descriptive statistics by group 
group: setosa
             vars  n mean   sd median trimmed  mad min max range skew kurtosis
Sepal.Length    1 50 5.01 0.35    5.0    5.00 0.30 4.3 5.8   1.5 0.11    -0.45
Sepal.Width     2 50 3.43 0.38    3.4    3.42 0.37 2.3 4.4   2.1 0.04     0.60
Petal.Length    3 50 1.46 0.17    1.5    1.46 0.15 1.0 1.9   0.9 0.10     0.65
Petal.Width     4 50 0.25 0.11    0.2    0.24 0.00 0.1 0.6   0.5 1.18     1.26
               se
Sepal.Length 0.05
Sepal.Width  0.05
Petal.Length 0.02
Petal.Width  0.01
----------------------------------------------------------- 
group: versicolor
             vars  n mean   sd median trimmed  mad min max range  skew
Sepal.Length    1 50 5.94 0.52   5.90    5.94 0.52 4.9 7.0   2.1  0.10
Sepal.Width     2 50 2.77 0.31   2.80    2.78 0.30 2.0 3.4   1.4 -0.34
Petal.Length    3 50 4.26 0.47   4.35    4.29 0.52 3.0 5.1   2.1 -0.57
Petal.Width     4 50 1.33 0.20   1.30    1.32 0.22 1.0 1.8   0.8 -0.03
             kurtosis   se
Sepal.Length    -0.69 0.07
Sepal.Width     -0.55 0.04
Petal.Length    -0.19 0.07
Petal.Width     -0.59 0.03
----------------------------------------------------------- 
group: virginica
             vars  n mean   sd median trimmed  mad min max range  skew
Sepal.Length    1 50 6.59 0.64   6.50    6.57 0.59 4.9 7.9   3.0  0.11
Sepal.Width     2 50 2.97 0.32   3.00    2.96 0.30 2.2 3.8   1.6  0.34
Petal.Length    3 50 5.55 0.55   5.55    5.51 0.67 4.5 6.9   2.4  0.52
Petal.Width     4 50 2.03 0.27   2.00    2.03 0.30 1.4 2.5   1.1 -0.12
             kurtosis   se
Sepal.Length    -0.20 0.09
Sepal.Width      0.38 0.05
Petal.Length    -0.37 0.08
Petal.Width     -0.75 0.04

2.頻數表與列聯表

利用table()生成一維列聯表：

> table(iris$Species)

    setosa versicolor  virginica 
        50         50         50

利用table()生成二維列聯表：

> table(mtcars$cyl,mtcars$am)
   
     0  1
  4  3  8
  6  4  3
  8 12  2

利用xtabs()生成列聯表：

> xtabs(~cyl,data=mtcars)
cyl
 4  6  8 
11  7 14 
> xtabs(~cyl+am,data=mtcars)
   am
cyl  0  1
  4  3  8
  6  4  3
  8 12  2

利用prop.table()將頻數轉化成比例值：

> prop.table(table(iris$Species))

    setosa versicolor  virginica 
 0.3333333  0.3333333  0.3333333 
> prop.table(table(mtcars$cyl,mtcars$am))
   
          0       1
  4 0.09375 0.25000
  6 0.12500 0.09375
  8 0.37500 0.06250

按行或列生成比例：

> #按行生成比例
> prop.table(table(mtcars$cyl,mtcars$am),1)
   
            0         1
  4 0.2727273 0.7272727
  6 0.5714286 0.4285714
  8 0.8571429 0.1428571
> #按列生成比例
> prop.table(table(mtcars$cyl,mtcars$am),2)
   
            0         1
  4 0.1578947 0.6153846
  6 0.2105263 0.2307692
  8 0.6315789 0.1538462

利用addmargins()爲列聯表添加邊際和：

> #添加邊際頻數和
> addmargins(table(mtcars$cyl,mtcars$am))
     
       0  1 Sum
  4    3  8  11
  6    4  3   7
  8   12  2  14
  Sum 19 13  32
> 
> #添加邊際頻率和
> addmargins(prop.table(table(mtcars$cyl,mtcars$am)))
     
            0       1     Sum
  4   0.09375 0.25000 0.34375
  6   0.12500 0.09375 0.21875
  8   0.37500 0.06250 0.43750
  Sum 0.59375 0.40625 1.00000

利用gmodels包中的CrossTable()函數來建立SPSS風格的二維列聯表：

CrossTable()函數有不少選項，能夠作許多事情：計算（行、列、單元格）的百分比；指定小數位數；進行卡方、Fisher和McNemar獨立性檢驗；計算指望和的殘差；將缺失值做爲一種有效值；進行行和列標題的標註；生成SPSS風格的輸出：

> library(gmodels)
> CrossTable(mtcars$cyl,mtcars$am)

 
   Cell Contents
|-------------------------|
|                       N |
| Chi-square contribution |
|           N / Row Total |
|           N / Col Total |
|         N / Table Total |
|-------------------------|

 
Total Observations in Table:  32 

 
             | mtcars$am 
  mtcars$cyl |         0 |         1 | Row Total | 
-------------|-----------|-----------|-----------|
           4 |         3 |         8 |        11 | 
             |     1.909 |     2.790 |           | 
             |     0.273 |     0.727 |     0.344 | 
             |     0.158 |     0.615 |           | 
             |     0.094 |     0.250 |           | 
-------------|-----------|-----------|-----------|
           6 |         4 |         3 |         7 | 
             |     0.006 |     0.009 |           | 
             |     0.571 |     0.429 |     0.219 | 
             |     0.211 |     0.231 |           | 
             |     0.125 |     0.094 |           | 
-------------|-----------|-----------|-----------|
           8 |        12 |         2 |        14 | 
             |     1.636 |     2.391 |           | 
             |     0.857 |     0.143 |     0.438 | 
             |     0.632 |     0.154 |           | 
             |     0.375 |     0.062 |           | 
-------------|-----------|-----------|-----------|
Column Total |        19 |        13 |        32 | 
             |     0.594 |     0.406 |           | 
-------------|-----------|-----------|-----------|

3.獨立性檢驗

3.1 卡方獨立性檢驗

利用chisq.test()來進行卡方獨立性檢驗，p值大於0.05（α）時，獨立：

> library(vcd)
> data("Arthritis")
> 
> #這裏利用with掛載mtcars以簡化變量名的輸入過程
> with(chisq.test(table(Improved,Sex)),data=Arthritis)

    Pearson's Chi-squared test

data:  table(Improved, Sex)
X-squared = 4.8407, df = 2, p-value = 0.08889

3.2 Fisher精確檢驗

利用fisher.test()函數進行Fisher精確檢驗。Fisher精確檢驗的原假設爲：邊界固定的列聯表中行和列相互獨立，p值小於α時拒絕原假設：

> with(fisher.test(table(Improved,Sex)),data=Arthritis)

    Fisher's Exact Test for Count Data

data:  table(Improved, Sex)
p-value = 0.1094
alternative hypothesis: two.sided

4.列聯表相關性的度量

利用vcd包中的assocstats()函數能夠用來計算二維列聯表的phi係數、列聯繫數和Cramer‘s V係數：

> with(assocstats(table(Improved,Treatment)),data=Arthritis)
                    X^2 df  P(> X^2)
Likelihood Ratio 13.530  2 0.0011536
Pearson          13.055  2 0.0014626

Phi-Coefficient   : NA 
Contingency Coeff.: 0.367 
Cramer's V        : 0.394

5.變量間相關性

5.1 Pearson、Spearman和Kendall相關係數

Pearson積差相關係數衡量了兩個定量變量之間的線性相關程度；Spearman等級相關係數則衡量分級定序變量之間的相關程度；Kendall’s Tau相關係數也是一種非參數的等級相關度量。

利用cor(x,use,method)函數來計算這三種相關係數，其中x爲矩陣或數據框；use指定確實數據的處理方式，可選的方式爲all.obs（假設不存在缺失數據，即遇到缺失數據時將報錯）、everything（遇到缺失數據時，相關係數的計算結果將被設爲missing）、complete.obs（行刪除）以及pairwise.complete.obs（成對刪除）;method指定相關係數的類型，可選類型有pearson、spearman或kendall；默認參數爲use='everything'和method='pearson':

> data(iris)
> cor(iris[,-5],method='pearson')
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length    1.0000000  -0.1175698    0.8717538   0.8179411
Sepal.Width    -0.1175698   1.0000000   -0.4284401  -0.3661259
Petal.Length    0.8717538  -0.4284401    1.0000000   0.9628654
Petal.Width     0.8179411  -0.3661259    0.9628654   1.0000000
> cor(iris[,-5],method='spearman')
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length    1.0000000  -0.1667777    0.8818981   0.8342888
Sepal.Width    -0.1667777   1.0000000   -0.3096351  -0.2890317
Petal.Length    0.8818981  -0.3096351    1.0000000   0.9376668
Petal.Width     0.8342888  -0.2890317    0.9376668   1.0000000
> cor(iris[,-5],method='kendall')
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length   1.00000000 -0.07699679    0.7185159   0.6553086
Sepal.Width   -0.07699679  1.00000000   -0.1859944  -0.1571257
Petal.Length   0.71851593 -0.18599442    1.0000000   0.8068907
Petal.Width    0.65530856 -0.15712566    0.8068907   1.0000000

5.2 偏相關係數

偏相關是指在控制一個或多個定量變量時，另外兩個定量變量之間的相互關係；

利用ggm包中的pcor(u,s)函數計算偏相關係數，其中u是一個數值向量，前兩個數值表示要計算相關係數的變量下標，其他的數值爲條件變量（即要排除影響的變量）的下標。s爲變量的協方差陣：

> library(ggm)
> 
> data(iris)
> 
> #以鳶尾花第2,4列數據做爲條件變量來計算第1,3列數據的偏相關係數
> pcor(c(1,3,2,4),cov(iris[,-5]))
[1] 0.7190656

5.3 相關性的顯著性檢驗

在變量間相關性假設檢驗中，經常使用的原假設爲變量間不相關（即整體的相關係數爲0）；

利用cor.test(x,y,alternative,method)對單個的相關係數進行檢驗，其中的x與y爲要檢驗相關性的變量，alternative用來指定進行雙側檢驗或單側檢驗（取值爲'two.side'、'less'、'greater'），而method用以指定要計算的相關類型（'pearson'、'kendall'、'spearman'）：

> #以鳶尾花第2,4列數據做爲條件變量來計算第1,3列數據的偏相關係數
> cor.test(iris[,1],iris[,3])

    Pearson's product-moment correlation

data:  iris[, 1] and iris[, 3]
t = 21.646, df = 148, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.8270363 0.9055080
sample estimates:
      cor 
0.8717538

p值小於0.05，拒絕原假設，斷定變量間不獨立。

利用psych中的corr.test()同時進行多變量間的相關顯著性檢驗：

> library(psych)
> data(iris)
> 
> corr.test(x=iris[,-5],use='complete',method='pearson')
Call:corr.test(x = iris[, -5], use = "complete", method = "pearson")
Correlation matrix 
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length         1.00       -0.12         0.87        0.82
Sepal.Width         -0.12        1.00        -0.43       -0.37
Petal.Length         0.87       -0.43         1.00        0.96
Petal.Width          0.82       -0.37         0.96        1.00
Sample Size 
[1] 150
Probability values (Entries above the diagonal are adjusted for multiple tests.) 
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length         0.00        0.15            0           0
Sepal.Width          0.15        0.00            0           0
Petal.Length         0.00        0.00            0           0
Petal.Width          0.00        0.00            0           0

 To see confidence intervals of the correlations, print with the short=FALSE option

6 t檢驗

檢驗兩組樣本是否有明顯區別，經常使用t檢驗：

6.1 獨立樣本的t檢驗

利用t.test(x,y)來進行t檢驗，原假設是兩個樣原本自同一整體，p值小於α時拒絕原假設：

> t.test(iris[,1],iris[,2])

    Welch Two Sample t-test

data:  iris[, 1] and iris[, 2]
t = 36.463, df = 225.68, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 2.63544 2.93656
sample estimates:
mean of x mean of y 
 5.843333  3.057333 

> 
> t.test(iris[,1],iris[,1]-0.01)

    Welch Two Sample t-test

data:  iris[, 1] and iris[, 1] - 0.01
t = 0.10458, df = 298, p-value = 0.9168
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.1781698  0.1981698
sample estimates:
mean of x mean of y 
 5.843333  5.833333