r - PIL 逊相关性 "Stratified"按分类变量

Question

我对 R 很陌生。

我有兴趣为我的数据计算 Pearson 相关性。我已经成功地想出了如何计算数据集中两个连续变量 x 和 y 的相关性；然而，我希望通过第三个分类变量来“分层”相关性:状态。我希望能够说“x 和 y 的相关系数/p 值是 [State] 中的 [Result]。”

我尝试了位于 dplyr 包中的 group_by 方法，该方法位于 cor.test 中(如下所示)。我需要系数和 p 值，所以我一直在尝试使用 cor.test 方法。我也尝试过使用矩阵方法，但也没有成功。

Data<-read.csv("PATHWAYNAME")
  library(dplyr)
  CCor<-cor.test(Data$x, Data$y,
          method=c("pearson"), group_by(State))
  CCor

我能够单独运行每个状态的每组值以获得系数和 p 值；但是，我确信有一种更有效的方法来完成这项任务。我的数据足够大，单独运行它们将不再乏味。

预先感谢您的帮助!

更新:
将此作为一个被极度截断的样本数据集，但同样代表了我自己的变量，我想知道平均收入是否与列出的每个州的访问次数相关；也就是说，平均收入与阿拉巴马州的访问次数呈正相关还是负相关？

>State  NumVis  AvgIncome
>IN       45        60000
>AL       100       56000
>AK       45        80000
>ME       89        54000
>NC       120       100000
>SC       356       43000
>ND       100       25000
>SD       63        20000
>MN       54        46000
>ID       85        55000

使用下面指示的代码运行此数据时，我的结果是这样的:

 CorrDat<-read.csv("File")
     CorrDat %>%
       group_by(State) %>%
        do(tidy(cor.test(CorrDat$NumVis, CorrDat$Income, method="pearson")))

Results

您能否帮助澄清我使用此代码做错了什么，或者我是否需要使用替代方法来完成此任务？

Answer 1

在 R 中有几种方法可以实现这一点。dplyr或更一般地tidyverse是一组能够达到预期结果的流行工具。这些工具的主要区别在于管道 %>%它提供了一种从左到右而不是从内到外编写代码的方法(或在环境中制作一堆中间对象)。尽管管道可以与基础 R 一起使用，但它的受欢迎程度来自 dplyr .

这里有几个关于 mtcars 数据集的例子。关键功能是do和 map这是非常多才多艺的。我建议运行 ?do和 ?map .

library(tidyverse)

mtcars %>%
  group_by(cyl) %>%
  summarize(cor = cor(mpg, disp))
#output
# A tibble: 3 x 2
    cyl correlation
  <dbl>       <dbl>
1     4  -0.8052361
2     6   0.1030827
3     8  -0.5197670

另一种方法是:

mtcars  %>% 
  group_by(cyl) %>%
  do(cor = cor(.$mpg, .$disp)) %>%
  unnest()

或更多变量:

mtcars  %>% 
  group_by(cyl) %>%
  do(cor = as.data.frame(cor(.[,-2]) )) %>%
  unnest() 

cor.test 的示例:

library(broom)

mtcars  %>% 
  group_by(cyl)  %>% 
  do(tidy(cor.test(.$mpg, .$disp))) 
#output
    cyl   estimate  statistic     p.value parameter   conf.low   conf.high                               method alternative
  <dbl>      <dbl>      <dbl>       <dbl>     <int>      <dbl>       <dbl>                               <fctr>      <fctr>
1     4 -0.8052361 -4.0740206 0.002782827         9 -0.9474526 -0.39724826 Pearson's product-moment correlation   two.sided
2     6  0.1030827  0.2317344 0.825929685         5 -0.7046776  0.79446840 Pearson's product-moment correlation   two.sided
3     8 -0.5197670 -2.1075838 0.056774876        12 -0.8232990  0.01492976 Pearson's product-moment correlation   two.sided

还有另一种使用 purrr::map 的方法:

mtcars  %>% 
  split(.$cyl)  %>% 
  map(~cor.test(.x$mpg, .x$disp))

它给出了一个列表，可以使用相同或另一个 map 函数对其进行操作:

mtcars  %>% 
  split(.$cyl)  %>% 
  map(~cor.test(.x$mpg, .x$disp)) %>%
  map_dbl("p.value")
#output:
          4           6           8 
0.002782827 0.825929685 0.056774876 

提取系数:

mtcars  %>% 
  split(.$cyl)  %>% 
  map(~cor.test(.x$mpg, .x$disp)) %>%
  map(~data.frame(cor = .x$estimate, p = .x$p.value)) #check also `map_dfr` and `map_dfc`

#output
$`4`
           cor           p
cor -0.8052361 0.002782827

$`6`
          cor         p
cor 0.1030827 0.8259297

$`8`
          cor          p
cor -0.519767 0.05677488

更新:回答更新的问题:

问题在于您如何指定 do称呼。这是对的:

df %>%
  group_by(State) %>%
  do(tidy(cor.test(.$NumVis, .$AvgIncome, method="pearson")))

哪里 .表示前一个管道传递的数据。在发布的示例中，结果为:

Error in cor.test.default(.$NumVis, .$AvgIncome, method = "pearson") : 
not enough finite observations

考虑到每组仅提供 1 个观察结果，这是合理的

你所做的是:

CorrDat<-read.csv("File")
     CorrDat %>%
       group_by(State) %>%
        do(tidy(cor.test(CorrDat$NumVis, CorrDat$Income, method="pearson")))

将整个 CorrDat 集传递给 do函数，因此它执行与组相同的操作次数。
%>%管道假定传递的数据将用作以下函数中的第一个参数，如果不应该，数据可以被引用为 . .您可以执行类似 .$column 的操作或 .[,2]等等。