R語言第七章基本統(tǒng)計(jì)分析

文 / 李 博

? 描述性統(tǒng)計(jì)分析
? 頻數(shù)表和列聯(lián)表
? 相關(guān)系數(shù)和協(xié)方差
? t檢驗(yàn)
? 非參數(shù)統(tǒng)計(jì)

7.1描述性統(tǒng)計(jì)分析

> myvars<-c("mpg","hp","wt")
> head(mtcars[myvars])#返回mtcars數(shù)據(jù)集中的三類
                   mpg  hp    wt
Mazda RX4         21.0 110 2.620
Mazda RX4 Wag     21.0 110 2.875
Datsun 710        22.8  93 2.320
Hornet 4 Drive    21.4 110 3.215
Hornet Sportabout 18.7 175 3.440
Valiant           18.1 105 3.460
7.1.1
  • 通過 summary() 計(jì)算描述性統(tǒng)計(jì)量。
> myvars<-c("mpg","hp","wt")
> summary(mtcars[myvars])#獲取描述下統(tǒng)計(jì)量:
      mpg              hp              wt       
 Min.   :10.40   Min.   : 52.0   Min.   :1.513  
 1st Qu.:15.43   1st Qu.: 96.5   1st Qu.:2.581  
 Median :19.20   Median :123.0   Median :3.325  
 Mean   :20.09   Mean   :146.7   Mean   :3.217  
 3rd Qu.:22.80   3rd Qu.:180.0   3rd Qu.:3.610  
 Max.   :33.90   Max.   :335.0   Max.   :5.424  
> #包括最小值、最大值、四分位數(shù)和數(shù)值型變量的均值,
> #以及因子向量和邏輯型向量的頻數(shù)統(tǒng)計(jì)。
  • 通過 sapply() 計(jì)算描述性統(tǒng)計(jì)量。
    sapply(x, FUN, options) 其中的 x 是你的數(shù)據(jù)框(或矩陣), FUN 為一個(gè)任意的函數(shù)。
> options(digits = 3)
> mystats<-function(x,na.omit=FALSE){
+   if(na.omit)
+     x<-x[!is.na(x)]#獲得除過缺失值的數(shù)值
+   m<-mean(x)#平均值
+   n<-length(x)#不為空的向量長度
+   s<-sd(x)#均值
+   skew<-sum((x-m)^3/s^3)/n#偏度
+   kurt<-sum((x-m)^4/s^4)/n-3#峰度
+   return(c(n=n,mean=m,stdev=s,skew=skew,kurtosis=kurt))
+ }
> myvars<-c("mpg","hp","wt")
> sapply(mtcars[myvars],mystats)
            mpg      hp      wt
n        32.000  32.000 32.0000
mean     20.091 146.688  3.2172
stdev     6.027  68.563  0.9785
skew      0.611   0.726  0.4231
kurtosis -0.373  -0.136 -0.0227
7.1.2
  • 通過 Hmisc 包中的 describe() 函數(shù)計(jì)算描述性統(tǒng)計(jì)量
> library(Hmisc)
> myvars<-c("mpg","hp","wt")
> describe(mtcars[myvars])#變量和觀測量、缺失值和唯一值得數(shù)目
#平均值、分位數(shù)和五個(gè)最大、五個(gè)最小的值
mtcars[myvars] 

 3  Variables #變量     32  Observations#觀測量
----------------------------------------------------
mpg 
       n  missing distinct     Info     Mean 
      32        0       25    0.999    20.09 
     Gmd      .05      .10      .25      .50 
   6.796    12.00    14.34    15.43    19.20 
     .75      .90      .95 
   22.80    30.09    31.30 

lowest : 10.4 13.3 14.3 14.7 15.0, highest: 26.0 27.3 30.4 32.4 33.9
----------------------------------------------------
hp 
       n  missing distinct     Info     Mean 
      32        0       22    0.997    146.7 
     Gmd      .05      .10      .25      .50 
   77.04    63.65    66.00    96.50   123.00 
     .75      .90      .95 
  180.00   243.50   253.55 

lowest :  52  62  65  66  91, highest: 215 230 245 264 335
----------------------------------------------------
wt 
       n  missing distinct     Info     Mean 
      32        0       29    0.999    3.217 
     Gmd      .05      .10      .25      .50 
   1.089    1.736    1.956    2.581    3.325 
     .75      .90      .95 
   3.610    4.048    5.293 

lowest : 1.513 1.615 1.835 1.935 2.140
highest: 3.845 4.070 5.250 5.345 5.424
----------------------------------------------------
  • 通過 pastecs 包中的 stat.desc() 函數(shù)計(jì)算描述性統(tǒng)計(jì)量。
> options(digits = 3)
> library(pastecs)
> myvars<-c("mpg","hp","wt")
> stat.desc(mtcars[myvars])
                    mpg       hp      wt
nbr.val#所有值       32.00   32.000  32.000
nbr.null#空值         0.00    0.000   0.000
nbr.na#缺失值得數(shù)量  0.00    0.000   0.000
min#最小值           10.40   52.000   1.513
max#最大值           33.90  335.000   5.424
range#值域          23.50  283.000   3.911
sum#總和           642.90 4694.000 102.952
median#中位數(shù)        19.20  123.000   3.325
mean#平均數(shù)          20.09  146.688   3.217
SE.mean#平均數(shù)的標(biāo)準(zhǔn)誤1.07   12.120   0.173
CI.mean.0.95平均數(shù)置信度95%的置信區(qū)間
 2.17   24.720   0.353
var#方差            36.32 4700.867   0.957
std.dev#標(biāo)準(zhǔn)差        6.03   68.563   0.978
coef.var#變異系數(shù)       0.30    0.467   0.304
  • 通過 psych 包中的 describe() 計(jì)算描述性統(tǒng)計(jì)量
> library(psych)
> myvars<-c("mpg","hp","wt")
> describe(mtcars[myvars])
    vars  n   mean    sd median trimmed   mad   min
mpg    1 32  20.09  6.03  19.20   19.70  5.41 10.40
hp     2 32 146.69 68.56 123.00  141.19 77.10 52.00
wt     3 32   3.22  0.98   3.33    3.15  0.77  1.51
       max  range skew kurtosis    se
mpg  33.90  23.50 0.61    -0.37  1.07
hp  335.00 283.00 0.73    -0.14 12.12
wt    5.42   3.91 0.42    -0.02  0.17

統(tǒng)計(jì)量結(jié)果包括:非缺失值的數(shù)量、平均數(shù)、標(biāo)準(zhǔn)差、中位數(shù)、截尾均值、絕對(duì)中位差、最小值、最大值、值域、偏度、峰度和平均值得標(biāo)準(zhǔn)誤。(忒豐富了,有點(diǎn)變態(tài)哈)

7.1.3 分組計(jì)算描述性統(tǒng)計(jì)量
  • 使用 aggregate() 分組獲取描述性統(tǒng)計(jì)量
> myvars<-c("mpg","hp","wt")
> aggregate(mtcars[myvars],by=list(am=mtcars$am),mean)#返回平均值
  am  mpg  hp   wt
1  0   17.1 160 3.77
2  1   24.4 127 2.41
> aggregate(mtcars[myvars],by=list(am=mtcars$am),sd)#返回標(biāo)準(zhǔn)差
  am  mpg   hp    wt
1  0   3.83 53.9 0.777
2  1   6.17 84.1 0.617

其中l(wèi)ist(am=mtcars$am)的作用是為第一列(am)指定一個(gè)標(biāo)簽,若是無指定,則是如下這樣:

> myvars<-c("mpg","hp","wt")
> aggregate(mtcars[myvars],by=list(mtcars$am),mean)
  Group.1  mpg  hp   wt
1       0  17.1 160 3.77
2       1  24.4 127 2.41
> aggregate(mtcars[myvars],by=list(mtcars$am),sd)
  Group.1  mpg   hp    wt
1       0  3.83 53.9 0.777
2       1  6.17 84.1 0.617

  • 使用 by() 分組計(jì)算描述性統(tǒng)計(jì)量
    by(data, INDICES, FUN),其中data是一個(gè)數(shù)據(jù)框或者矩陣,INDICSE是一個(gè)因子或因子組成的列表,定義了分組。FUN是任意函數(shù)。
> mystats <- function(x, na.omit=FALSE){
+   if (na.omit)
+     x <- x[!is.na(x)]
+   m <- mean(x)
+   n <- length(x)
+   s <- sd(x)
+   skew <- sum((x-m)^3/s^3)/n
+   kurt <- sum((x-m)^4/s^4)/n - 3
+   return(c(n=n, mean=m, stdev=s, skew=skew, kurtosis=kurt))
+ }
> dstats<-function(x)sapply(x,mystats)
> myvars<-c("mpg","hp","wt")
> by(mtcars[myvars],mtcars$am,dstats)
mtcars$am: 0
            mpg       hp     wt
n        19.000  19.0000 19.000
mean     17.147 160.2632  3.769
stdev     3.834  53.9082  0.777
skew      0.014  -0.0142  0.976
kurtosis -0.803  -1.2097  0.142
--------------------------------------------------------------------- 
mtcars$am: 1
             mpg      hp     wt
n        13.0000  13.000 13.000
mean     24.3923 126.846  2.411
stdev     6.1665  84.062  0.617
skew      0.0526   1.360  0.210
kurtosis -1.4554   0.563 -1.174

7.1.4 分組計(jì)算的擴(kuò)展
  • 使用 doBy 包中的 summaryBy() 分組計(jì)算概述統(tǒng)計(jì)量
> library(doBy)
> summaryBy(mpg+hp+wt~am,data=mtcars,FUN=mystats)
  am mpg.n mpg.mean mpg.stdev mpg.skew mpg.kurtosis hp.n hp.mean hp.stdev hp.skew
1  0    19     17.1      3.83   0.0140       -0.803   19     160     53.9 -0.0142
2  1    13     24.4      6.17   0.0526       -1.455   13     127     84.1  1.3599
  hp.kurtosis wt.n wt.mean wt.stdev wt.skew wt.kurtosis
1      -1.210   19    3.77    0.777   0.976       0.142
2       0.563   13    2.41    0.617   0.210      -1.174

其中左側(cè)的變量是需要分析的數(shù)值型變量(mpg+hp+wt),右側(cè)(am)是類別型的分組變量。

  • 使用 psych 包中的 describeBy() 分組計(jì)算概述統(tǒng)計(jì)量
> library(psych)
> myvars<-c("mpg","hp","wt")
> describeBy(mtcars[myvars],list(am=mtcars$am))

 Descriptive statistics by group 
am: 0
    vars  n   mean    sd median trimmed   mad   min    max  range  skew kurtosis    se
mpg    1 19  17.15  3.83  17.30   17.12  3.11 10.40  24.40  14.00  0.01    -0.80  0.88
hp     2 19 160.26 53.91 175.00  161.06 77.10 62.00 245.00 183.00 -0.01    -1.21 12.37
wt     3 19   3.77  0.78   3.52    3.75  0.45  2.46   5.42   2.96  0.98     0.14  0.18
--------------------------------------------------------------------- 
am: 1
    vars  n   mean    sd median trimmed   mad   min    max  range skew kurtosis    se
mpg    1 13  24.39  6.17  22.80   24.38  6.67 15.00  33.90  18.90 0.05    -1.46  1.71
hp     2 13 126.85 84.06 109.00  114.73 63.75 52.00 335.00 283.00 1.36     0.56 23.31
wt     3 13   2.41  0.62   2.32    2.39  0.68  1.51   3.57   2.06 0.21    -1.17  0.17

其中describeBy() 函數(shù)不允許指定任意函數(shù)。

7.2 頻數(shù)表和列聯(lián)表(主要是類比型變量,而非描述性統(tǒng)計(jì)量)

> library(vcd)
載入需要的程輯包:grid
Warning message:
程輯包‘vcd’是用R版本3.3.3 來建造的 
> head(Arthritis)
  ID Treatment  Sex Age Improved
1 57   Treated Male  27     Some
2 46   Treated Male  29     None
3 77   Treated Male  30     None
4 17   Treated Male  32   Marked
5 36   Treated Male  46   Marked
6 23   Treated Male  58   Marked

本節(jié)中的數(shù)據(jù)來自 vcd 包中的 Arthritis 數(shù)據(jù)集,上例為前幾個(gè)觀測。

7.2.1 生成頻數(shù)表

table(var1, var2, ..., varN)使用 N 個(gè)類別型變量(因子)創(chuàng)建一個(gè) N 維列聯(lián)表
xtabs(formula, data)根據(jù)一個(gè)公式和一個(gè)矩陣或數(shù)據(jù)框創(chuàng)建一個(gè) N 維列聯(lián)表
prop.table(table, margins)依 margins 定義的邊際列表將表中條目表示為分?jǐn)?shù)形式
margin.table(table, margins)依 margins 定義的邊際列表計(jì)算表中條目的和
addmargins(table, margins)將概述邊 margins (默認(rèn)是求和結(jié)果)放入表中
ftable(table)創(chuàng)建一個(gè)緊湊的“平鋪”式列聯(lián)表

1. 一維列聯(lián)表

> #一維列聯(lián)表
> mytable<-with(Arthritis,table(Improved))
> mytable
Improved
  None   Some Marked 
    42     14     28 

轉(zhuǎn)化為比例值和百分比:

> prop.table(mytable)
Improved
  None   Some Marked 
 0.500  0.167  0.333 
> prop.table(mytable)*100
Improved
  None   Some Marked 
  50.0   16.7   33.3 

2.二維列聯(lián)表
對(duì)于二維列聯(lián)表,table()函數(shù)的使用格式為:mytable<-table(A,B),其中的 A 是行變量, B 是列變量。
使用mytable <- xtabs(~ A + B, data=mydata)時(shí),mydata是一個(gè)矩陣或數(shù)據(jù)框。

> mytable<-xtabs(~Treatment+Improved,data=Arthritis)
> mytable#對(duì)于Arthitis數(shù)據(jù)集,創(chuàng)建一個(gè)自己想要內(nèi)容的表,然后再利用其它函數(shù)進(jìn)行處理
         Improved
Treatment None Some Marked
  Placebo   29    7      7
  Treated   13    7     21

生成邊際頻數(shù)和比例。

>margin.table(mytable,1)#下標(biāo)1指代table()語句中的第一個(gè)變量Treatment
Treatment
Placebo Treated 
     43      41 
> prop.table(mytable,1)# 下標(biāo) 1 指代 table() 語句中的第一個(gè)變量Treatment
         Improved
Treatment  None  Some Marked
  Placebo 0.674 0.163  0.163
  Treated 0.317 0.171  0.512
> margin.table(mytable,2)
Improved
  None   Some Marked 
    42     14     28 
> prop.table(mytable,2)#下標(biāo)2指的是table中的Improved
         Improved
Treatment None Some Marked
  Placebo 0.69 0.50   0.25
  Treated 0.31 0.50   0.75

各單元所占比例如下:

> prop.table(mytable)
         Improved
Treatment   None   Some Marked
  Placebo 0.3452 0.0833 0.0833
  Treated 0.1548 0.0833 0.2500

使用 addmargins() 函數(shù)為這些表格添加邊際和。默認(rèn)情況是為表中所有變量創(chuàng)建邊際和。

> addmargins(mytable)
         Improved
Treatment None Some Marked Sum
  Placebo   29    7      7  43
  Treated   13    7     21  41
  Sum       42   14     28  84
> addmargins(prop.table(mytable))
         Improved
Treatment   None   Some Marked    Sum
  Placebo 0.3452 0.0833 0.0833 0.5119
  Treated 0.1548 0.0833 0.2500 0.4881
  Sum     0.5000 0.1667 0.3333 1.0000

僅添加各行/各列的和:

> addmargins(prop.table(mytable,1),2)
         Improved
Treatment  None  Some Marked   Sum
  Placebo 0.674 0.163  0.163 1.000
  Treated 0.317 0.171  0.512 1.000
> addmargins(prop.table(mytable,2),1)
         Improved
Treatment None Some Marked
  Placebo 0.69 0.50   0.25
  Treated 0.31 0.50   0.75
  Sum     1.00 1.00   1.00

使用 CrossTable 生成二維列聯(lián)表

> library(gmodels)
> CrossTable(Arthritis$Treatment,Arthritis$Improved)

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

 
Total Observations in Table:  84 

 
                    | Arthritis$Improved 
Arthritis$Treatment |      None |      Some |    Marked | Row Total | 
--------------------|-----------|-----------|-----------|-----------|
            Placebo |        29 |         7 |         7 |        43 | 
                    |     2.616 |     0.004 |     3.752 |           | 
                    |     0.674 |     0.163 |     0.163 |     0.512 | 
                    |     0.690 |     0.500 |     0.250 |           | 
                    |     0.345 |     0.083 |     0.083 |           | 
--------------------|-----------|-----------|-----------|-----------|
            Treated |        13 |         7 |        21 |        41 | 
                    |     2.744 |     0.004 |     3.935 |           | 
                    |     0.317 |     0.171 |     0.512 |     0.488 | 
                    |     0.310 |     0.500 |     0.750 |           | 
                    |     0.155 |     0.083 |     0.250 |           | 
--------------------|-----------|-----------|-----------|-----------|
       Column Total |        42 |        14 |        28 |        84 | 
                    |     0.500 |     0.167 |     0.333 |           | 
--------------------|-----------|-----------|-----------|-----------|

3.多維列聯(lián)表
table()和xtabs()都可以基于三個(gè)或更多的類別型變量生成多維列聯(lián)表。
margin.table()、prop.table()和addmargins()函數(shù)可以自然地推廣到高于二維的情況。

> mytable<-xtabs(~Treatment+Sex+Improved,data = Arthritis)
> mytable#生成各單元格的頻數(shù)
, , Improved = None
         Sex
Treatment Female Male
  Placebo     19   10
  Treated      6    7
, , Improved = Some
         Sex
Treatment Female Male
  Placebo      7    0
  Treated      5    2
, , Improved = Marked
         Sex
Treatment Female Male
  Placebo      6    1
  Treated     16    5

邊際頻數(shù)。

> margin.table(mytable,1)#對(duì)使用公式 ~Treatement+Sex+Improve 創(chuàng)建的表,
#Treatment需要下標(biāo)1來引用
Treatment
Placebo Treated 
     43      41 
> margin.table(mytable,2)#Sex需要下標(biāo)2來引用
Sex
Female   Male 
    59     25 
> margin.table(mytable,3)#Improved需要下標(biāo)3來引用
Improved
  None   Some Marked 
    42     14     28 

治療情況( Treatment ) × 改善情況( Improved )的邊際頻數(shù)

> margin.table(mytable,c(1,3))
         Improved
Treatment None Some Marked
  Placebo   29    7      7
  Treated   13    7     21

治療情況( Treatment ) × 性別( Sex )的各類改善情況比例

> ftable(prop.table(mytable,c(1,2)))
                         Improved   None   Some Marked
Treatment Sex                                 
Placebo   Female          0.5938 0.2188 0.1875
          Male            0.9091 0.0000 0.0909
Treated   Female          0.2222 0.1852 0.5926
          Male            0.5000 0.1429 0.3571
> ftable(addmargins(prop.table(mytable,c(1,2)),3))#以較為緊湊的方式輸出多維列聯(lián)表,
#并添加來邊際和。
                 Improved   None   Some Marked    Sum
Treatment Sex                                        
Placebo   Female          0.5938 0.2188 0.1875 1.0000
          Male            0.9091 0.0000 0.0909 1.0000
Treated   Female          0.2222 0.1852 0.5926 1.0000
          Male            0.5000 0.1429 0.3571 1.0000
7.2.2 獨(dú)立性檢驗(yàn)(對(duì)生成列聯(lián)表中的變量是否相關(guān)或獨(dú)立進(jìn)行檢驗(yàn))

** 1. 卡方獨(dú)立性檢驗(yàn)**(卡方檢驗(yàn)就是統(tǒng)計(jì)樣本的實(shí)際觀測值與理論推斷值之間的偏離程度,實(shí)際觀測值與理論推斷值之間的偏離程度就決定卡方值的大小,卡方值越大,越不符合;卡方值越小,偏差越小,越趨于符合,若兩個(gè)值完全相等時(shí),卡方值就為0,表明理論值完全符合。)

> #卡方獨(dú)立性檢驗(yàn)
> library(vcd)
> mytable<-xtabs(~Treatment+Improved,data = Arthritis)
> chisq.test(mytable)#對(duì)二維列表的行變量和列變量進(jìn)行檢驗(yàn)

    Pearson's Chi-squared test

data:  mytable
X-squared = 13.055, df = 2, p-value = 0.001463
> #本例中對(duì)治療情況和改善情況進(jìn)行檢驗(yàn),P<0.01說明存在一定的相關(guān)性。
> mytable<-xtabs(~Improved+Sex,data=Arthritis)
> chisq.test(mytable)

    Pearson's Chi-squared test
data:  mytable
X-squared = 4.8407, df = 2, p-value = 0.08889
Warning message:
In chisq.test(mytable) : Chi-squared近似算法有可能不準(zhǔn)(信息量太少,且有一個(gè)小于5的值)

由于P=0.089>0.01,故說明治療結(jié)果和性別之間是不獨(dú)立的。

2. Fisher精確檢驗(yàn)(邊界固定的列聯(lián)表中行和列是相互獨(dú)立的)

> mytable<-xtabs(~Treatment+Improved,data=Arthritis)
> fisher.test(mytable)

    Fisher's Exact Test for Count Data

data:  mytable
p-value = 0.001393
alternative hypothesis: two.sided#替代假設(shè)

3. Cochran-Mantel-Haenszel檢驗(yàn)(兩個(gè)名義變量在第三個(gè)變量的每一層中都是條件獨(dú)立的)

> #兩個(gè)名義變量在第三個(gè)變量的每一層中都是條件獨(dú)立的
> mytable<-xtabs(~Treatment+Improved+Sex,data = Arthritis)
> mantelhaen.test(mytable)

    Cochran-Mantel-Haenszel test

data:  mytable
Cochran-Mantel-Haenszel M^2 = 14.632, df = 2,
p-value = 0.0006647
#表明:分性別來看,用藥治療的患者較接受安慰劑的患者有了更多的改善
7.2.3 相關(guān)性的度量(數(shù)值越大表明相關(guān)性越強(qiáng))

二維列聯(lián)表的相關(guān)性度量(很是抽象)

> #二維列聯(lián)表的相關(guān)性度量
> library(vcd)
> mytable<-xtabs(~Treatment+Improved,data = Arthritis)
> assocstats(mytable)
                    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 

7.3 相關(guān)

7.3.1 相關(guān)的類型

1. Pearson、Spearman和Kendall相關(guān)
使用的函數(shù)為:cor(x, use= , method= ),以及cov()函數(shù)。
協(xié)方差和相關(guān)系數(shù)

> #協(xié)方差和相關(guān)系數(shù)
> options(digits = 3)
> states<-state.x77[,1:6]#使用R基礎(chǔ)安裝中的數(shù)據(jù)集,提供了大量的關(guān)于人口、收入、
#文盲率、預(yù)期壽命、
#謀殺率和高中畢業(yè)率數(shù)據(jù)。
> cov(states)#計(jì)算樣本協(xié)方差(所有變量之間兩兩計(jì)算相關(guān))
           Population Income Illiteracy Life Exp  Murder  HS Grad
Population   19931684 571230    292.868 -407.842 5663.52 -3551.51
Income         571230 377573   -163.702  280.663 -521.89  3076.77
Illiteracy        293   -164      0.372   -0.482    1.58    -3.24
Life Exp         -408    281     -0.482    1.802   -3.87     6.31
Murder           5664   -522      1.582   -3.869   13.63   -14.55
HS Grad         -3552   3077     -3.235    6.313  -14.55    65.24
> cor(states)# 計(jì)算了Pearson積差相關(guān)系數(shù)
           Population Income Illiteracy Life Exp Murder HS Grad
Population     1.0000  0.208      0.108  -0.0681  0.344 -0.0985
Income         0.2082  1.000     -0.437   0.3403 -0.230  0.6199
Illiteracy     0.1076 -0.437      1.000  -0.5885  0.703 -0.6572
Life Exp      -0.0681  0.340     -0.588   1.0000 -0.781  0.5822
Murder         0.3436 -0.230      0.703  -0.7808  1.000 -0.4880
HS Grad       -0.0985  0.620     -0.657   0.5822 -0.488  1.0000
> cor(states,method = "spearman")# 計(jì)算了Spearman等級(jí)相關(guān)系數(shù)
           Population Income Illiteracy Life Exp Murder HS Grad
Population      1.000  0.125      0.313   -0.104  0.346  -0.383
Income          0.125  1.000     -0.315    0.324 -0.217   0.510
Illiteracy      0.313 -0.315      1.000   -0.555  0.672  -0.655
Life Exp       -0.104  0.324     -0.555    1.000 -0.780   0.524
Murder          0.346 -0.217      0.672   -0.780  1.000  -0.437
HS Grad        -0.383  0.510     -0.655    0.524 -0.437   1.000

從表中可以看到收入和高中畢業(yè)率之間存在很強(qiáng)的正相關(guān)。Pearson相關(guān)關(guān)系為0.6199,Spearman相關(guān)關(guān)系為0.510,兩者均為極強(qiáng)的正相關(guān),但預(yù)期壽命和文盲率存在很強(qiáng)的負(fù)相關(guān)。Pearson相關(guān)關(guān)系為-0.588,Spearman相關(guān)關(guān)系為-0.555。

2. 偏相關(guān)(偏相關(guān)是指在控制一個(gè)或多個(gè)定量變量時(shí),另外兩個(gè)定量變量之間的相互關(guān)系)
使用pcor(u, S)函數(shù)。

> library(ggm)
> colnames(states)
[1] "Population" "Income"     "Illiteracy" "Life Exp"   "Murder"     "HS Grad"   
> pcor(c(1,5,2,3,6),cov(states))#前兩個(gè)1,5表示要計(jì)算相關(guān)系數(shù)
#的變量下標(biāo),其余數(shù)值為條件變量的下標(biāo)。
[1] 0.346#人口和謀殺率之間的相關(guān)系數(shù)0.346
7.3.2 相關(guān)性的顯著性檢驗(yàn)

使用cor.test(x, y, alternative = , method = )函數(shù),其中的 x 和 y 為要檢驗(yàn)相關(guān)性的變量, alternative 則用來指定進(jìn)行雙側(cè)檢驗(yàn)或單側(cè)檢驗(yàn)(取值為 "two.side" 、 "less" 或 "greater" ),而 method 用以指定要計(jì)算的相關(guān)類型( "pearson" 、"kendall" 或 "spearman" )。

  • 檢驗(yàn)?zāi)撤N相關(guān)系數(shù)的顯著性
> #檢驗(yàn)?zāi)撤N相關(guān)系數(shù)的顯著性
> cor.test(states[,3],states[,5])

    Pearson's product-moment correlation

data:  states[, 3] and states[, 5]
t = 7, df = 50, p-value = 1e-08
alternative hypothesis: true correlation is not equal to 0#替代假設(shè),
95 percent confidence interval:#95%的置信區(qū)間
 0.528 0.821
sample estimates:#樣本估計(jì)
  cor 
0.703 
  • 通過 corr.test 計(jì)算相關(guān)矩陣并進(jìn)行顯著性檢驗(yàn)
> library(psych)
> corr.test(states,use = "complete")
Call:corr.test(x = states, use = "complete")
Correlation matrix 
           Population Income Illiteracy Life Exp Murder HS Grad
Population       1.00   0.21       0.11    -0.07   0.34   -0.10
Income           0.21   1.00      -0.44     0.34  -0.23    0.62
Illiteracy       0.11  -0.44       1.00    -0.59   0.70   -0.66
Life Exp        -0.07   0.34      -0.59     1.00  -0.78    0.58
Murder           0.34  -0.23       0.70    -0.78   1.00   -0.49
HS Grad         -0.10   0.62      -0.66     0.58  -0.49    1.00
Sample Size 
[1] 50
Probability values (Entries above the diagonal are adjusted for multiple tests.) 
           Population Income Illiteracy Life Exp Murder HS Grad
Population       0.00   0.59       1.00      1.0   0.10       1
Income           0.15   0.00       0.01      0.1   0.54       0
Illiteracy       0.46   0.00       0.00      0.0   0.00       0
Life Exp         0.64   0.02       0.00      0.0   0.00       0
Murder           0.01   0.11       0.00      0.0   0.00       0
HS Grad          0.50   0.00       0.00      0.0   0.00       0

 To see confidence intervals of the correlations, print with the short=FALSE option
最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點(diǎn),簡書系信息發(fā)布平臺(tái),僅提供信息存儲(chǔ)服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容