仅比较某些组组合的均值

Question

数据如下：

> data <- read.csv("data.csv")
> head(data)
  ï..class.1    rev.1 class.2    rev.2
1          7 136.9900    1318  31.9900
2       1223  24.0984    1001   0.0000
3       1318  61.9900    6851 104.2655
4       1014  39.9800    1318  29.9800
5          7  32.9800    7     52.9900
6        291 107.6674     797  31.2741

我想执行显着性检验来比较 rev.1 和 rev.2 的均值，仅在分组 class.1=class.2 的情况下。例如，我试图比较所有“7”类，然后比较所有 1318 个类。我尝试使用 ANOVA 和 TukeyHSD 进行此操作，但仅比较我想要的组时遇到问题。任何指导将不胜感激！

Answer 1

看起来行数据未对齐。如果是这样，那么数据可以分成两组并堆叠在一起。列名是：class、class_id、rev。从那里您可以过滤感兴趣的 class_id，然后继续您的分析。

library(dplyr)
library(tidyr)

# create some data
rev.1 <- rnorm(100, 200,50)
rev.2 <- rnorm(100, 180,35)

class.1 <- seq.int(from = 1000, by = 10, length.out = 100)
class.2 <- seq.int(from = 1000, by = 20, length.out = 100)

df <- tibble(class.1 = class.1, rev.1 = rev.1, class.2 = class.2, rev.2 = rev.2)


# split the data and stack
group_1 <- df %>%
  select(class.1, rev.1) %>%
  gather(key = class,
         value = class_id,
         -rev.1) %>%
  rename(rev = rev.1)

group_2 <- df %>%
  select(class.2, rev.2) %>%
  gather(key = class,
         value = class_id,
         -rev.2) %>%
  rename(rev = rev.2)

df_stacked <- rbind(group_1, group_2)

# filter for the class_id of interest
df_filtered <- df_stacked %>%
  filter(class_id == 1020)

Answer 2

如果您想比较两组的均值，在我看来，t 检验是一个不错的选择。这是使用tidyverse 的选项。首先，我创建了一个名为 dat 的示例数据框。

# Load package
library(tidyverse)

# Set seed
set.seed(12345)

# Create example data frame
dat <- expand.grid(class1 = 1:5, class2 = 1:5) %>%
  slice(rep(1:n(), 5)) %>%
  mutate(rev1 = rnorm(n()), rev2 = rnorm(n())) %>%
  mutate(rev2 = sample(rev2, size = n(), replace = TRUE))
# View the head of data frame
dat
# # A tibble: 125 x 4
#    class1 class2   rev1    rev2
#     <int>  <int>  <dbl>   <dbl>
#  1      1      1  0.586  0.548 
#  2      2      1  0.709  0.868 
#  3      3      1 -0.109  0.0784
#  4      4      1 -0.453 -0.567 
#  5      5      1  0.606 -0.0767
#  6      1      2 -1.82   0.167 
#  7      2      2  0.630  2.66  
#  8      3      2 -0.276  0.831 
#  9      4      2 -0.284 -1.70  
# 10      5      2 -0.919 -2.13  
# # ... with 115 more rows

之后，我过滤了class1 == class2时的数据框，将数据按class1分组，然后使用do函数进行t检验。最后map_dbl可以得到每个t.test的p.value到一个新的数据框。

dat2 <- dat %>%
  filter(class1 == class2) %>%
  group_by(class1) %>%
  do(data_frame(class = .$class1[1],
                TTest = list(t.test(.$rev1, .$rev2)))) %>%
  mutate(PValue = map_dbl(TTest, "p.value"))
dat2
# # A tibble: 5 x 4
# # Groups: class1 [5]
#   class1 class TTest       PValue
#    <int> <int> <list>       <dbl>
# 1      1     1 <S3: htest> 0.700 
# 2      2     2 <S3: htest> 0.381 
# 3      3     3 <S3: htest> 0.859 
# 4      4     4 <S3: htest> 0.0580
# 5      5     5 <S3: htest> 0.206 

如果您想访问特定班级的测试结果，您可以执行以下操作。

# Get the result of the first class
dat2$TTest[dat2$class == 1]
# [[1]]
# 
# Welch Two Sample t-test
# 
# data:  .$rev1 and .$rev2
# t = 0.40118, df = 7.3956, p-value = 0.6996
# alternative hypothesis: true difference in means is not equal to 0
# 95 percent confidence interval:
#   -0.9379329  1.3262368
# sample estimates:
# mean of x mean of y 
# 0.6033533 0.4092013 

这里是另一种选择，我们也可以将数据框拆分为一个列表，并通过列表应用t检验。

# Split the data frame and conduct T-test
dat_list <- dat %>%
  filter(class1 == class2) %>%
  split(.$class1) %>%
  map(~t.test(.$rev1, .$rev2))

# Get the result of the first class
dat_list$`1`

# Welch Two Sample t-test
# 
# data:  .$rev1 and .$rev2
# t = 0.40118, df = 7.3956, p-value = 0.6996
# alternative hypothesis: true difference in means is not equal to 0
# 95 percent confidence interval:
#   -0.9379329  1.3262368
# sample estimates:
# mean of x mean of y 
# 0.6033533 0.4092013